Merge branch 'main' into flexible_report

Author: mwprestonjr

README.rst

@@ -82,6 +82,7 @@ There are also optional dependencies, which are not required for model fitting i
 
 - `matplotlib <https://github.com/matplotlib/matplotlib>`_ is needed to visualize data and model fits
 - `tqdm <https://github.com/tqdm/tqdm>`_ is needed to print progress bars when fitting many models
+- `pandas <https://github.com/pandas-dev/pandas>`_ is needed to for exporting model fit results to dataframes
 - `pytest <https://github.com/pytest-dev/pytest>`_ is needed to run the test suite locally
 
 We recommend using the `Anaconda <https://www.anaconda.com/distribution/>`_ distribution to manage these requirements.

examples/analyses/plot_mne_example.py

@@ -5,7 +5,7 @@
 Parameterizing neural power spectra with MNE, doing a topographical analysis.
 
 This tutorial requires that you have `MNE <https://mne-tools.github.io/>`_
-installed.
+installed. This tutorial needs mne >= 1.2.
 
 If you don't already have MNE, you can follow instructions to get it
 `here <https://mne-tools.github.io/stable/getting_started.html>`_.
@@ -23,10 +23,7 @@
 
 # Import MNE, as well as the MNE sample dataset
 import mne
-from mne import io
 from mne.datasets import sample
-from mne.viz import plot_topomap
-from mne.time_frequency import psd_welch
 
 # FOOOF imports
 from fooof import FOOOFGroup
@@ -52,16 +49,15 @@
 ###################################################################################################
 
 # Get the data path for the MNE example data
-raw_fname = sample.data_path() + '/MEG/sample/sample_audvis_filt-0-40_raw.fif'
-event_fname = sample.data_path() + '/MEG/sample/sample_audvis_filt-0-40_raw-eve.fif'
+raw_fname = sample.data_path() / 'MEG' / 'sample' / 'sample_audvis_filt-0-40_raw.fif'
 
 # Load the example MNE data
 raw = mne.io.read_raw_fif(raw_fname, preload=True, verbose=False)
 
 ###################################################################################################
 
 # Select EEG channels from the dataset
-raw = raw.pick_types(meg=False, eeg=True, eog=False, exclude='bads')
+raw = raw.pick(['eeg'], exclude='bads')
 
 ###################################################################################################
 
@@ -110,15 +106,16 @@ def check_nans(data, nan_policy='zero'):
 # frequency representations - meaning we have to calculate power spectra.
 #
 # To do so, we will leverage the time frequency tools available with MNE,
-# in the `time_frequency` module. In particular, we can use the ``psd_welch``
-# function, that takes in MNE data objects and calculates and returns power spectra.
+# in the `time_frequency` module. In particular, we can use the ``compute_psd``
+# method, that takes in MNE data objects and calculates and returns power spectra.
 #
 
 ###################################################################################################
 
-# Calculate power spectra across the the continuous data
-spectra, freqs = psd_welch(raw, fmin=1, fmax=40, tmin=0, tmax=250,
-                           n_overlap=150, n_fft=300)
+# Calculate power spectra across the continuous data
+psd = raw.compute_psd(method="welch", fmin=1, fmax=40, tmin=0, tmax=250,
+                      n_overlap=150, n_fft=300)
+spectra, freqs = psd.get_data(return_freqs=True)
 
 ###################################################################################################
 # Fitting Power Spectrum Models
@@ -193,7 +190,7 @@ def check_nans(data, nan_policy='zero'):
 ###################################################################################################
 
 # Plot the topography of alpha power
-plot_topomap(alpha_pw, raw.info, cmap=cm.viridis, contours=0);
+mne.viz.plot_topomap(alpha_pw, raw.info, cmap=cm.viridis, contours=0, size=4)
 
 ###################################################################################################
 #
@@ -214,8 +211,7 @@ def check_nans(data, nan_policy='zero'):
     band_power = check_nans(get_band_peak_fg(fg, band_def)[:, 1])
 
     # Create a topomap for the current oscillation band
-    mne.viz.plot_topomap(band_power, raw.info, cmap=cm.viridis, contours=0,
-                         axes=axes[ind], show=False);
+    mne.viz.plot_topomap(band_power, raw.info, cmap=cm.viridis, contours=0, axes=axes[ind])
 
     # Set the plot title
     axes[ind].set_title(label + ' power', {'fontsize' : 20})
@@ -268,7 +264,7 @@ def check_nans(data, nan_policy='zero'):
 ###################################################################################################
 
 # Plot the topography of aperiodic exponents
-plot_topomap(exps, raw.info, cmap=cm.viridis, contours=0)
+mne.viz.plot_topomap(exps, raw.info, cmap=cm.viridis, contours=0, size=4)
 
 ###################################################################################################
 #
@@ -297,6 +293,3 @@ def check_nans(data, nan_policy='zero'):
 # In this example, we have seen how to apply power spectrum models to data that is
 # managed and processed with MNE.
 #
-
-###################################################################################################
-#

fooof/core/reports.py

@@ -22,7 +22,7 @@
 ###################################################################################################
 
 @check_dependency(plt, 'matplotlib')
-def save_report_fm(fm, file_name, file_path=None, plt_log=True, **plot_kwargs):
+def save_report_fm(fm, file_name, file_path=None, plt_log=False, add_settings=True, **plot_kwargs):
     """Generate and save out a PDF report for a power spectrum model fit.
 
     Parameters
@@ -35,41 +35,46 @@ def save_report_fm(fm, file_name, file_path=None, plt_log=True, **plot_kwargs):
         Path to directory to save to. If None, saves to current directory.
     plt_log : bool, optional, default: False
         Whether or not to plot the frequency axis in log space.
+    add_settings : bool, optional, default: True
+        Whether to add a print out of the model settings to the end of the report.
     plot_kwargs : keyword arguments
         Keyword arguments to pass into the plot method.
     """
 
+    # Define grid settings based on what is to be plotted
+    n_rows = 3 if add_settings else 2
+    height_ratios = [0.5, 1.0, 0.25] if add_settings else [0.45, 1.0]
+
     # Set up outline figure, using gridspec
     _ = plt.figure(figsize=REPORT_FIGSIZE)
-    grid = gridspec.GridSpec(3, 1, height_ratios=[0.45, 1.0, 0.25])
+    grid = gridspec.GridSpec(n_rows, 1, hspace=0.25, height_ratios=height_ratios)
 
     # First - text results
     ax0 = plt.subplot(grid[0])
     results_str = gen_results_fm_str(fm)
     ax0.text(0.5, 0.7, results_str, REPORT_FONT, ha='center', va='center')
     ax0.set_frame_on(False)
-    ax0.set_xticks([])
-    ax0.set_yticks([])
+    ax0.set(xticks=[], yticks=[])
 
     # Second - data plot
     ax1 = plt.subplot(grid[1])
     fm.plot(plt_log=plt_log, ax=ax1, **plot_kwargs)
 
     # Third - FOOOF settings
-    ax2 = plt.subplot(grid[2])
-    settings_str = gen_settings_str(fm, False)
-    ax2.text(0.5, 0.1, settings_str, REPORT_FONT, ha='center', va='center')
-    ax2.set_frame_on(False)
-    ax2.set_xticks([])
-    ax2.set_yticks([])
+    if add_settings:
+        ax2 = plt.subplot(grid[2])
+        settings_str = gen_settings_str(fm, False)
+        ax2.text(0.5, 0.1, settings_str, REPORT_FONT, ha='center', va='center')
+        ax2.set_frame_on(False)
+        ax2.set(xticks=[], yticks=[])
 
     # Save out the report
     plt.savefig(fpath(file_path, fname(file_name, SAVE_FORMAT)))
     plt.close()
 
 
 @check_dependency(plt, 'matplotlib')
-def save_report_fg(fg, file_name, file_path=None):
+def save_report_fg(fg, file_name, file_path=None, add_settings=True):
     """Generate and save out a PDF report for a group of power spectrum models.
 
     Parameters
@@ -80,19 +85,26 @@ def save_report_fg(fg, file_name, file_path=None):
         Name to give the saved out file.
     file_path : str, optional
         Path to directory to save to. If None, saves to current directory.
+    add_settings : bool, optional, default: True
+        Whether to add a print out of the model settings to the end of the report.
     """
 
+    # Define grid settings based on what is to be plotted
+    n_rows = 4 if add_settings else 3
+    height_ratios = [1.0, 1.0, 1.0, 0.5] if add_settings else [0.8, 1.0, 1.0]
+
     # Initialize figure
     _ = plt.figure(figsize=REPORT_FIGSIZE)
-    grid = gridspec.GridSpec(3, 2, wspace=0.4, hspace=0.25, height_ratios=[0.8, 1.0, 1.0])
+    grid = gridspec.GridSpec(n_rows, 2, wspace=0.4, hspace=0.25, height_ratios=height_ratios)
 
     # First / top: text results
     ax0 = plt.subplot(grid[0, :])
     results_str = gen_results_fg_str(fg)
     ax0.text(0.5, 0.7, results_str, REPORT_FONT, ha='center', va='center')
     ax0.set_frame_on(False)
-    ax0.set_xticks([])
-    ax0.set_yticks([])
+    ax0.set(xticks=[], yticks=[])
+
+    # Second - data plots
 
     # Aperiodic parameters plot
     ax1 = plt.subplot(grid[1, 0])
@@ -106,6 +118,14 @@ def save_report_fg(fg, file_name, file_path=None):
     ax3 = plt.subplot(grid[2, :])
     plot_fg_peak_cens(fg, ax3)
 
+    # Third - Model settings
+    if add_settings:
+        ax4 = plt.subplot(grid[3, :])
+        settings_str = gen_settings_str(fg, False)
+        ax4.text(0.5, 0.1, settings_str, REPORT_FONT, ha='center', va='center')
+        ax4.set_frame_on(False)
+        ax4.set(xticks=[], yticks=[])
+
     # Save out the report
     plt.savefig(fpath(file_path, fname(file_name, SAVE_FORMAT)))
     plt.close()

fooof/core/strings.py

@@ -36,7 +36,7 @@ def gen_width_warning_str(freq_res, bwl):
     output = '\n'.join([
         '',
         'FOOOF WARNING: Lower-bound peak width limit is < or ~= the frequency resolution: ' + \
-            '{:1.2f} <= {:1.2f}'.format(freq_res, bwl),
+            '{:1.2f} <= {:1.2f}'.format(bwl, freq_res),
         '\tLower bounds below frequency-resolution have no effect ' + \
         '(effective lower bound is the frequency resolution).',
         '\tToo low a limit may lead to overfitting noise as small bandwidth peaks.',

fooof/data/conversions.py

@@ -0,0 +1,106 @@
+"""Conversion functions for organizing model results into alternate representations."""
+
+import numpy as np
+
+from fooof import Bands
+from fooof.core.funcs import infer_ap_func
+from fooof.core.info import get_ap_indices, get_peak_indices
+from fooof.core.modutils import safe_import, check_dependency
+from fooof.analysis.periodic import get_band_peak
+
+pd = safe_import('pandas')
+
+###################################################################################################
+###################################################################################################
+
+def model_to_dict(fit_results, peak_org):
+    """Convert model fit results to a dictionary.
+
+    Parameters
+    ----------
+    fit_results : FOOOFResults
+        Results of a model fit.
+    peak_org : int or Bands
+        How to organize peaks.
+        If int, extracts the first n peaks.
+        If Bands, extracts peaks based on band definitions.
+
+    Returns
+    -------
+    dict
+        Model results organized into a dictionary.
+    """
+
+    fr_dict = {}
+
+    # aperiodic parameters
+    for label, param in zip(get_ap_indices(infer_ap_func(fit_results.aperiodic_params)),
+                            fit_results.aperiodic_params):
+        fr_dict[label] = param
+
+    # periodic parameters
+    peaks = fit_results.peak_params
+
+    if isinstance(peak_org, int):
+
+        if len(peaks) < peak_org:
+            nans = [np.array([np.nan] * 3) for ind in range(peak_org-len(peaks))]
+            peaks = np.vstack((peaks, nans))
+
+        for ind, peak in enumerate(peaks[:peak_org, :]):
+            for pe_label, pe_param in zip(get_peak_indices(), peak):
+                fr_dict[pe_label.lower() + '_' + str(ind)] = pe_param
+
+    elif isinstance(peak_org, Bands):
+        for band, f_range in peak_org:
+            for label, param in zip(get_peak_indices(), get_band_peak(peaks, f_range)):
+                fr_dict[band + '_' + label.lower()] = param
+
+    # goodness-of-fit metrics
+    fr_dict['error'] = fit_results.error
+    fr_dict['r_squared'] = fit_results.r_squared
+
+    return fr_dict
+
+@check_dependency(pd, 'pandas')
+def model_to_dataframe(fit_results, peak_org):
+    """Convert model fit results to a dataframe.
+
+    Parameters
+    ----------
+    fit_results : FOOOFResults
+        Results of a model fit.
+    peak_org : int or Bands
+        How to organize peaks.
+        If int, extracts the first n peaks.
+        If Bands, extracts peaks based on band definitions.
+
+    Returns
+    -------
+    pd.Series
+        Model results organized into a dataframe.
+    """
+
+    return pd.Series(model_to_dict(fit_results, peak_org))
+
+
+@check_dependency(pd, 'pandas')
+def group_to_dataframe(fit_results, peak_org):
+    """Convert a group of model fit results into a dataframe.
+
+    Parameters
+    ----------
+    fit_results : list of FOOOFResults
+        List of FOOOFResults objects.
+    peak_org : int or Bands
+        How to organize peaks.
+        If int, extracts the first n peaks.
+        If Bands, extracts peaks based on band definitions.
+
+    Returns
+    -------
+    pd.DataFrame
+        Model results organized into a dataframe.
+    """
+
+    return pd.DataFrame([model_to_dataframe(f_res, peak_org) for f_res in fit_results])