ADA-SVR (3/4) PR example of models

doc_conf/api.rst

@@ -16,6 +16,7 @@ Functions
    :toctree: generated/
 
    ada_svr
+   ada_svr_pvalue
    aggregate_quantiles
    clustered_inference
    data_simulation

doc_conf/conf.py

@@ -90,7 +90,7 @@
 # built documents.
 #
 # The short X.Y version.
-from hidimstat._version import __version__  # noqa
+from hidimstat import __version__
 
 # The full version, including alpha/beta/rc tags.
 release = __version__
@@ -217,9 +217,9 @@
     "python": ("https://docs.python.org/3", None),
     "numpy": ("https://numpy.org/devdocs", None),
     "scipy": ("https://scipy.github.io/devdocs", None),
-    "matplotlib": ("https://matplotlib.org", None),
+    "matplotlib": ("https://matplotlib.org/stable/", None),
     "sklearn": ("https://scikit-learn.org/stable", None),
-    "numba": ("https://numba.pydata.org/numba-doc/latest", None),
+    "numba": ("https://numba.readthedocs.io/en/stable/", None),
     "joblib": ("https://joblib.readthedocs.io/en/latest", None),
     "pandas": ("https://pandas.pydata.org/pandas-docs/stable", None),
     "seaborn": ("https://seaborn.pydata.org/", None),
@@ -228,7 +228,6 @@
 
 examples_dirs = ["../examples"]
 gallery_dirs = ["auto_examples"]
-import mne
 
 scrapers = ("matplotlib",)
 try:
@@ -240,6 +239,7 @@
     pass
 if any(x in scrapers for x in ("pyvista")):
     from traits.api import push_exception_handler
+    import mne
 
     push_exception_handler(reraise_exceptions=True)
     report_scraper = mne.report._ReportScraper()
@@ -259,6 +259,8 @@
     "abort_on_example_error": False,
     "image_scrapers": scrapers,
     "show_memory": True,
+    'filename_pattern': r'\.py',
+    'ignore_pattern': r'__init__\.py',
     # 'reference_url': {
     #     'numpy': 'http://docs.scipy.org/doc/numpy-1.9.1',
     #     'scipy': 'http://docs.scipy.org/doc/scipy-0.17.0/reference',

doc_conf/index.rst

@@ -46,8 +46,7 @@ is also needed to install ``pytest``.
 Documentation & Examples
 ------------------------
 
-Documentation about the main HiDimStat functions is available
-`here <api.html>`_ and examples are available `here <auto_examples/index.html>`_.
+Documentation of HiDimStat is composed of an `API <api.html>`_ and `examples <auto_examples/index.html>`_.
 
 As of now, there are three different examples (Python scripts) that
 illustrate how to use the main HiDimStat functions.
@@ -118,15 +117,15 @@ Application to source localization (MEG/EEG data):
 
 Single/Group statistically validated importance using conditional permutations:
 
-* Chamma, A., Thirion, B., & Engemann, D. (2024). **Variable importance in
-high-dimensional settings requires grouping**. In Proceedings of the 38th
-Conference of the Association for the Advancement of Artificial
-Intelligence(AAAI 2024), Vancouver, Canada.
+* Chamma, A., Thirion, B., & Engemann, D. (2024). Variable importance in
+  high-dimensional settings requires grouping. In Proceedings of the 38th
+  Conference of the Association for the Advancement of Artificial
+  Intelligence(AAAI 2024), Vancouver, Canada.
 
-* Chamma, A., Engemann, D., & Thirion, B. (2023). **Statistically Valid Variable
-Importance Assessment through Conditional Permutations**. In Proceedings of the
-37th Conference on Neural Information Processing Systems (NeurIPS 2023), New
-Orleans, USA.
+* Chamma, A., Engemann, D., & Thirion, B. (2023). Statistically Valid Variable
+  Importance Assessment through Conditional Permutations. In Proceedings of the
+  37th Conference on Neural Information Processing Systems (NeurIPS 2023), New
+  Orleans, USA.
 
 If you use our packages, we would appreciate citations to the relevant
 aforementioned papers.

doc_conf/references.bib

@@ -177,4 +177,28 @@ @article{liuFastPowerfulConditional2021
   archiveprefix = {arxiv},
   keywords = {Statistics - Methodology},
   file = {/home/ahmad/Zotero/storage/8HRQZX3H/Liu et al. - 2021 - Fast and Powerful Conditional Randomization Testin.pdf;/home/ahmad/Zotero/storage/YFNDKN2B/2006.html}
-}
+}
+
+@article{gaonkar_deriving_2012,
+	title = {Deriving statistical significance maps for {SVM} based image classification and group comparisons},
+	volume = {15},
+	url = {https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3703958/},
+	abstract = {Population based pattern analysis and classification for quantifying structural and functional differences between diverse groups has been shown to be a powerful tool for the study of a number of diseases, and is quite commonly used especially in neuroimaging. The alternative to these pattern analysis methods, namely mass univariate methods such as voxel based analysis and all related methods, cannot detect multivariate patterns associated with group differences, and are not particularly suitable for developing individual-based diagnostic and prognostic biomarkers. A commonly used pattern analysis tool is the support vector machine ({SVM}). Unlike univariate statistical frameworks for morphometry, analytical tools for statistical inference are unavailable for the {SVM}. In this paper, we show that null distributions ordinarily obtained by permutation tests using {SVMs} can be analytically approximated from the data. The analytical computation takes a small fraction of the time it takes to do an actual permutation test, thereby rendering it possible to quickly create statistical significance maps derived from {SVMs}. Such maps are critical for understanding imaging patterns of group differences and interpreting which anatomical regions are important in determining the classifier's decision.},
+	pages = {723--730},
+	number = {0},
+	journaltitle = {Medical image computing and computer-assisted intervention : {MICCAI} ... International Conference on Medical Image Computing and Computer-Assisted Intervention},
+	journal = {Med Image Comput Comput Assist Interv},
+	author = {Gaonkar, Bilwaj and Davatzikos, Christos},
+	urldate = {2024-12-16},
+	year = {2012},
+	pmid = {23285616},
+	pmcid = {PMC3703958},
+	file = {PubMed Central Full Text PDF:/home/likusch/Zotero/storage/DX8QQAF5/Gaonkar and Davatzikos - 2012 - Deriving statistical significance maps for SVM based image classification and group comparisons.pdf:application/pdf},
+}
+
+@book{molnar2020interpretable,
+  title={Interpretable machine learning},
+  author={Molnar, Christoph},
+  year={2020},
+  publisher={Lulu. com}
+}

examples/README.txt

@@ -5,4 +5,4 @@ Examples Gallery
 
 .. contents:: Contents
    :local:
-   :depth: 3
+   :depth: 0

examples/_utils/plot_dataset.py

@@ -0,0 +1,246 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from matplotlib.colors import LogNorm
+from scipy import stats
+
+
+def plot_dataset1D(X, y, beta, title="Toy dataset"):
+    """
+    Plot a 1D toy dataset with the true regression line.
+    Parameters:
+    - X: Input features (n_samples, n_features)
+    - y: Vectors of the regression  (n_samples,)
+    - beta: Coefficients of the true variable of importance (n_features,)
+    - title: Title of the plot (str)
+
+    Returns:
+    - a figure with 3 subplots
+    """
+    # Create a figure and a set of subplots
+    fig, ([ax11, ax12], [ax21, ax22]) = plt.subplots(
+        2, 2, width_ratios=[0.9, 0.01], height_ratios=[0.9, 0.1], figsize=(6.4, 4.8)
+    )
+    # plot data
+    im_X = ax11.imshow(X, aspect="auto", interpolation="nearest")
+    cbaxes_X = fig.add_axes([0.05, 0.275, 0.03, 0.6])
+    col_X = plt.colorbar(im_X, cax=cbaxes_X, location="left")
+    col_X.ax.set_xlabel("X values", loc="center", labelpad=10)
+    col_X.ax.xaxis.set_label_position("top")
+    ax11.set_ylabel("n samples")
+    ax11.set_xlabel("n features")
+    ax11.set_title("X:data", fontdict={"fontweight": "bold"})
+    # plot regression
+    im_Y = ax12.imshow(np.expand_dims(y, 1), aspect="auto", interpolation="nearest")
+    ax12.set_ylabel("y:regression        ", fontdict={"fontweight": "bold"})
+    ax12.yaxis.tick_right()
+    ax12.set_xticks([])
+    cbaxes_Y = fig.add_axes([0.95, 0.275, 0.03, 0.6])
+    col_Y = plt.colorbar(im_Y, cax=cbaxes_Y, location="left")
+    col_Y.ax.set_xlabel("Y values   ", loc="center", labelpad=10)
+    col_Y.ax.xaxis.set_label_position("top")
+    # plot beta
+    ax21.imshow(np.expand_dims(beta, 0), aspect="auto", interpolation="nearest")
+    ax21.set_xlabel("beta:variable of importance", fontdict={"fontweight": "bold"})
+    ax21.set_yticks([])
+    ax22.axis("off")
+
+    plt.suptitle(title)
+    plt.subplots_adjust(hspace=0.3, left=0.15, right=0.85)
+
+
+def plot_validate_variable_importance(beta, beta_hat, vmin=0.0, vmax=1.0):
+    """
+    Plot for validating of the variable importance estimation.
+    Parameters:
+    - beta: Coefficients of the true variable of importance (n_features,)
+    - beta_hat: Coefficients of the estimated variable of importance (n_features,)
+    - vmin: Minimum value of the colorbar (float)
+    - vmax: Maximum value of the colorbar (float)
+
+    Returns:
+    - a figure with 2 subplots
+    """
+    fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(6.4, 2.8))
+    # plot beta
+    ax1.imshow(np.expand_dims(beta, 0), vmin=vmin, vmax=vmax)
+    ax1.set_xlabel("beta:variable of importance", fontdict={"fontweight": "bold"})
+    ax1.set_yticks([])
+    # plot beta_hat
+    im = ax2.imshow(np.expand_dims(beta_hat, 0), vmin=vmin, vmax=vmax)
+    ax2.set_xlabel("beta hat:variable of importance", fontdict={"fontweight": "bold"})
+    ax2.set_yticks([])
+    plt.colorbar(
+        im, ax=ax2, orientation="horizontal", label="Variable importance", pad=0.5
+    )
+
+    plt.subplots_adjust(top=1.0, bottom=0.2)
+
+
+def plot_pvalue_H0(
+    beta_hat,
+    pvalue,
+    pvalue_corrected,
+    one_minus_pvalue,
+    one_minus_pvalue_corrected,
+    vmin=0.0,
+    vmax=1.0,
+):
+    """
+    Plot for the confidence in the hipothse that the variables are important.
+    Parameters:
+    - beta_hat: Coefficients of the estimated variable of importance (n_features,)
+    - pvalue: pvalue of each variable of importance (n_features,)
+    - pvalue_corrected: corrected pvalue of each variable of importance (n_features,)
+    - one_minus_pvalue: 1 - pvalue of each variable of importance (n_features,)
+    - one_minus_pvalue_corrected: 1 - corrected pvalue of each variable of importance (n_features,)
+    - vmin: Minimum value of the colorbar (float)
+    - vmax: Maximum value of the colorbar (float)
+
+    Returns:
+    - a figure with 3 subplots
+    """
+    fig, (ax1, ax2, ax3) = plt.subplots(3, 1, figsize=(6.4, 4.8))
+    # plot beta_hat
+    im_beta_hat = ax1.imshow(np.expand_dims(beta_hat, 0), vmin=vmin, vmax=vmax)
+    ax1.set_title("beta hat:variable of importance", fontdict={"fontweight": "bold"})
+    ax1.set_yticks([])
+    colbar_beta_hat = plt.colorbar(
+        im_beta_hat,
+        ax=ax1,
+        orientation="horizontal",
+        label="Variable importance",
+        pad=0.2,
+    )
+    colbar_beta_hat.ax.xaxis.labelpad = 0
+    # plot pvalue
+    im_pvalue = ax2.imshow(
+        np.expand_dims(pvalue, 0),
+        norm=LogNorm(vmin=np.min(pvalue), vmax=np.max(pvalue)),
+        cmap=plt.cm.viridis.reversed(),
+    )
+    ax2.set_title(
+        "pvalue of each variable of importance", fontdict={"fontweight": "bold"}
+    )
+    ax2.set_yticks([])
+    colbar_pvalue = plt.colorbar(
+        im_pvalue, ax=ax2, orientation="horizontal", label="pvalue", pad=0.2
+    )
+    colbar_pvalue.ax.xaxis.labelpad = 0
+    # plot pvalue_corrected
+    im_pvalue_corrected = ax3.imshow(
+        np.expand_dims(pvalue_corrected, 0),
+        norm=LogNorm(vmin=np.min(pvalue_corrected), vmax=np.max(pvalue_corrected)),
+        cmap=plt.cm.viridis.reversed(),
+    )
+    ax3.set_title(
+        "corrected pvalue of each variable of importance",
+        fontdict={"fontweight": "bold"},
+    )
+    ax3.set_yticks([])
+    colbar_pvalue_corrected = plt.colorbar(
+        im_pvalue_corrected,
+        ax=ax3,
+        orientation="horizontal",
+        label="corrected pvalue",
+        pad=0.2,
+    )
+    colbar_pvalue_corrected.ax.xaxis.labelpad = 0
+
+    plt.subplots_adjust(top=1.0, hspace=0.3)
+
+
+def plot_pvalue_H1(
+    beta_hat,
+    pvalue,
+    pvalue_corrected,
+    one_minus_pvalue,
+    one_minus_pvalue_corrected,
+    vmin=0.0,
+    vmax=1.0,
+):
+    """
+    Plot for the confidence in the hipothse that the variables are not important.
+    Parameters:
+    - beta_hat: Coefficients of the estimated variable of importance (n_features,)
+    - pvalue: pvalue of each variable of importance (n_features,)
+    - pvalue_corrected: corrected pvalue of each variable of importance (n_features,)
+    - one_minus_pvalue: 1 - pvalue of each variable of importance (n_features,)
+    - one_minus_pvalue_corrected: 1 - corrected pvalue of each variable of importance (n_features,)
+    - vmin: Minimum value of the colorbar (float)
+    - vmax: Maximum value of the colorbar (float)
+
+    Returns:
+    - a figure with 3 subplots
+    """
+    fig, (ax1, ax2, ax3) = plt.subplots(3, 1, figsize=(6.4, 4.8))
+    im_beta_hat = ax1.imshow(np.expand_dims(beta_hat, 0), vmin=vmin, vmax=vmax)
+    ax1.set_title("beta hat:variable of importance", fontdict={"fontweight": "bold"})
+    ax1.set_yticks([])
+    colbar_beta_hat = plt.colorbar(
+        im_beta_hat,
+        ax=ax1,
+        orientation="horizontal",
+        label="Variable importance",
+        pad=0.2,
+    )
+    colbar_beta_hat.ax.xaxis.labelpad = 0
+    im_pvalue = ax2.imshow(
+        np.expand_dims(one_minus_pvalue_corrected, 0),
+        norm=LogNorm(
+            vmin=np.min(one_minus_pvalue), vmax=np.max(one_minus_pvalue_corrected)
+        ),
+        cmap=plt.cm.viridis.reversed(),
+    )
+    ax2.set_title(
+        "pvalue of each variable of importance", fontdict={"fontweight": "bold"}
+    )
+    ax2.set_yticks([])
+    colbar_pvalue = plt.colorbar(
+        im_pvalue, ax=ax2, orientation="horizontal", label="pvalue", pad=0.2
+    )
+    colbar_pvalue.ax.xaxis.labelpad = 0
+    im_pvalue_corrected = ax3.imshow(
+        np.expand_dims(one_minus_pvalue_corrected, 0),
+        norm=LogNorm(
+            vmin=np.min(one_minus_pvalue_corrected),
+            vmax=np.max(one_minus_pvalue_corrected),
+        ),
+        cmap=plt.cm.viridis.reversed(),
+    )
+    ax3.set_title(
+        "corrected pvalue of each variable of importance",
+        fontdict={"fontweight": "bold"},
+    )
+    ax3.set_yticks([])
+    colbar_pvalue_corrected = plt.colorbar(
+        im_pvalue_corrected,
+        ax=ax3,
+        orientation="horizontal",
+        label="corrected pvalue",
+        pad=0.2,
+    )
+    colbar_pvalue_corrected.ax.xaxis.labelpad = 0
+    plt.subplots_adjust(top=1.0, hspace=0.3)
+
+
+def plot_compare_proba_estimated(proba, beta_hat, scale):
+    """
+    Plot the comparison of the estimated probability with the normal distribution.
+    Parameters:
+    - proba: Probability of each coefficent (n_features, n_permutations)
+    - beta_hat: Coefficients of the estimated variable of importance (n_features,)
+    - scale: Standard deviation of the distribution of the coefficients (n_features,)
+
+    Returns:
+    - a figure with 5*5 subplots
+    """
+    fig, axs = plt.subplots(5, 5, figsize=(6.4, 4.8))
+    for index, ax in enumerate(axs.flat):
+        # plot the histogram of the proba
+        ax.hist(
+            proba[:, index], bins=100, density=True, alpha=0.5, color="b", label="proba"
+        )
+        # plot the normal distribution
+        x = np.linspace(-3 * scale[index], +3 * scale[index], 100)
+        ax.plot(x, stats.norm.pdf(x, 0.0, scale[index]), "r-", lw=2, label="norm pdf")
+        ax.set_title(f"beta_hat[{index}]")

examples/methods/README.txt

@@ -0,0 +1,7 @@
+Description of the methods of the package
+=========================================
+The package contains the following methods:
+
+.. contents:: Contents
+   :local:
+   :depth: 0