ENH: Add show hide menu with freesurfer atlas

[alexrockhill]

Fixes #26

[alexrockhill]

@larsoner do you have a minute to point me in the right direction for how to fix these CI errors? looks like they were solved on mne-python recently

[larsoner]

I would suggest to update pinning like mne-tools/mne-python#12580

[alexrockhill]

Ok for notes about this PR: I went to add show/hide for the head and I felt like there were a lot of buttons in a lot of different places and I looked at the text and a lot of those buttons were show/hide (max intensity projection, local maxima, brain slices) so I consolidated them into a show/hide menu while adding the show/hide head functionality.

That will play really well with #10 which adds a lot of other elements that it will be nice to show/hide separately.

Lastly, I wrapped into this PR showing the anatomical segmentation atlas because that's also very helpful to see relative to the contacts and also needs show/hide functionality. It looks really good, super useful I think, especially if you want to make a figure of one contact's location or really understand the positioning of that contact.

[alexrockhill]

Also is nice for showing/hiding surfaces for the volume source estimate viewer.

[alexrockhill]

Ok, just had to fix the style everything else passed. If you want to give it a quick spin, here's some scripts that aren't too heavy

import numpy as np
import matplotlib.pyplot as plt

import nibabel as nib
import nilearn.plotting
from dipy.align import resample

import mne
import mne_gui_addons as mne_gui
from mne.datasets import fetch_fsaverage
import mne_bids

# paths to mne datasets: sample sEEG and FreeSurfer's fsaverage subject,
# which is in MNI space
misc_path = mne.datasets.misc.data_path()
sample_path = mne.datasets.sample.data_path()
subjects_dir = sample_path / "subjects"

# use mne-python's fsaverage data
fetch_fsaverage(subjects_dir=subjects_dir, verbose=True)  # downloads if needed

# GUI requires pyvista backend
mne.viz.set_3d_backend("pyvistaqt")

T1 = nib.load(misc_path / "seeg" / "sample_seeg" / "mri" / "T1.mgz")
CT_orig = nib.load(misc_path / "seeg" / "sample_seeg_CT.mgz")

reg_affine = np.array(
    [
        [0.99270756, -0.03243313, 0.11610254, -133.094156],
        [0.04374389, 0.99439665, -0.09623816, -97.58320673],
        [-0.11233068, 0.10061512, 0.98856381, -84.45551601],
        [0.0, 0.0, 0.0, 1.0],
    ]
)

# use a cval='1%' here to make the values outside the domain of the CT
# the same as the background level during interpolation
CT_aligned = mne.transforms.apply_volume_registration(
    CT_orig, T1, reg_affine, cval="1%"
)
subj_trans = mne.coreg.estimate_head_mri_t("sample_seeg", misc_path / "seeg")

raw = mne.io.read_raw(misc_path / "seeg" / "sample_seeg_ieeg.fif")

# you may want to add `block=True` to halt execution until you have interacted
# with the GUI to find the channel positions, that way the raw object can
# be used later in the script (e.g. saved with channel positions)
self = mne_gui.locate_ieeg(
    raw.info,
    subj_trans,
    CT_aligned,
    subject="sample_seeg",
    subjects_dir=misc_path / "seeg",
)

import numpy as np

import mne_gui_addons

import mne
from mne.cov import compute_covariance
from mne.datasets import somato
from mne.time_frequency import csd_tfr
from mne.beamformer import make_dics, apply_dics_csd, make_lcmv, apply_lcmv_cov
from mne.minimum_norm import make_inverse_operator, apply_inverse_cov

print(__doc__)

# %%
# Reading the raw data and creating epochs:

data_path = somato.data_path()
subject = "01"
subjects_dir = data_path / "derivatives" / "freesurfer" / "subjects"
task = "somato"
raw_fname = (
    data_path
    / "sub-{}".format(subject)
    / "meg"
    / "sub-{}_task-{}_meg.fif".format(subject, task)
)

# crop to 5 minutes to save memory
raw = mne.io.read_raw_fif(raw_fname).crop(0, 300)

# We are interested in the beta band (12-30 Hz)
raw.load_data().filter(12, 30)

# The DICS beamformer currently only supports a single sensor type.
# We'll use the gradiometers in this example.
picks = mne.pick_types(raw.info, meg="grad", exclude="bads")

# Read epochs
events = mne.find_events(raw)
epochs = mne.Epochs(
    raw, events, event_id=1, tmin=-1.5, tmax=2, picks=picks, preload=True, decim=3
)

# Read forward operator and point to freesurfer subject directory
fwd_fname = (
    data_path
    / "derivatives"
    / "sub-{}".format(subject)
    / "sub-{}_task-{}-fwd.fif".format(subject, task)
)
fwd = mne.read_forward_solution(fwd_fname)

# %%
# Compute covariances and cross-spectral density
# ----------------------------------------------
# ERS activity starts at 0.5 seconds after stimulus onset. Because these
# data have been processed by MaxFilter directly (rather than MNE-Python's
# version), we have to be careful to compute the rank with a more conservative
# threshold in order to get the correct data rank (64). Once this is used in
# combination with an advanced covariance estimator like "shrunk", the rank
# will be correctly preserved.

rank = mne.compute_rank(epochs, tol=1e-6, tol_kind="relative")
win_active = (0.5, 1.5)
win_baseline = (-1, 0)
cov_baseline = compute_covariance(
    epochs,
    tmin=win_baseline[0],
    tmax=win_baseline[1],
    method="shrunk",
    rank=rank,
    verbose=True,
)
cov_active = compute_covariance(
    epochs,
    tmin=win_active[0],
    tmax=win_active[1],
    method="shrunk",
    rank=rank,
    verbose=True,
)

# when the covariance objects are added together, they are scaled by the size
# of the window used to create them so that the average is properly weighted
cov_common = cov_baseline + cov_active
cov_baseline.plot(epochs.info)

freqs = np.logspace(np.log10(12), np.log10(30), 9)

# time-frequency decomposition
epochs_tfr = mne.time_frequency.tfr_morlet(
    epochs,
    freqs=freqs,
    n_cycles=freqs / 2,
    return_itc=False,
    average=False,
    output="complex",
)
epochs_tfr.decimate(20)  # decimate for speed

surface = subjects_dir / subject / "bem" / "inner_skull.surf"
vol_src = mne.setup_volume_source_space(
    subject=subject,
    subjects_dir=subjects_dir,
    surface=surface,
    pos=10,
    add_interpolator=False,
)  # just for speed!

conductivity = (0.3,)  # one layer for MEG
model = mne.make_bem_model(
    subject=subject,
    ico=3,  # just for speed
    conductivity=conductivity,
    subjects_dir=subjects_dir,
)
bem = mne.make_bem_solution(model)

trans = fwd["info"]["mri_head_t"]
vol_fwd = mne.make_forward_solution(
    raw.info,
    trans=trans,
    src=vol_src,
    bem=bem,
    meg=True,
    eeg=True,
    mindist=5.0,
    n_jobs=1,
    verbose=True,
)

# Compute source estimate using MNE solver
snr = 3.0
lambda2 = 1.0 / snr**2
method = "MNE"  # use MNE method (could also be dSPM or sLORETA)

# make a different inverse operator for each frequency so as to properly
# whiten the sensor data
inverse_operator = list()
for freq_idx in range(epochs_tfr.freqs.size):
    # for each frequency, compute a separate covariance matrix
    cov_baseline = csd_baseline.get_data(index=freq_idx, as_cov=True)
    cov_baseline["data"] = cov_baseline["data"].real  # only normalize by real
    # then use that covariance matrix as normalization for the inverse
    # operator
    inverse_operator.append(
        mne.minimum_norm.make_inverse_operator(epochs.info, vol_fwd, cov_baseline)
    )

# finally, compute the stcs for each epoch and frequency
stcs = mne.minimum_norm.apply_inverse_tfr_epochs(
    epochs_tfr, inverse_operator, lambda2, method=method, pick_ori="vector"
)

# %%
# Plot volume source estimates
# ----------------------------

viewer = mne_gui_addons.view_vol_stc(
    stcs, subject=subject, subjects_dir=subjects_dir, src=vol_src, inst=epochs_tfr
)
viewer.go_to_extreme()  # show the maximum intensity source vertex
viewer.set_cmap(vmin=0.25, vmid=0.8)
viewer.set_3d_view(azimuth=40, elevation=35, distance=350)

[alexrockhill]

Ok to merge?