Add IBMAWorkflow

docs/api.rst

@@ -336,6 +336,7 @@ For more information about fetching data from the internet, see :ref:`fetching t
    workflows.macm_workflow
    workflows.cbma.CBMAWorkflow
    workflows.cbma.PairwiseCBMAWorkflow
+   workflows.ibma.IBMAWorkflow
 
 :mod:`nimare.reports`: NiMARE report
 --------------------------------------------------

examples/02_meta-analyses/12_plot_ibma_workflow.py

@@ -0,0 +1,114 @@
+"""
+
+.. _ibma_workflow:
+
+====================================================
+Run an image-based meta-analysis (IBMA) workflow
+====================================================
+
+NiMARE provides a plethora of tools for performing meta-analyses on neuroimaging data.
+Sometimes it's difficult to know where to start, especially if you're new to meta-analysis.
+This tutorial will walk you through using a IBMA workflow function which puts together
+the fundamental steps of a IBMA meta-analysis.
+"""
+import os
+from pathlib import Path
+
+import matplotlib.pyplot as plt
+from nilearn.plotting import plot_stat_map
+
+from nimare.extract import download_nidm_pain
+
+###############################################################################
+# Download data
+# -----------------------------------------------------------------------------
+
+dset_dir = download_nidm_pain()
+
+###############################################################################
+# Load Dataset
+# -----------------------------------------------------------------------------
+import nibabel as nib
+from nilearn.image import resample_to_img
+
+from nimare.dataset import Dataset
+from nimare.transforms import ImageTransformer
+from nimare.utils import get_resource_path
+
+dset_file = os.path.join(get_resource_path(), "nidm_pain_dset.json")
+dset = Dataset(dset_file)
+dset.update_path(dset_dir)
+
+# Resample images
+imgs = dset.get({"z_maps": ("image", "z")})["z_maps"]
+mask_img = dset.masker.mask_img
+for img in imgs:
+    res_img = resample_to_img(nib.load(img), mask_img, clip=True, interpolation="linear")
+    res_img.to_filename(img)
+
+# Calculate missing images
+xformer = ImageTransformer(target=["varcope", "z"])
+dset = xformer.transform(dset)
+
+###############################################################################
+# Run IBMA Workflow
+# -----------------------------------------------------------------------------
+# The fit method of a IBMA workflow class runs the following steps:
+#
+# 1. Runs a meta-analysis using the specified method (default: Stouffers)
+# 2. Applies a corrector to the meta-analysis results (default: FDRCorrector, indep)
+# 3. Generates cluster tables and runs diagnostics on the corrected results (default: Jackknife)
+#
+# All in one call!
+from nimare.workflows.ibma import IBMAWorkflow
+
+workflow = IBMAWorkflow()
+result = workflow.fit(dset)
+
+###############################################################################
+# Plot Results
+# -----------------------------------------------------------------------------
+# The fit method of the IBMA workflow class returns a :class:`~nimare.results.MetaResult` object,
+# where you can access the corrected results of the meta-analysis and diagnostics tables.
+#
+# Corrected map:
+img = result.get_map("z_corr-FDR_method-indep")
+plot_stat_map(
+    img,
+    cut_coords=4,
+    display_mode="z",
+    threshold=1.65,  # voxel_thresh p < .05, one-tailed
+    cmap="RdBu_r",
+    vmax=4,
+)
+plt.show()
+
+###############################################################################
+# Clusters table
+# ``````````````````````````````````````````````````````````````````````````````
+result.tables["z_corr-FDR_method-indep_tab-clust"]
+
+###############################################################################
+# Contribution table
+# ``````````````````````````````````````````````````````````````````````````````
+result.tables["z_corr-FDR_method-indep_diag-Jackknife_tab-counts_tail-positive"]
+
+###############################################################################
+# Report
+# -----------------------------------------------------------------------------
+# Finally, a NiMARE report is generated from the MetaResult.
+from nimare.reports.base import run_reports
+
+# root_dir = Path(os.getcwd()).parents[1] / "docs" / "_build"
+# Use the previous root to run the documentation locally.
+root_dir = Path(os.getcwd()).parents[1] / "_readthedocs"
+html_dir = root_dir / "html" / "auto_examples" / "02_meta-analyses" / "12_plot_ibma_workflow"
+html_dir.mkdir(parents=True, exist_ok=True)
+
+run_reports(result, html_dir)
+
+####################################
+# .. raw:: html
+#
+#     <iframe src="./12_plot_ibma_workflow/report.html" style="border:none;" seamless="seamless"\
+#        width="100%" height="1000px"></iframe>

nimare/reports/base.py

@@ -30,14 +30,14 @@
 from pkg_resources import resource_filename as pkgrf
 
 from nimare.meta.cbma.base import CBMAEstimator, PairwiseCBMAEstimator
-from nimare.meta.ibma import IBMAEstimator
 from nimare.reports.figures import (
     gen_table,
     plot_clusters,
     plot_coordinates,
     plot_heatmap,
     plot_interactive_brain,
     plot_mask,
+    plot_ridgeplot,
     plot_static_brain,
 )
 
@@ -54,6 +54,12 @@
     "method": "Method",
     "alpha": "Alpha",
     "prior": "Prior",
+    "use_sample_size": "Use sample size for weights",
+    "beta": "Parameter estimate",
+    "se": "Standard error of the parameter estimate",
+    "varcope": "Variance of the parameter estimate",
+    "t": "T-statistic",
+    "z": "Z-statistic",
 }
 
 PNG_SNIPPET = """\
@@ -71,11 +77,7 @@
 
 SUMMARY_TEMPLATE = """\
 <ul class="elem-desc">
-<li>Number of studies: {n_studies:d}</li>
-<li>Number of experiments: {n_exps:d}</li>
-<li>Number of experiments included: {n_exps_sel:d}</li>
-<li>Number of foci: {n_foci:d} </li>
-<li>Number of foci outside the mask: {n_foci_nonbrain:d} </li>
+{meta_text}
 </ul>
 <details>
 <summary>Experiments excluded</summary><br />
@@ -86,7 +88,7 @@
 ESTIMATOR_TEMPLATE = """\
 <ul class="elem-desc">
 <li>Estimator: {est_name}</li>
-<li>Kernel Transformer: {kernel_transformer}{ker_params_text}</li>
+{ker_text}
 {est_params_text}
 </ul>
 """
@@ -109,26 +111,100 @@
 """
 
 
-def _gen_est_summary(obj, out_filename):
-    """Generate html with parameter use in obj (e.g., estimator)."""
-    params_dict = obj.get_params()
-    est_params = {k: v for k, v in params_dict.items() if not k.startswith("kernel_transformer")}
-    ker_params = {k: v for k, v in params_dict.items() if k.startswith("kernel_transformer__")}
+def _get_cbma_summary(dset, sel_ids):
+    n_studies = len(dset.coordinates["study_id"].unique())
+
+    mask = dset.masker.mask_img
+    sel_ids = dset.get_studies_by_mask(mask)
+    sel_dset = dset.slice(sel_ids)
+
+    n_foci = dset.coordinates.shape[0]
+    n_foci_sel = sel_dset.coordinates.shape[0]
+    n_foci_nonbrain = n_foci - n_foci_sel
+
+    n_exps = len(dset.ids)
+    n_exps_sel = len(sel_ids)
+
+    cbma_text = [
+        f"<li>Number of studies: {n_studies:d}</li>",
+        f"<li>Number of experiments: {n_exps:d}</li>",
+        f"<li>Number of experiments included: {n_exps_sel:d}</li>",
+        f"<li>Number of foci: {n_foci:d} </li>",
+        f"<li>Number of foci outside the mask: {n_foci_nonbrain:d} </li>",
+    ]
+
+    return " ".join(cbma_text)
+
+
+def _get_ibma_summary(dset, sel_ids):
+    img_df = dset.images
+    n_studies = len(img_df["study_id"].unique())
+
+    ignore_columns = ["id", "study_id", "contrast_id"]
+    map_type = [c for c in img_df if not c.endswith("__relative") and c not in ignore_columns]
+
+    n_imgs = len(dset.ids)
+    n_sel_ids = len(sel_ids)
+
+    ibma_text = [
+        f"<li>Number of studies: {n_studies:d}</li>",
+        f"<li>Number of images: {n_imgs:d}</li>",
+        f"<li>Number of images included: {n_sel_ids:d}</li>",
+    ]
+
+    maptype_text = ["<li>Available maps: ", "<ul>"]
+    maptype_text.extend(f"<li>{PARAMETERS_DICT[m]} ({m})</li>" for m in map_type)
+    maptype_text.extend(["</ul>", "</li>"])
+
+    ibma_text.extend(maptype_text)
+    return " ".join(ibma_text)
+
+
+def _gen_summary(dset, sel_ids, meta_type, out_filename):
+    """Generate preliminary checks from dataset for the report."""
+    exc_ids = list(set(dset.ids) - set(sel_ids))
+    exc_ids_str = ", ".join(exc_ids)
+
+    meta_text = (
+        _get_cbma_summary(dset, sel_ids)
+        if meta_type == "CBMA"
+        else _get_ibma_summary(dset, sel_ids)
+    )
+
+    summary_text = SUMMARY_TEMPLATE.format(
+        meta_text=meta_text,
+        exc_ids=exc_ids_str,
+    )
+    (out_filename).write_text(summary_text, encoding="UTF-8")
+
 
+def _get_kernel_summary(params_dict):
+    kernel_transformer = str(params_dict["kernel_transformer"])
+    ker_params = {k: v for k, v in params_dict.items() if k.startswith("kernel_transformer__")}
     ker_params_text = ["<ul>"]
     ker_params_text.extend(f"<li>{PARAMETERS_DICT[k]}: {v}</li>" for k, v in ker_params.items())
     ker_params_text.append("</ul>")
     ker_params_text = "".join(ker_params_text)
 
+    return f"<li>Kernel Transformer: {kernel_transformer}{ker_params_text}</li>"
+
+
+def _gen_est_summary(obj, out_filename):
+    """Generate html with parameter use in obj (e.g., estimator)."""
+    params_dict = obj.get_params()
+
+    # Add kernel transformer parameters to summary if obj is a CBMAEstimator
+    ker_text = _get_kernel_summary(params_dict) if isinstance(obj, CBMAEstimator) else ""
+
+    est_params = {k: v for k, v in params_dict.items() if not k.startswith("kernel_transformer")}
     est_params_text = [f"<li>{PARAMETERS_DICT[k]}: {v}</li>" for k, v in est_params.items()]
     est_params_text = "".join(est_params_text)
 
     est_name = obj.__class__.__name__
 
     summary_text = ESTIMATOR_TEMPLATE.format(
         est_name=est_name,
-        kernel_transformer=str(params_dict["kernel_transformer"]),
-        ker_params_text=ker_params_text,
+        ker_text=ker_text,
         est_params_text=est_params_text,
     )
     (out_filename).write_text(summary_text, encoding="UTF-8")
@@ -188,33 +264,6 @@ def _gen_fig_summary(img_key, threshold, out_filename):
     (out_filename).write_text(summary_text, encoding="UTF-8")
 
 
-def _gen_summary(dset, out_filename):
-    """Generate preliminary checks from dataset for the report."""
-    mask = dset.masker.mask_img
-    n_studies = len(dset.coordinates["study_id"].unique())
-    sel_ids = dset.get_studies_by_mask(mask)
-    sel_dset = dset.slice(sel_ids)
-
-    n_foci = dset.coordinates.shape[0]
-    n_foci_sel = sel_dset.coordinates.shape[0]
-    n_foci_nonbrain = n_foci - n_foci_sel
-
-    n_exps = len(dset.ids)
-    n_exps_sel = len(sel_dset.ids)
-    exc_ids = list(set(dset.ids) - set(sel_dset.ids))
-    exc_ids_str = ", ".join(exc_ids)
-
-    summary_text = SUMMARY_TEMPLATE.format(
-        n_studies=n_studies,
-        n_exps=n_exps,
-        n_exps_sel=n_exps_sel,
-        n_foci=n_foci,
-        n_foci_nonbrain=n_foci_nonbrain,
-        exc_ids=exc_ids_str,
-    )
-    (out_filename).write_text(summary_text, encoding="UTF-8")
-
-
 def _gen_figures(results, img_key, diag_name, threshold, fig_dir):
     """Generate html and png objects for the report."""
     # Plot brain images if not empty
@@ -362,6 +411,7 @@ def __init__(
         out_filename="report.html",
     ):
         self.results = results
+        meta_type = "CBMA" if issubclass(type(self.results.estimator), CBMAEstimator) else "IBMA"
         self._is_pairwise_estimator = issubclass(
             type(self.results.estimator), PairwiseCBMAEstimator
         )
@@ -374,31 +424,51 @@ def __init__(
         self.fig_dir = self.out_dir / "figures"
         self.fig_dir.mkdir(parents=True, exist_ok=True)
 
-        datasets = (
-            [self.results.estimator.dataset1, self.results.estimator.dataset2]
-            if self._is_pairwise_estimator
-            else [self.results.estimator.dataset]
-        )
-        for dset_i, dataset in enumerate(datasets):
+        if self._is_pairwise_estimator:
+            datasets = [self.results.estimator.dataset1, self.results.estimator.dataset2]
+            sel_ids = [
+                self.results.estimator.inputs_["id1"],
+                self.results.estimator.inputs_["id2"],
+            ]
+        else:
+            datasets = [self.results.estimator.dataset]
+            sel_ids = [self.results.estimator.inputs_["id"]]
+
+        for dset_i, (dataset, sel_id) in enumerate(zip(datasets, sel_ids)):
             # Generate summary text
-            _gen_summary(dataset, self.fig_dir / f"preliminary_dset-{dset_i+1}_summary.html")
+            _gen_summary(
+                dataset,
+                sel_id,
+                meta_type,
+                self.fig_dir / f"preliminary_dset-{dset_i+1}_summary.html",
+            )
 
             # Plot mask
             plot_mask(
                 dataset.masker.mask_img,
                 self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-mask.png",
             )
 
-            if issubclass(type(self.results.estimator), CBMAEstimator):
+            if meta_type == "CBMA":
                 # Plot coordinates for CBMA estimators
                 plot_coordinates(
                     dataset.coordinates,
                     self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-static.png",
                     self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-interactive.html",
                     self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-legend.png",
                 )
-            elif issubclass(type(self.results.estimator), IBMAEstimator):
-                raise NotImplementedError
+            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]
+                ids_ = self.results.estimator.inputs_["id"]
+                x_label = "Z" if key_maps == "z_maps" else "Beta"
+                plot_ridgeplot(
+                    maps_arr,
+                    ids_,
+                    x_label,
+                    self.fig_dir / f"preliminary_dset-{dset_i+1}_figure-ridgeplot.png",
+                )
 
         _gen_est_summary(self.results.estimator, self.fig_dir / "estimator_summary.html")
         _gen_cor_summary(self.results.corrector, self.fig_dir / "corrector_summary.html")

nimare/reports/default.yml

@@ -51,6 +51,8 @@ sections:
   - bids: {value: preliminary, dset: 1, suffix: figure-interactive}
     title: Explorer
     iframe: True
+  - bids: {value: preliminary, dset: 1, suffix: figure-ridgeplot}
+    title: Ridge Plot
   - bids: {value: preliminary, dset: 2, suffix: summary}
     title: Dataset 2
   - bids: {value: preliminary, dset: 2, suffix: figure-mask}