Add advanced plots to IBMA report (#864)

* Add Summary Stats and N True Voxels advanced plots

* Save raw_data to use in Reports

* Add dropdown advanced menu for percentage of voxels and summary stats

* Update figures.py

* Update base.py

* Update figures.py

* Update ibma.py

* Update figures.py

* Update figures.py

* Update figures.py

* Add y-axis ticks label

* save inputs_

* Scipy latest only support Python >=3.9. Build the documentation with 3.9

* [skip ci] test documentation build in python 3.9

* Update figures.py

Author: JulioAPeraza

.readthedocs.yaml

@@ -8,7 +8,7 @@ version: 2
 build:
   os: "ubuntu-22.04"
   tools:
-    python: "3.8"
+    python: "3.9"
 
 # Build documentation in the docs/ directory with Sphinx
 sphinx:

nimare/meta/ibma.py

@@ -114,8 +114,7 @@ def _preprocess_input(self, dataset):
                 # Mask required input images using either the dataset's mask or the estimator's.
                 temp_arr = masker.transform(img4d)
 
-                data = masker.transform(img4d)
-                temp_image_inputs[name] = data
+                temp_image_inputs[name] = temp_arr
                 if self.aggressive_mask:
                     # An intermediate step to mask out bad voxels.
                     # Can be dropped once PyMARE is able to handle masked arrays or missing data.
@@ -132,8 +131,8 @@ def _preprocess_input(self, dataset):
                             good_voxels_bool,
                         )
                 else:
-                    self.inputs_[name] = data  # This data is saved only to use in Reports
-                    data_bags = zip(*_apply_liberal_mask(data))
+                    self.inputs_[name] = temp_arr
+                    data_bags = zip(*_apply_liberal_mask(temp_arr))
 
                     keys = ["values", "voxel_mask", "study_mask"]
                     self.inputs_["data_bags"][name] = [dict(zip(keys, bag)) for bag in data_bags]
@@ -153,6 +152,8 @@ def _preprocess_input(self, dataset):
             for name, raw_masked_data in temp_image_inputs.items():
                 self.inputs_[name] = raw_masked_data[:, self.inputs_["aggressive_mask"]]
 
+        self.inputs_["raw_data"] = temp_image_inputs  # This data is saved only to use in Reports
+
 
 class Fishers(IBMAEstimator):
     """An image-based meta-analytic test using t- or z-statistic images.

nimare/reports/base.py

@@ -35,6 +35,8 @@
     _plot_dof_map,
     _plot_relcov_map,
     _plot_ridgeplot,
+    _plot_sumstats,
+    _plot_true_voxels,
     gen_table,
     plot_clusters,
     plot_coordinates,
@@ -378,11 +380,17 @@ def __init__(self, out_dir, config=None):
             ext = "".join(src.suffixes)
             desc_text = config.get("caption")
             iframe = config.get("iframe", False)
+            dropdown = config.get("dropdown", False)
 
             contents = None
             html_anchor = src.relative_to(out_dir)
             if ext == ".html":
                 contents = IFRAME_SNIPPET.format(html_anchor) if iframe else src.read_text()
+                if dropdown:
+                    contents = (
+                        f"<details><summary>Advanced ({self.title})</summary>{contents}</details>"
+                    )
+                    self.title = ""
             elif ext == ".png":
                 contents = PNG_SNIPPET.format(html_anchor)
 
@@ -475,16 +483,19 @@ def __init__(
                 )
             elif meta_type == "IBMA":
                 # Use "z_maps", for Fishers, and Stouffers; otherwise use "beta_maps".
-                key_maps = "z_maps" if "z_maps" in self.results.estimator.inputs_ else "beta_maps"
-                maps_arr = self.results.estimator.inputs_[key_maps]
+                key_maps = (
+                    "z_maps"
+                    if "z_maps" in self.results.estimator.inputs_["raw_data"]
+                    else "beta_maps"
+                )
+                maps_arr = self.results.estimator.inputs_["raw_data"][key_maps]
                 ids_ = self.results.estimator.inputs_["id"]
                 x_label = "Z" if key_maps == "z_maps" else "Beta"
 
                 if self.results.estimator.aggressive_mask:
                     _plot_relcov_map(
                         maps_arr,
                         self.results.estimator.masker,
-                        self.results.estimator.inputs_["aggressive_mask"],
                         self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-relcov.png",
                     )
                 else:
@@ -494,13 +505,25 @@ def __init__(
                         self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-dof.png",
                     )
 
+                _plot_true_voxels(
+                    maps_arr,
+                    ids_,
+                    self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-truevoxels.html",
+                )
+
                 _plot_ridgeplot(
                     maps_arr,
                     ids_,
                     x_label,
                     self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-ridgeplot.html",
                 )
 
+                _plot_sumstats(
+                    maps_arr,
+                    ids_,
+                    self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-summarystats.html",
+                )
+
                 similarity_table = _compute_similarities(maps_arr, ids_)
 
                 plot_heatmap(

nimare/reports/default.yml

@@ -55,8 +55,14 @@ sections:
     title: Relative Coverage Map
   - bids: {value: preliminary, dset: 1, suffix: figure-dof}
     title: DoF Map
+  - bids: {value: preliminary, dset: 1, suffix: figure-truevoxels}
+    title: Percentage of Valid Voxels
+    dropdown: True
   - bids: {value: preliminary, dset: 1, suffix: figure-ridgeplot}
     title: Ridge Plot
+  - bids: {value: preliminary, dset: 1, suffix: figure-summarystats}
+    title: Summary Statistics
+    dropdown: True
   - bids: {value: preliminary, dset: 1, suffix: figure-similarity}
     title: Similarity of Input Maps
   - bids: {value: preliminary, dset: 2, suffix: summary}

nimare/reports/figures.py

@@ -15,10 +15,9 @@
     view_img,
 )
 from ridgeplot import ridgeplot
+from scipy import stats
 from scipy.cluster.hierarchy import leaves_list, linkage, optimal_leaf_ordering
 
-from nimare.utils import _boolean_unmask
-
 TABLE_STYLE = [
     dict(
         selector="th, td",
@@ -402,6 +401,60 @@ def plot_clusters(img, out_filename):
     fig.close()
 
 
+def _plot_true_voxels(maps_arr, ids_, out_filename):
+    """Plot percentage of valid voxels.
+
+    .. versionadded:: 0.2.2
+
+    """
+    n_studies, n_voxels = maps_arr.shape
+    mask = ~np.isnan(maps_arr) & (maps_arr != 0)
+
+    x_label, y_label = "Voxels Included", "ID"
+    perc_voxs = mask.sum(axis=1) / n_voxels
+    valid_df = pd.DataFrame({y_label: ids_, x_label: perc_voxs})
+    valid_sorted_df = valid_df.sort_values(x_label, ascending=True)
+
+    fig = px.bar(
+        valid_sorted_df,
+        x=x_label,
+        y=y_label,
+        orientation="h",
+        color=x_label,
+        color_continuous_scale="blues",
+        range_color=(0, 1),
+    )
+
+    fig.update_xaxes(
+        showline=True,
+        linewidth=2,
+        linecolor="black",
+        visible=True,
+        showticklabels=False,
+        title=None,
+    )
+    fig.update_yaxes(
+        showline=True,
+        linewidth=2,
+        linecolor="black",
+        visible=True,
+        ticktext=valid_sorted_df[y_label].str.slice(0, MAX_CHARS).tolist(),
+    )
+
+    height = n_studies * PXS_PER_STD
+    fig.update_layout(
+        height=height,
+        autosize=True,
+        font_size=14,
+        plot_bgcolor="white",
+        xaxis_gridcolor="white",
+        yaxis_gridcolor="white",
+        xaxis_gridwidth=2,
+        showlegend=False,
+    )
+    fig.write_html(out_filename, full_html=True, include_plotlyjs=True)
+
+
 def _plot_ridgeplot(maps_arr, ids_, x_label, out_filename):
     """Plot histograms of the images.
 
@@ -446,7 +499,88 @@ def _plot_ridgeplot(maps_arr, ids_, x_label, out_filename):
     fig.write_html(out_filename, full_html=True, include_plotlyjs=True)
 
 
-def _plot_relcov_map(maps_arr, masker, aggressive_mask, out_filename):
+def _plot_sumstats(maps_arr, ids_, out_filename):
+    """Plot summary statistics of the images.
+
+    .. versionadded:: 0.2.2
+
+    """
+    n_studies = len(ids_)
+    mask = ~np.isnan(maps_arr) & (maps_arr != 0)
+    maps_lst = [maps_arr[i][mask[i]] for i in range(n_studies)]
+
+    stats_lbls = [
+        "Mean",
+        "STD",
+        "Var",
+        "Median",
+        "Mode",
+        "Min",
+        "Max",
+        "Skew",
+        "Kurtosis",
+        "Range",
+        "Moment",
+        "IQR",
+    ]
+    scores, id_lst = [], []
+    for id_, map_ in zip(ids_, maps_lst):
+        scores.extend(
+            [
+                np.mean(map_),
+                np.std(map_),
+                np.var(map_),
+                np.median(map_),
+                stats.mode(map_)[0],
+                np.min(map_),
+                np.max(map_),
+                stats.skew(map_),
+                stats.kurtosis(map_),
+                np.max(map_) - np.min(map_),
+                stats.moment(map_, moment=4),
+                stats.iqr(map_),
+            ]
+        )
+        id_lst.extend([id_] * len(stats_lbls))
+
+    stats_labels = stats_lbls * n_studies
+    data_df = pd.DataFrame({"ID": id_lst, "Score": scores, "Stat": stats_labels})
+
+    fig = px.strip(
+        data_df,
+        y="Score",
+        color="ID",
+        facet_col="Stat",
+        stripmode="group",
+        facet_col_wrap=4,
+        facet_col_spacing=0.08,
+    )
+
+    fig.update_xaxes(showline=True, linewidth=2, linecolor="black", mirror=True)
+    fig.update_yaxes(
+        constrain="domain",
+        matches=None,
+        showline=True,
+        linewidth=2,
+        linecolor="black",
+        mirror=True,
+        title=None,
+    )
+    fig.update_layout(
+        height=900,
+        autosize=True,
+        font_size=14,
+        plot_bgcolor="white",
+        xaxis_gridcolor="white",
+        yaxis_gridcolor="white",
+        xaxis_gridwidth=2,
+        showlegend=False,
+    )
+    fig.for_each_yaxis(lambda yaxis: yaxis.update(showticklabels=True))
+    fig.write_html(out_filename, full_html=True, include_plotlyjs=True)
+
+
+def _plot_relcov_map(maps_arr, masker, out_filename):
     """Plot relative coverage map.
 
     .. versionadded:: 0.2.0
@@ -460,8 +594,6 @@ def _plot_relcov_map(maps_arr, masker, aggressive_mask, out_filename):
     binary_maps_arr = np.where((-epsilon > maps_arr) | (maps_arr > epsilon), 1, 0)
     coverage_arr = np.sum(binary_maps_arr, axis=0) / binary_maps_arr.shape[0]
 
-    # Add bad voxels back to the arr to transform it back to an image
-    coverage_arr = _boolean_unmask(coverage_arr, aggressive_mask)
     coverage_img = masker.inverse_transform(coverage_arr)
 
     # Plot coverage map

setup.cfg

@@ -77,7 +77,7 @@ doc =
     sphinx>=3.5
     sphinx-argparse
     sphinx-copybutton
-    sphinx_gallery==0.10.1
+    sphinx-gallery
     sphinx_rtd_theme>=1.3.0
     sphinxcontrib-bibtex
 tests =