coh_all.m

Find all cohere method

function [PDC, GPDC, PDCF] = coh_all (TS, N, Fs, P, title)
% Brain Dynamics Lab's coherence analysis tool
% 捞 橇肺弊伐篮 Biosig 4 啊 鞘荐利栏肺 鞘夸钦聪促.
% Biosig4 狼 plota function 阑 捞侩窍咯 累己茄 巴涝聪促.
% --------------------------------------------------------------
% Examples:
% [PDC, GPDC, GGC, pCOH2] = find_coh (EEG, 250, 100, 5);
%   => plot 4 figures
% --------------------------------------------------------------
% INPUT : 
% -------------------------------------------------------------
% A, B	multivariate polynomials defining the transfer function
%
%    a0*Y(n) = b0*X(n) + b1*X(n-1) + ... + bq*X(n-q)
%                          - a1*Y(n-1) - ... - ap*Y(:,n-p)
%
%  A=[a0,a1,a2,...,ap] and B=[b0,b1,b2,...,bq] must be matrices of
%  size  Mx((p+1)*M) and Mx((q+1)*M), respectively. 
%
%  C is the covariance of the input noise X (i.e. D'*D if D is the mixing matrix)
%  N if scalar, N is the number of frequencies 
%    if N is a vector, N are the designated frequencies. 
%  Fs sampling rate [default 2*pi]
% --------------------------------------------------------------
% OUTPUT : 
% --------------------------------------------------------------
% PDC 	partial directed coherence [2]
% pCOH2 partial coherence - alternative method 
% GGC	a modified version of Geweke's Granger Causality [Geweke 1982]
%	   !!! it uses a Multivariate AR model, and computes the bivariate GGC as in [Bressler et al 2007]. 
%	   This is not the same as using bivariate AR models and GGC as in [Bressler et al 2007]
% GPDC 	Generalized Partial Directed Coherence [9,10]
%
% see also: FREQZ, MVFILTER, MVAR
% 
% REFERENCE(S):
% [1] H. Liang et al. Neurocomputing, 32-33, pp.891-896, 2000. 
% [2] L.A. Baccala and K. Samashima, Biol. Cybern. 84,463-474, 2001. 
% [3] A. Korzeniewska, et al. Journal of Neuroscience Methods, 125, 195-207, 2003. 
% [4] Piotr J. Franaszczuk, Ph.D. and Gregory K. Bergey, M.D.
% 	Fast Algorithm for Computation of Partial Coherences From Vector Autoregressive Model Coefficients
%	World Congress 2000, Chicago. 
% [5] Nolte G, Bai O, Wheaton L, Mari Z, Vorbach S, Hallett M.
%	Identifying true brain interaction from EEG data using the imaginary part of coherency.
%	Clin Neurophysiol. 2004 Oct;115(10):2292-307. 
% [6] Schlogl A., Supp G.
%       Analyzing event-related EEG data with multivariate autoregressive parameters.
%       (Eds.) C. Neuper and W. Klimesch, 
%       Progress in Brain Research: Event-related Dynamics of Brain Oscillations. 
%       Analysis of dynamics of brain oscillations: methodological advances. Elsevier. 
% [7] Bressler S.L., Richter C.G., Chen Y., Ding M. (2007)
%	Cortical fuctional network organization from autoregressive modelling of loal field potential oscillations.
%	Statistics in Medicine, doi: 10.1002/sim.2935 
% [8] Geweke J., 1982	
%	J.Am.Stat.Assoc., 77, 304-313.
% [9] L.A. Baccala, D.Y. Takahashi, K. Sameshima. (2006) 
% 	Generalized Partial Directed Coherence. 
%	Submitted to XVI Congresso Brasileiro de Automatica, Salvador, Bahia.  
% [10] L.A. Baccala, D.Y. Takahashi, K. Sameshima. 
% 	Computer Intensive Testing for the Influence Between Time Series, 
%	Eds. B. Schelter, M. Winterhalder, J. Timmer: 
%	Handbook of Time Series Analysis - Recent Theoretical Developments and Applications
%	Wiley, p.413, 2006.
% [11] M. Eichler
%	On the evaluation of informatino flow in multivariate systems by the directed transfer function
%	Biol. Cybern. 94: 469-482, 2006. 	

% -----------------------------------------------------------------
% Contact : KAIST, Brain Dynamics LAB , Chae Yong Wook
%           [e-mail]

      [AR,RC,PE] = mvar(TS,P);
      M = size(AR,1); % number of channels       

    X.A = [eye(M),-AR];
    X.B = eye(M); 
    X.C = PE(:,M*P+1:M*(P+1)); 
    X.datatype = 'MVAR';

    [S,h,PDC,COH,DTF,DC,pCOH,dDTF,ffDTF, pCOH2, PDCF, coh,GGC,Af,GPDC,GGC2]=mvfreqz(X.B,X.A,X.C,N,Fs);

%     F = 1;
% %     Display PDC
%     figure(F); F = F+1; 
% %     figure;
% % %     figure(
%     [H, PDC] = plota(X,'PDC',N,Fs);         
%     tmp = sprintf('PDC: Figure ');
% %     tmp = title;
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
% % %     
%     figure(F); F = F+1; 
%     [H, GPDC] = plota(X,'GPDC',N,Fs);         
%     tmp = sprintf('GPDC: Figure ');
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
%     Display DTF
%     figure(F); F = F+1; 
%     plota(X,'DTF',N,1);         
%     tmp = sprintf('DTF: Figure ');
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
%     % Display PC
%     figure(F); F = F+1; 
%     [H, pCOH2] = plota(X,'pCOH2',N,Fs);         
%     tmp = sprintf('pCOH2: Figure ');
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
%     figure(F); F = F+1; 
%     plota(X,'COH',N,1);         
%     tmp = sprintf('COH: Figure ');
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
%     figure(F); F = F+1; 
%     [H, GGC] = plota(X,'GGC',N,Fs);         
%     tmp = sprintf('GGC: Figure ');
%     if exist('suptitle','file')
%             suptitle(tmp);        
%     end;
end