Skip to content

Commit

Permalink
small updates
Browse files Browse the repository at this point in the history
  • Loading branch information
robertoostenveld committed Apr 26, 2021
1 parent 18c0b25 commit d84c83e
Show file tree
Hide file tree
Showing 2 changed files with 19 additions and 13 deletions.
24 changes: 15 additions & 9 deletions exercises_code/ni2_inverse5_scanning.m
Original file line number Diff line number Diff line change
Expand Up @@ -97,6 +97,7 @@
cfg.headmodel = headmodel;
% cfg.normalize = 'yes';
sourcemodel = ft_prepare_leadfield(cfg);

L = cat(2,sourcemodel.leadfield{sourcemodel.inside});

% compute the covariance
Expand Down Expand Up @@ -236,41 +237,45 @@
cfg.funparameter = 'avg.pow';
cfg.method = 'slice';
cfg.nslices = 10;
cfg.funcolorlim=[0 0.2];
cfg.funcolorlim = [0 0.2];
ft_sourceplot(cfg,source);

%% contrast between 2 conditions

% create the sensor data for the second condition
sensordata2 = 1.25.*leadfield1*s1 + 0.8.*leadfield2*s2 + 0.8.*leadfield3*s3 + 1.25.*leadfield4*s4 + randn(301,1000)*0.04e-8;

% compute the covariance
C2 = cov([sensordata sensordata2]');
iC2 = pinv(C2); % so that we compute it only once
% compute the covariance and its inverse for both conditions combined
C2 = cov([sensordata sensordata2]');
iC2 = inv(C2);
iCr2 = inv(C2+eye(301)*1e-19);

% compute the beamformer spatial filter
% compute the beamformer spatial filter for condition 2
for ii = 1:size(L,2)/3
indx=(ii-1)*3+(1:3);
Lr = L(:,indx); % Lr is the leadfield for source r
wbfr2(indx,:)=pinv(Lr'*iCr2*Lr)*Lr'*iCr2;
end
sbfr2 = wbfr2*sensordata2;
pbfr2 = var(sbfr2,[],2);
pbfr2 = sum(reshape(pbfr2,3,[]));

sbfr1 = wbfr2*sensordata;
pbfr1 = var(sbfr1,[],2);
pbfr1 = sum(reshape(pbfr1,3,[]));

sbfr2 = wbfr2*sensordata2;
pbfr2 = var(sbfr2,[],2);
pbfr2 = sum(reshape(pbfr2,3,[]));

source.avg.pow(source.inside)=(pbfr1-pbfr2)./(pbfr1+pbfr2);

cfg = [];
cfg.funparameter='avg.pow';
cfg.method='slice';
cfg.method='ortho';
cfg.nslices = 10;
cfg.funcolorlim=[-.2 .2];
ft_sourceplot(cfg,source);

%% correlated sources

sens = ni2_sensors('type','meg');
headmodel = ni2_headmodel('type','spherical','nshell',1);

Expand All @@ -295,6 +300,7 @@
cfg.vol = headmodel;
% cfg.normalize = 'yes';
sourcemodel = ft_prepare_leadfield(cfg);

L = cat(2,sourcemodel.leadfield{sourcemodel.inside});

% compute the covariance
Expand Down
8 changes: 4 additions & 4 deletions exercises_text/ni2_inverse5_scanning.md
Original file line number Diff line number Diff line change
Expand Up @@ -304,10 +304,10 @@ We will use the variables `sensordata`, `sourcemodel` and `L` that we also used

Now we will simulate data from a 'second’ condition (compared to the original variable sensordata), where the sources have the exact same locations and time courses, but the amplitude of two sources is decreased, and the amplitude of the other sources is increased, relative to the 'first’ condition.

sensordata2 = 1.25.*leadfield1_s1 + ...
0.80.*leadfield2_s2 + ...
0.80.*leadfield3_s3 + ...
1.25.*leadfield4_s4 + ...
sensordata2 = 1.25 .* leadfield1*s1 + ...
0.80 .* leadfield2*s2 + ...
0.80 .* leadfield3*s3 + ...
1.25 .* leadfield4*s4 + ...
randn(301, 1000)*0.04e-8;

We will now compute the spatial filters using the covariance estimated from the data combined across the two conditions. In this way we will end up with a single set of spatial filters, and thus, in comparing across the two conditions, we ensure that the depth bias is the same for condition 1 and 2. The covariance of the data combined across the 2 conditions can be computed in several ways, e.g. by averaging the single condition covariances. Here, we first concatenate the data and subsequently compute the covariance. Now we will use the MATLAB `cov` function, rather than computing the covariance 'by hand’ (i.e. by first computing the mean across time etc.).
Expand Down

0 comments on commit d84c83e

Please sign in to comment.