mne-python: Size of resolution matrix (make_inverse_resolution_matrix) inconsistent with documentation
Describe the bug
It is unclear what the correct way to compute the resolution matrix is when using mne.minimum_norm.make_inverse_resolution_matrix with an inverse operator that has loose orientations. The function _get_matrix_from_inverse_operator is called (in line 59) to return a matrix that can be multiplied with the leadfield, and the documentation there states:
For loose orientation constraint, the CTFs are computed for the normal component (pick_ori=‘normal’).
However, it does not seem like pick_ori='normal' is actually used.
Steps to reproduce
import mne
data_path = mne.datasets.sample.data_path()
fname_inv = data_path + '/MEG/sample/sample_audvis-eeg-oct-6-eeg-inv.fif' # free source orientation
inverse_operator = mne.minimum_norm.read_inverse_operator(fname_inv)
fname_fwd = data_path + '/MEG/sample/sample_audvis-eeg-oct-6-fwd.fif'
forward = mne.read_forward_solution(fname_fwd)
method = "sLORETA"
snr = 3.
lambda2 = 1. / snr ** 2
res_matrix = mne.minimum_norm.make_inverse_resolution_matrix(forward,inverse_operator,method=method,lambda2=lambda2)
Expected results
I expect the res_matrix to have dimensions nsources x nsources, where nsources is the number of sources in the source space.
Actual results
The res_matrix has dimensions (3nsources) x (3nsources).
Additional information
Output of mne.sys_info():
Qt WebEngine seems to be initialized from a plugin. Please set Qt::AA_ShareOpenGLContexts using QCoreApplication::setAttribute before constructing QGuiApplication.
Platform: Linux-5.11.0-44-generic-x86_64-with-debian-bullseye-sid
Python: 3.7.9 (default, Aug 31 2020, 12:42:55) [GCC 7.3.0]
Executable: /home/katharina/snap/anaconda3/envs/py37cmp-gui/bin/python
CPU: x86_64: 16 cores
Memory: 31.3 GB
mne: 0.20.7
numpy: 1.18.5 {blas=mkl_rt, lapack=mkl_rt}
scipy: 1.5.2
matplotlib: 3.2.2 {backend=Qt5Agg}
sklearn: 1.0.1
numba: Not found
nibabel: 3.2.1
cupy: Not found
pandas: 1.3.4
dipy: 1.1.1
mayavi: Not found
pyvista: Not found
vtk: Not found
About this issue
- Original URL
- State: closed
- Created 2 years ago
- Comments: 15 (7 by maintainers)
Hello!
Sorry for the long silence - I’m still on this! Thank you for the discussion, it really helped me to understand that a resolution matrix of (3ndipoles) x (3ndipoles) is the correct way to go when handling free orientations.
I tried the solution in pull request #8639, but I think it doesn’t entirely solve the problem. Yes, in
_get_psf_ctf, there is a bit that reduces the vector to a scalar (starting at line 150 inresolution_matrix.py), but it only does that in one direction and I still end up with a resolution matrix of ndipoles x (3*ndipoles). Unless I’m missing something?The solution I am using right now is to compute the full matrix and then take the spectral norm of each 3x3 matrix, as suggested by @larsoner and @pettetmw It’s quite slow and not very memory efficient, but I don’t understand that einsum business well enough to make a better suggestion 😬 and it works for now.
So, thanks for your help! I’ll keep an eye on this issue, maybe the
combine_xyz='spectral'argument will find its way into a new version 😁… sorry I meant try #8639