Connectivity analysisBrain functioning is based on the dynamic interactions between functional specialized areas. Many connectivity measures have been proposed to quantify the interactions occurring in the brain, and each metric has its own advantages and disadvantages. An important caveat to consider is that connectivity estimates between sensor recordings – namely between signals picked up by the electrodes – are difficult to interpret because of the volume conduction problem. Therefore, we strongly advise to study neuronal interactions on the level of the reconstructed sources. Another important caveat is that, even when working at the source level, different problems other than the volume conduction may affect the interpretation of the estimated brain interactions (see Bastos and Schoffelen, 2016). In the second part of the previous section (Source analysis), we have computed the source activity from signals recorded during resting state. We will recompute the variable ‘srcELoreta’ as a starting point to estimate source connectivity. 
Return to the wiki

%% Loading relevant variables to recompute srcELoreta load('workspaceSrcELoreta.mat') % Recomputing source reconstruction using eLORETA cfg = []; cfg.method = 'eloreta'; % eLORETA method cfg.eloreta.lambda = 0.05; % parameter for the regularization cfg.eloreta.keepfilter = 'yes'; % to keep the spatial filters in the output cfg.headmodel = headModel; cfg.sourcemodel = leadField; srcELoreta = ft_sourceanalysis(cfg,freq); srcELoreta.cfg = []; % to save memory srcELoreta.avg.momdimord = '{pos}_ori_rpttap_freq'; srcELoreta.avg.oridimord = '{pos}_ori_freq'; clear cfg 
↑ Go up 
Notice that the field ‘srcELoreta.avg.mom’ contains three channels at each of the 8004 positions on the cortical surface. The three channels represent the x, the y and the zcomponent of the Fourier transform of the equivalent current dipole source. Computing connectivity between triplets of channels makes the interpretation of results much more difficult. Therefore, we will project the frequency representations of the source activity along the dipole orientation that explains most variance. This projection can be performed using the singular value decomposition that allow us determining the largest eigenvector. 

%% Projecting along the strongest dipole direction srcELoretaProj = srcELoreta; inside = srcELoreta.inside; tmpMom = cell(1,length(srcELoreta.avg.mom)); for i = 1:length(srcELoreta.avg.mom) if inside(i) for j = 1:length(srcELoreta.freq) % Using singular value decomposition [u,s,~] = svd(real(srcELoreta.avg.mom{i}(:,:,j)),'econ'); tmpMom{i}(1,:,j) = u(:,1).' * srcELoreta.avg.filter{i} * ... freq.fourierspctrm(:,:,j).'; clear u s end end end srcELoretaProj.avg.mom = tmpMom; clear tmpMom i j inside 
↑ Go up 
We will use the FieldTrip function ‘ft_connectivityanalysis’ for the actual computation of connectivity. In particular, we will estimate the imaginary part of the coherency (Nolte et al., 2004), by specifying cfg.method = ‘coh’ and cfg.complex = ‘absimag’. 

%% Connectivity analysis % Computing connectivity for each frequency using the imaginary part of coherency tmpConn = cell(1,length(srcELoretaProj.freq)); inside = leadField.inside; cfg = []; cfg.method = 'coh'; cfg.complex = 'absimag'; for i = 1:length(srcELoretaProj.freq) tmpMom = srcELoretaProj.avg.mom; for j = 1:length(tmpMom) if inside(j) tmpMom{j} = tmpMom{j}(:,:,i); end end tmpSrc = srcELoretaProj; tmpSrc.avg.mom = tmpMom; tmpConn{i} = ft_connectivityanalysis(cfg,tmpSrc); end % Collapsing the cell array tmpConn in a single structure srcConn = tmpConn{1}; tmp = nan(size(srcConn.cohspctrm,1),size(srcConn.cohspctrm,2), ... length(srcELoretaProj.freq)); for i = 1:length(srcELoretaProj.freq) tmp(:,:,i) = tmpConn{i}.cohspctrm; end srcConn.cohspctrm = tmp; % Computing connectivity collapsing all frequencies tmpSrc = srcELoretaProj; for j = 1:length(tmpSrc.avg.mom) if inside(j) tmpSrc.avg.mom{j} = mean(tmpSrc.avg.mom{j},3); end end srcConnAvg = ft_connectivityanalysis(cfg,tmpSrc); srcConnAvg.freq = mean(srcConnAvg.freq); srcConnAvg.dimord = 'pos_pos'; clear i j cfg tmp tmpMom tmpSrc tmpConn inside 
↑ Go up 
The field ‘srcConn.cohspctrm’ contains the connectivity matrix estimated for each frequency. Each connectivity matrix contains about 64 million (= 8004 x 8004) elements and, because of the low spatial resolution of EEG, it is convenient to reduce the matrix dimensions using a parcellation approach (see also the section ‘Source analysis’). 

%% Loading atlas and applying a parcellation scheme to reduce data dimensionality load('atlas_MMP1.0_4k.mat') atlasMMP1 = atlas; atlasMMP1 = ft_convert_units(atlasMMP1,leadField.unit); clear atlas atlasMMP1.pos = leadField.pos; cfg = []; cfg.parcellation = 'parcellation'; cfg.method = 'mean'; cfg.parameter = 'cohspctrm'; srcConnPrcl = ft_sourceparcellate(cfg,srcConn,atlasMMP1); srcConnAvgPrcl = ft_sourceparcellate(cfg,srcConnAvg,atlasMMP1); clear cfg 
↑ Go up 
Now we can plot the results. 

%% Plotting the connectivity matrix between all pairs of parcels for i = 1:length(srcConnPrcl.freq)+1 if i <= length(srcConnPrcl.freq) figure; imagesc(srcConnPrcl.cohspctrm(:,:,i)); title([num2str(srcConnPrcl.freq(i)) ' Hz']); else figure; imagesc(srcConnAvgPrcl.cohspctrm); title('Average across frequencies') end colorbar; caxis([0 0.25]); end clear i % Plotting the connectivity between a pair of parcels as a function of frequency ch1 = 10; ch2 = 222; figure; plot(srcConnPrcl.freq,squeeze(srcConnPrcl.cohspctrm(ch1,ch2,:))) title(sprintf('Connectivity between %s and %s',srcConnPrcl.label{ch1},... srcConnPrcl.label{ch2}),'Interpreter','none'); xlabel('freq (Hz)') ylabel('Im(coherence)') clear ch1 ch2 
↑ Go up 
The imaginary part of coherency is a symmetric measure that does not take into account the direction of information flow between any pair of signals. To study directionality in connectivity, FieldTrip provides methods based on Granger causality. The interested reader can straightforwardly compute Granger causality by specifying cfg.method = ‘granger’. Before concluding this section, we will save the relevant variables. 

%% Saving relevant variables to disk save('workspaceSrcConn','srcConnPrcl','srcConnAvgPrcl') 
↑ Go up 
References

↑ Go up 