OPM array simulation and source reconstruction

[qishengjie9]

hello, I find the OPM sensor number is different from yours in your paper titled Pragmatic spatial sampling for wearable MEG arrays, such as 940 for 10mm space vs 784 in your paper. The codes I use are as follows. In addition, I got a problem when I try to use the sensor number as the spatial number in the module of source inversion iterative with the error information of not many this lead field.

for i = 1:Narray
S =[];
S.space = space(i);
S.meshres = 3;
S.offset = 6.5;
S.wholehead = 0;
S.fname = sprintf('sim_opm_%dmm',S.space);
D = spm_opm_sim(S);
newfile=fullfile(out_path, sprintf('sim_opm_%dmm.mat',S.space));
save(newfile,'D');
end

[tierneytim]

Hi, the algorithm for sensor placement is not deterministic and has changed since the publication so the number of sensors will change when you run it. As for the iterative source inversion, I did not write that code and am not familiar with the issue you have described. It might be best to contact the author directly

[qishengjie9]

Thanks for your kind help and advice. I will continue to try to solve the problem!

[georgeoneill]

Hello @qishengjie9 if you could copy in the code you use to perform the source inversion, and copy the error messages you get, we should be able to identify what the problem is.

[qishengjie9]

Hello @georgeoneill I have copied the code and error in my email. I guess the problem is that the rank of the lead field matrix decomposed by svd is less than the sensor number when the number is 940 in my situation. You can check the code to give me inspiration and suggestions. Thanks for your kind help.

[qishengjie9]

Hello，these are the snippet of my code where the error occurs and the error message. I find the error occured after performing the svd of the lead field matrix.

[georgeoneill]

I think this hsould be easy to fix.

You want to avoid the spatial decomposotion step entirely,

To do this, if n_channels is the number of channels of the dataset, set the code do be;

matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.nsmodes = n_channels;
matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.umodes= eye(n_channels);

the umodes option may need to be specified as a cell e.g. umodes= {eye(n_channels)};

[qishengjie9]

The umodes I have uesd is the same as yours. And I found I can not avoid the spatial decomposition. The error message is the same as before. The spatial number available is 596 which is less than the sensor number of 940 in my simulation.

[georgeoneill]

Okay, can I think we need to check that everything communicated from the batch script is getting to the function.

1) can you show me the entire section of the matlabbatch you are using to call the function?
2) can you put a breakpoint at the line of the error. Run the code and it should pause at just before it throws an error. When it does, can you type inverse and report back what the command line returns?

[qishengjie9]

Hello，I have done as you suggested. The matlabbatch code are as follows.

matlabbatch = []; matlabbatch{1}.spm.meeg.source.invertiter.D = {simfile}; matlabbatch{1}.spm.meeg.source.invertiter.val = 1; matlabbatch{1}.spm.meeg.source.invertiter.whatconditions.all = 1; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.invfunc = 'Classic'; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.invtype = algorithm; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.woi = [0 1000]; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.foi = [0 40]; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.hanning = 1; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.isfixedpatch.randpatch.npatches = 512; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.isfixedpatch.randpatch.niter = 1; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.patchfwhm = 0.6; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.mselect = 0; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.nsmodes = Nmodes; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.umodes = spatialmodesname; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.ntmodes = []; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.priors.priorsmask = {''}; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.priors.space = 1; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.restrict.locs = zeros(0, 3); matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.restrict.radius = 32; matlabbatch{1}.spm.meeg.source.invertiter.isstandard.custom.outinv = ''; matlabbatch{1}.spm.meeg.source.invertiter.modality = {'All'}; matlabbatch{1}.spm.meeg.source.invertiter.crossval = [0 1]; [a2,~]=spm_jobman('run', matlabbatch);

And the report is shown below after adding a breakpoint and type inverse.
`Contains the following fields struct:

   priors: []
     type: 'EBB'
      woi: [0 999]
      Han: 1
      lpf: 0
      hpf: 40
mergeflag: 0
       Np: 512
       Nm: 940
       Nt: 4
   smooth: 0.6000
 modality: 'MEG'
     allF: 0
  PostMax: [20484×1 double]
        A: {}
       Ip: []

`

[georgeoneill]

It looks like the spatial modes are not being transferred to the function (it would be in inverse.A), the question I have is what is spatialmodesname? Is it a string? I beleive it should be either an array or an array ecapsulated in a cell. (e.g. umodes= {eye(n_channels)};) can you check please?

[qishengjie9]

Thanks for your suggestions. And the spatialmodesname is a string pointing to the spatial modes which contains the {eye(n_channels)}. The inverse.A is also empty after specifying umodes= {eye(n_channels)}. I think the inverse.A is dispensable if I choose all the sensor data(ie the number of spatial modes is equal to the number of sensors). Additionally, the inverse process will run successfully when the number of sensors is not that many in my simulation(eg about 200 sensors). I don't know if I understand correctly. Could you please give me some advice.