Maxwell Filter and Projections

Hi all,

Hopefully a quick question: Is there any reason why projections should not
be applied in advance of using the Maxwell filter? I understand that doing
so results in error for the Neuromag MaxFilter program, as well as the
newly implemented maxwell_filter function in mne-python. Any reason for
this?

Thanks!
-Sam

Hi Sam,

Maybe the place to start is why would you want to do this?

d

Absolutely. The short of it is that we were running into some issues with
the Maxwell filter function, such that the SVD will not converge for
certain buffers if there is an abundance of noise. As such, we were
performing epoch rejections on the raw to determine time periods with high
amounts of noise; removing these time periods from the raw; and passing the
resulting raw object to Maxwell filter.

The process above would yield the following processing steps:

   1. Perform epoch rejections
   2. Crop raw data and re-stitch good
   3. Maxwell filter data
   4. Save filtered raw data
   5. Make EOG/ECG projections
   6. Perform epoch rejections
   7. Save epochs

We were wondering if we could skip a few steps by making the EOG/ECG
projections ahead of time. That way we wouldn't need to re-epoch and save
out an additional raw file. Instead we could make and apply the projections
before Step (1). Does that make sense?

Thanks!
-Sam

The short of it is that we were running into some issues with the Maxwell
filter function, such that the SVD will not converge for certain buffers if
there is an abundance of noise.

I have observed seen this "failure to converge" problem with SVD before.
Are you running Anaconda? If so, could you look here, see if you can
reproduce the error with the .npy file from the post, and also do the `cat
/proc/cpuinfo` step? It will help isolate the problem, and keep attention
on it:

If you're not on Anaconda, I'm curious what Python setup you're using.

Such "high noise" segments really shouldn't be a problem for SVD. If the
Anaconda folks don't fix the issue, we're going to try to put a workaround
solution in upstream scipy:

Eric
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.nmr.mgh.harvard.edu/pipermail/mne_analysis/attachments/20160325/4154cdfb/attachment.html

Absolutely. The short of it is that we were running into some issues with
the Maxwell filter function, such that the SVD will not converge for certain
buffers if there is an abundance of noise. As such, we were performing epoch
rejections on the raw to determine time periods with high amounts of noise;
removing these time periods from the raw; and passing the resulting raw
object to Maxwell filter.

Sam, the easier way to do this, would be to find windows of noise, and feed
those to the call to maxfilter (using -skip, or the mne-python equivalent).

Honestly, there might be a much bigger issue lurking somewhere... Your
data shouldn't have noise that causes maxfilter not to work (we have an
incredible shielded room, and I don't think you should ever have noise issues
that maxfilter can't handle.

hth
d

@ Dan: We were trying to move away from the command line version for a
number of reasons. As far as I know, there is no -skip option in the python
implementation.

@ Eric: We are using the Anaconda distribution, and I was not able to
reproduce the error from the thread you linked. When we are getting an
error, it looks like this:

RuntimeWarning: invalid value encountered in greater

  num = np.sum(s > tol, dtype=int)

Traceback (most recent call last):

  File "motion_correction.py", line 32, in <module>
    raw = maxwell_filter(raw, st_duration=st_duration,
destination=trans.mean(axis=0), bad_condition='warning')

File "<string>", line 2, in maxwell_filter

File "/homes/5/szoro/.local/lib/python2.7/site-packages/mne/utils.py",
line 551, in verbose
    return function(*args, **kwargs)

File
"/homes/5/szoro/.local/lib/python2.7/site-packages/mne/preprocessing/maxwell.py",
line 404, in maxwell_filter
    t_proj = _overlap_projector(orig_in_data, resid, st_correlation)

File
"/homes/5/szoro/.local/lib/python2.7/site-packages/mne/preprocessing/maxwell.py",
line 1237, in _overlap_projector
    overwrite_a=True, mode='economic', **check_disable)[0].T

File
"/homes/5/szoro/.local/lib/python2.7/site-packages/scipy/linalg/decomp_qr.py",
line 142, in qr
    overwrite_a=overwrite_a)

File
"/homes/5/szoro/.local/lib/python2.7/site-packages/scipy/linalg/decomp_qr.py",
line 20, in safecall
    ret = f(*args, **kwargs)

ValueError: failed to create intent(cache|hide)|optional array-- must have

defined dimensions but got (0,)

Unfortunately I can't quite parse what's going on here in the error
message.

That error is different from the SVD one, but I suspect it has the same
underlying cause (Anaconda using outdated build tools). To avoid derailing
the original topic too much, can you open an mne-python issue for this
specific error scipy.linalg error? We can tackle debugging steps there. If
you can anonymize the file, upload it somewhere, give the minimal steps to
reproduce it, and paste the output of `mne.sys_info()` (assuming you're on
`master`) that would help.

Eric

What does the data look like where you need to "fix it"? It still
isn't clear what problems the data are causing.

Alternatively, you could just set the raw to 0s there (the same
equivalent as skip in maxfilter).

@ Dan: Ahh, we were wondering how the -skip function worked. (We looked
through the manual but did not specify.)

@ Eric: Absolutely will do.

All of which brings me back to the original question: Does anyone know a
reason why one shouldn't apply projections before Maxwell filtering?

Note, I am not a physicist, but I think this is roughly helpful:
(if a real physicist reads this, please correct me)

short version it doesn't make sense, both ssp and sss use the signal
space to remove artifacts
but sss pushes converts the data into spherical harmonics to do this,
but if you reduce the rank
with an ssp vector i don't think going into spherical harmonics would work.
I think this is also why you have to have a whole head type system to
work in maxfilter in general.

hth
d

My couple of cents

SSP and SSS (maxfilter) are both spatial filter.

SSS is a single point multipolar expansion using spherical harmonics basis
functions, assumes data to be represented in a particular spatial
orientation.

SSP rotates the spatial orientation hence basic assumption of the SSS
basis function no longer remain valid.

You can always do SSP after SSS

Sheraz

In other words, in order for the SSS to work if SSP has been already applied, the same SSP projector should be applied to the SSS multipolar basis vectors as well.

- Matti

And this is why it is pragmatically easier to do it the other way around
(SSP after SSS). Even though the rank of the data is not full after SSS,
SSP does a PCA/SVD so there will be a bunch of eigenvectors with near 0
eigenvalues, which does not cause problems.

In fact, you can use the SVD to estimate the true rank of your data after
SSS, which can be helpful. For example if you want to do multivariate
classification of the MEG data in sensor space, you might be better off
doing it in the eigenspace on only the eigenvectors that do not have near
0 eigenvalues associated with them rather than the full sensor matrix and
the eigenvectors will be orthogonal, so you can use a relatively simpler
classifier that assumes independence, like a Na?ve Bayes classifier.

This was all very helpful. Thank you all for your responses.