Skip to content

Commit

Permalink
MAINT: Update roadmap (mne-tools#12202)
Browse files Browse the repository at this point in the history
Co-authored-by: Daniel McCloy <[email protected]>
  • Loading branch information
larsoner and drammock authored Nov 13, 2023
1 parent 26a0cdc commit fa0e8cf
Show file tree
Hide file tree
Showing 9 changed files with 116 additions and 113 deletions.
3 changes: 0 additions & 3 deletions .pre-commit-config.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,3 @@ repos:
hooks:
- id: rstcheck
files: ^doc/.*\.(rst|inc)$
# https://github.com/rstcheck/rstcheck/issues/199
# https://github.com/rstcheck/rstcheck/issues/200
exclude: ^doc/(help/faq|install/manual_install|install/mne_c|install/advanced|install/updating|_includes/channel_interpolation|_includes/inverse|_includes/ssp)\.rst$
11 changes: 5 additions & 6 deletions doc/_includes/channel_interpolation.rst
Original file line number Diff line number Diff line change
Expand Up @@ -22,12 +22,11 @@ In short, data repair using spherical spline interpolation :footcite:`PerrinEtAl
Spherical splines assume that the potential :math:`V(\boldsymbol{r_i})` at any point :math:`\boldsymbol{r_i}` on the surface of the sphere can be represented by:

.. math:: V(\boldsymbol{r_i}) = c_0 + \sum_{j=1}^{N}c_{i}g_{m}(cos(\boldsymbol{r_i}, \boldsymbol{r_{j}}))
:label: model
:name: model

where the :math:`C = (c_{1}, ..., c_{N})^{T}` are constants which must be estimated. The function :math:`g_{m}(\cdot)` of order :math:`m` is given by:

.. math:: g_{m}(x) = \frac{1}{4 \pi}\sum_{n=1}^{\infty} \frac{2n + 1}{(n(n + 1))^m}P_{n}(x)
:label: legendre

where :math:`P_{n}(x)` are `Legendre polynomials`_ of order :math:`n`.

Expand All @@ -36,22 +35,22 @@ where :math:`P_{n}(x)` are `Legendre polynomials`_ of order :math:`n`.
To estimate the constants :math:`C`, we must solve the following two equations simultaneously:

.. math:: G_{ss}C + T_{s}c_0 = X
:label: matrix_form
:name: matrix_form

.. math:: {T_s}^{T}C = 0
:label: constraint
:name: constraint

where :math:`G_{ss} \in R^{N \times N}` is a matrix whose entries are :math:`G_{ss}[i, j] = g_{m}(cos(\boldsymbol{r_i}, \boldsymbol{r_j}))` and :math:`X \in R^{N \times 1}` are the potentials :math:`V(\boldsymbol{r_i})` measured at the good channels. :math:`T_{s} = (1, 1, ..., 1)^\top` is a column vector of dimension :math:`N`. Equation :eq:`matrix_form` is the matrix formulation of Equation :eq:`model` and equation :eq:`constraint` is like applying an average reference to the data. From equation :eq:`matrix_form` and :eq:`constraint`, we get:

.. math:: \begin{bmatrix} c_0 \\ C \end{bmatrix} = {\begin{bmatrix} {T_s}^{T} && 0 \\ T_s && G_{ss} \end{bmatrix}}^{-1} \begin{bmatrix} 0 \\ X \end{bmatrix} = C_{i}X
:label: estimate_constant
:name: estimate_constant

:math:`C_{i}` is the same as matrix :math:`{\begin{bmatrix} {T_s}^{T} && 0 \\ T_s && G_{ss} \end{bmatrix}}^{-1}` but with its first column deleted, therefore giving a matrix of dimension :math:`(N + 1) \times N`.

Now, to estimate the potentials :math:`\hat{X} \in R^{M \times 1}` at the bad channels, we have to do:

.. math:: \hat{X} = G_{ds}C + T_{d}c_0
:label: estimate_data
:name: estimate_data

where :math:`G_{ds} \in R^{M \times N}` computes :math:`g_{m}(\boldsymbol{r_i}, \boldsymbol{r_j})` between the bad and good channels. :math:`T_{d} = (1, 1, ..., 1)^\top` is a column vector of dimension :math:`M`. Plugging in equation :eq:`estimate_constant` in :eq:`estimate_data`, we get

Expand Down
24 changes: 12 additions & 12 deletions doc/_includes/inverse.rst
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ this by writing :math:`R' = R/ \lambda^2 = R \lambda^{-2}`, which yields the
inverse operator

.. math::
:label: inv_m
:name: inv_m
M &= R' G^\top (G R' G^\top + C)^{-1} \\
&= R \lambda^{-2} G^\top (G R \lambda^{-2} G^\top + C)^{-1} \\
Expand Down Expand Up @@ -106,12 +106,12 @@ The MNE software employs data whitening so that a 'whitened' inverse operator
assumes the form

.. math:: \tilde{M} = M C^{^1/_2} = R \tilde{G}^\top (\tilde{G} R \tilde{G}^\top + \lambda^2 I)^{-1}\ ,
:label: inv_m_tilde
:name: inv_m_tilde

where

.. math:: \tilde{G} = C^{-^1/_2}G
:label: inv_g_tilde
:name: inv_g_tilde

is the spatially whitened gain matrix. We arrive at the whitened inverse
operator equation :eq:`inv_m_tilde` by making the substitution for
Expand All @@ -128,7 +128,7 @@ operator equation :eq:`inv_m_tilde` by making the substitution for
The expected current values are

.. math::
:label: inv_j_hat_t
:name: inv_j_hat_t
\hat{j}(t) &= Mx(t) \\
&= M C^{^1/_2} C^{-^1/_2} x(t) \\
Expand All @@ -137,7 +137,7 @@ The expected current values are
knowing :eq:`inv_m_tilde` and taking

.. math::
:label: inv_tilde_x_t
:name: inv_tilde_x_t
\tilde{x}(t) = C^{-^1/_2}x(t)
Expand All @@ -151,7 +151,7 @@ to raw data. To reflect the decrease of noise due to averaging, this matrix,
C_0 / L`.

.. note::
When EEG data are included, the gain matrix :math:`G` needs to be average referenced when computing the linear inverse operator :math:`M`. This is incorporated during creating the spatial whitening operator :math:`C^{-^1/_2}`, which includes any projectors on the data. EEG data average reference (using a projector) is mandatory for source modeling and is checked when calculating the inverse operator.
When EEG data are included, the gain matrix :math:`G` needs to be average referenced when computing the linear inverse operator :math:`M`. This is incorporated during creating the spatial whitening operator :math:`C^{-^1/_2}`, which includes any projectors on the data. EEG data average reference (using a projector) is mandatory for source modeling and is checked when calculating the inverse operator.

As shown above, regularization of the inverse solution is equivalent to a
change in the variance of the current amplitudes in the Bayesian *a priori*
Expand Down Expand Up @@ -224,7 +224,7 @@ computational convenience we prefer to take another route, which employs the
singular-value decomposition (SVD) of the matrix

.. math::
:label: inv_a
:name: inv_a
A &= \tilde{G} R^{^1/_2} \\
&= U \Lambda V^\top
Expand All @@ -238,7 +238,7 @@ Combining the SVD from :eq:`inv_a` with the inverse equation :eq:`inv_m` it is
easy to show that

.. math::
:label: inv_m_tilde_svd
:name: inv_m_tilde_svd
\tilde{M} &= R \tilde{G}^\top (\tilde{G} R \tilde{G}^\top + \lambda^2 I)^{-1} \\
&= R^{^1/_2} A^\top (A A^\top + \lambda^2 I)^{-1} \\
Expand All @@ -253,23 +253,23 @@ where the elements of the diagonal matrix :math:`\Gamma` are simply
.. `reginv` in our code:
.. math::
:label: inv_gamma_k
:name: inv_gamma_k
\gamma_k = \frac{\lambda_k}{\lambda_k^2 + \lambda^2}\ .
From our expected current equation :eq:`inv_j_hat_t` and our whitened
measurement equation :eq:`inv_tilde_x_t`, if we take

.. math::
:label: inv_w_t
:name: inv_w_t
w(t) &= U^\top \tilde{x}(t) \\
&= U^\top C^{-^1/_2} x(t)\ ,
we can see that the expression for the expected current is just

.. math::
:label: inv_j_hat_t_svd
:name: inv_j_hat_t_svd
\hat{j}(t) &= R^{^1/_2} V \Gamma w(t) \\
&= \sum_k {\bar{v_k} \gamma_k w_k(t)}\ ,
Expand Down Expand Up @@ -314,7 +314,7 @@ normalization factors, it's convenient to reuse our "weighted eigenleads"
definition from equation :eq:`inv_j_hat_t` in matrix form as

.. math::
:label: inv_eigenleads_weighted
:name: inv_eigenleads_weighted
\bar{V} = R^{^1/_2} V\ .
Expand Down
12 changes: 6 additions & 6 deletions doc/_includes/ssp.rst
Original file line number Diff line number Diff line change
Expand Up @@ -30,13 +30,13 @@ Without loss of generality we can always decompose any :math:`n`-channel
measurement :math:`b(t)` into its signal and noise components as

.. math:: b(t) = b_s(t) + b_n(t)
:label: additive_model
:name: additive_model

Further, if we know that :math:`b_n(t)` is well characterized by a few field
patterns :math:`b_1 \dotso b_m`, we can express the disturbance as

.. math:: b_n(t) = Uc_n(t) + e(t)\ ,
:label: pca
:name: pca

where the columns of :math:`U` constitute an orthonormal basis for :math:`b_1
\dotso b_m`, :math:`c_n(t)` is an :math:`m`-component column vector, and the
Expand All @@ -48,12 +48,12 @@ such that the conditions described above are satisfied. We can now construct
the orthogonal complement operator

.. math:: P_{\perp} = I - UU^\top
:label: projector
:name: projector

and apply it to :math:`b(t)` in Equation :eq:`additive_model` yielding

.. math:: b_{s}(t) \approx P_{\perp}b(t)\ ,
:label: result
:name: result

since :math:`P_{\perp}b_n(t) = P_{\perp}(Uc_n(t) + e(t)) \approx 0` and
:math:`P_{\perp}b_{s}(t) \approx b_{s}(t)`. The projection operator
Expand Down Expand Up @@ -102,12 +102,12 @@ typical in EEG analysis to subtract the average reference from all the sensor
signals :math:`b^{1}(t), ..., b^{n}(t)`. That is:

.. math:: {b}^{j}_{s}(t) = b^{j}(t) - \frac{1}{n}\sum_{k}{b^k(t)}
:label: eeg_proj
:name: eeg_proj

where the noise term :math:`b_{n}^{j}(t)` is given by

.. math:: b_{n}^{j}(t) = \frac{1}{n}\sum_{k}{b^k(t)}
:label: noise_term
:name: noise_term

Thus, the projector vector :math:`P_{\perp}` will be given by
:math:`P_{\perp}=\frac{1}{n}[1, 1, ..., 1]`
Expand Down
2 changes: 2 additions & 0 deletions doc/conf.py
Original file line number Diff line number Diff line change
Expand Up @@ -1324,6 +1324,8 @@ def reset_warnings(gallery_conf, fname):
# nilearn
"pkg_resources is deprecated as an API",
r"The .* was deprecated in Matplotlib 3\.7",
# scipy
r"scipy.signal.morlet2 is deprecated in SciPy 1\.12",
):
warnings.filterwarnings( # deal with other modules having bad imports
"ignore", message=".*%s.*" % key, category=DeprecationWarning
Expand Down
Loading

0 comments on commit fa0e8cf

Please sign in to comment.