Back to the main page.

Bug 3298 - Correspondence channels-coils in HCP files

Status CLOSED FIXED
Reported 2017-05-23 10:25:00 +0200
Modified 2018-01-26 09:50:21 +0100
Product: FieldTrip
Component: documentation
Version: unspecified
Hardware: All
Operating System: All
Importance: P5 trivial
Assigned to: Jan-Mathijs Schoffelen
URL:
Tags:
Depends on:
Blocks:
See also:

François Tadel - 2017-05-23 10:25:09 +0200

Hello, I tried on the mailing list, but didn't get any reply: https://mailman.science.ru.nl/pipermail/fieldtrip/2017-May/011546.html I wrote a FieldTrip>Brainstorm import function based on Elekta files, and I was using the .tra matrix to get the correspondence between the channels and the coils: https://github.com/brainstorm-tools/brainstorm3/blob/master/toolbox/io/in_data_fieldtrip.m#L107 I based this code on your description of the .tra field: "NxM matrix with the weight of each coil into each channel": http://www.fieldtriptoolbox.org/faq/how_are_electrodes_magnetometers_or_gradiometers_described This doesn't work with the processed 4D recordings from the HCP dataset because the .tra matrix is dense. Is this because the 4D compensation matrix is integrated in this .tra matrix? Or because there is some additional PCA/ICA projection done during the preprocessing? How do I recover the sparse .tra information that simply contains the weights channels-coils from the HCP processed files? Thanks, François


Robert Oostenveld - 2017-05-23 11:40:43 +0200

how about ft_read_sens (or ft_read_header) on the original 4D files?


François Tadel - 2017-05-23 13:38:47 +0200

Yes, I guess we could read directly the positions of the coils and noise compensation matrix from the original file. But this is a completely different processing pipeline, and we lose the possible ICA/SSP projectors that were applied (which we would like apply on the forward model as well). What I was primarily trying to write here is a conversion function from FieldTrip to Brainstorm. Is there any way to recover the sparse channels-coils weight matrix from a saved FieldTrip structure (such as the ones distributed in the HCP)?


Jan-Mathijs Schoffelen - 2017-06-17 18:05:22 +0200

you could try ft_apply_montage while specifying 'inverse','yes' as a key-value pair, with the appropriate sequence of transformations applied to the tra-matrix (in reverse order).


François Tadel - 2017-08-19 15:22:19 +0200

Hi Jan-Mathijs, Sorry I'm getting back to this topic only now. I'm not sure I understand what you mean. Is there a history of the modifications of the .tra matrix saved in the .mat files? Could you give a precise sequence of Matlab/FieldTrip commands I could use to get the sparse coil-channel correspondence matrix from the pre-processed MEG recordings in the HCP database? If I could get the two other matrices (4D noise compensation matrix, SSP/ICA projectors) separately, this would be even better. Thanks Francois


François Tadel - 2017-08-24 18:23:24 +0200

Looking again at the pre-processed HCP files, I got what you meant. The balance field contains the history of all the transformations, and could think about inverting them one after the other. I was hoping to get something like this to work: % Get the list of montages to undo montageList = {grad.balance.current}; if ~isempty(grad.balance.previous) montageList = cat(2, montageList, grad.balance.previous{:}); end % Undo montages one by one for iMontage = 1:length(montageList) % Make sure the montage is defined in the structure mont = montageList{iMontage}; if ~isfield(grad.balance, mont) continue; end % Reverse transformation grad = ft_apply_montage(grad, grad.balance.(mont), 'inverse', 'yes'); end % Remove small values (to keep only the ones) grad.tra(abs(grad.tra) < 0.1) = 0; for i = 1:nChannels % Find channel index ichan = find(strcmpi(ftMat.label{i}, grad.label)); % Find corresponding coils icoils = find(grad.tra(ichan,:)); [...] end I have many problems though: 1) The output channel units in "pca" are set to "unknown". Error in ft_scalingfactor (line 181) ft_error('cannot convert %s to %s', old, new); (And it looks like some montages include the field "chanunitorg", and others the field "chanunitold"). So in my loop I added the following code: % Fix units and types if isfield(grad.balance.(mont), 'chanunitnew') && strcmpi(grad.balance.(mont).chanunitnew{1}, 'unknown') if isfield(grad.balance.(mont), 'chanunitorg') grad.balance.(mont).chanunitnew = grad.balance.(mont).chanunitorg; elseif isfield(grad.balance.(mont), 'chanunitold') grad.balance.(mont).chanunitnew = grad.balance.(mont).chanunitold; end end if isfield(grad.balance.(mont), 'chantypenew') && strcmpi(grad.balance.(mont).chantypenew{1}, 'unknown') if isfield(grad.balance.(mont), 'chantypeorg') grad.balance.(mont).chantypenew = grad.balance.(mont).chantypeorg; elseif isfield(grad.balance.(mont), 'chantypeold') grad.balance.(mont).chantypenew = grad.balance.(mont).chantypeold; end end 2) The channel types across the various transformations do not match. In some the sensors are called "meg", in others "megmag". Error in ft_apply_montage (line 333) ft_error('inconsistent chantype in data and montage'); How am I supposed to fix these inconsistencies? To reproduce these errors: my current test dataset is from the HCP-MEG release, subject 891667, MEG/Restin/rmegpreproc/891667_MEG_3-Restin_rmegpreproc.mat Could you please help me with this?


Jan-Mathijs Schoffelen - 2017-08-24 21:59:26 +0200

Typically, such glitches might be dealt with properly by calling: grad = ft_datatype_sens(grad), prior to start the inverse balancing procedure. This should at least deal with the old/org inconsistentcies at the highest level of the grad, but I am not sure whether this will solve it recursively in the balance structures. It's worth a shot though. As a side note, I'd use ft_apply_montage with the key-value pair 'keepunused', 'yes'. Although I don't think that this is causing you problems.


François Tadel - 2017-08-27 11:14:25 +0200

Thanks, but these two suggestions do not help. The function ft_datatype_sens do not fix the structures in the .balance field. If this is not too complicated to do, could you add it to this function? Given the inconsistencies in the structures of the HCP pre-processed files, it sounds like an interesting feature. Or do you have any alternative to suggest? (I'm trying to address requests from people who contacted me because they are trying to re-use these files for various purposes)


Jan-Mathijs Schoffelen - 2017-08-27 20:40:09 +0200

Hi Francois, Unfortunately, I don't have time to solve this for you right now. If you want a speedy resolution, I suggest that you make a proposal for some code improvement through github: github.com/fieldtrip/fieldtrip.git The following link: http://www.fieldtriptoolbox.org/development/git explains a relatively painfree way to contribute code. If you want me to look and think along, it would help if you upload an example of such a faulty grad structure. This would save me a couple of hours to lay my hands on some HCP data myself, of which, surprisingly, I don't have a copy lying around.


François Tadel - 2017-08-28 16:03:44 +0200

I posted an example as a github issue (easier to handle than bugzilla): https://github.com/fieldtrip/fieldtrip/issues/507 You just have to download the .zip, run the .m file, and tell me what is wrong with the structure of the HCP files, so I can fix it in my code. Thanks


Jan-Mathijs Schoffelen - 2017-11-28 12:15:54 +0100

I believe that this has been addressed