[FieldTrip] Inquiry About Default Model in Granger Causality Analysis
Jeff
wcy131608 at gmail.com
Sun Jul 7 08:34:00 CEST 2024
Dear Jan-Mathijs,
Thank you very much for your detailed response to my previous query. I
apologize for the delay in getting back to you as I have been supplementing
my understanding of Granger causality in this field.
After further study, I found that my task-related data is non-stationary,
as determined by the Augmented Dickey-Fuller unit-root test and the KPSS
unit-root test. I have come across literature suggesting that using wavelet
transform to compute the cross-power spectrum at each time-frequency point,
followed by non-parametric Granger causality analysis, is suitable for
non-stationary data (Dhamala, M., Rangarajan, G., & Ding, M. (2008)). Based
on this, and after reviewing past community emails, I decided to use the
following code for my task-related data non-parametric Granger causality
analysis:
cfg = [];
cfg.channel = 'all';
cfg.method = 'tfr';
cfg.output = 'powandcsd';
cfg.foi = 3:1:250; %alternative method:
cfg.foi = 0.02:1:250;
cfg.toi = [-0.6:0.02:-0.4,0:0.02:2.7]; % epoch time range [-1 3]
cfg.width = 3;
cfg.pad = 'nextpow2';
subj_data_tf = ft_freqanalysis(cfg, ft_data);
subj_data_tf.freq= round(subj_data_tf.freq);
grangercfg = [];
grangercfg.method = 'granger';
grangercfg.granger.conditional = 'no';
grangercfg.granger.sfmethod = 'bivariate';
gdata = ft_connectivityanalysis(grangercfg, subj_data_tf);
However, I learned from your previous response that the frequency axis
should start from 0 and consist of integer elements. Given my sampling rate
is 500 Hz, ideally, my frequency axis should be set to 0:1:250 because I
need to integrate over the frequency axis to obtain time-domain Granger
causality. I found that cfg.method = 'tfr' does not support setting the
lowest frequency to 0 Hz. I attempted to set the frequency axis to
0.02:1:250.02 and used subj_data_tf.freq = round(subj_data_tf.freq) after
wavelet transformation to approximate the ideal 0:1:250 frequency axis.
This rough approximation required zero-padding the epoch at both ends to
estimate the cross-power spectrum at the boundary time-frequency points.
However, I discovered that zero-padding at the epoch ends affects the
estimates for all time-frequency points within the epoch, regardless of
their proximity to the boundaries. Therefore, I am uncertain whether this
approach is valid, and this forms my first question.
My second question pertains to the frequency axis configuration. I
understand that if I set the frequency axis to cfg.foi = 3:1:250, the
ft_connectivityanalysis function will automatically pad the frequency axis
to 0:1:250 and add zero matrices for the cross-power spectra at 0 Hz, 1 Hz,
and 2 Hz. Which method is more accurate: this automatic padding method or
the method involving zero-padding the epoch data and setting the frequency
axis to 0.02:1:250.02 for wavelet transformation?
Thank you once again for your time and assistance.
Best regards,
Chengyuan Wu
Schoffelen, J.M. (Jan Mathijs) via fieldtrip <fieldtrip at science.ru.nl>
于2024年6月15日周六 02:28写道：
> Hi Jeff,
>
> I agree that the documentation regarding several ways to compute
> transfer-function derived connectivity metrics is a bit sparse.
>
> On 13 Jun 2024, at 15:30, Jeff via fieldtrip <fieldtrip at science.ru.nl>
> wrote:
>
> To whom it may concern,
>
> I hope this message finds you well. I am writing to seek assistance
> regarding the computation of frequency-domain Granger causality using
> non-parametric methods as described in the FieldTrip tutorial
> documentation. Specifically, I am interested in the section on conditional
> Granger causality (
> https://www.fieldtriptoolbox.org/example/connectivity_conditional_granger/
> <https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.fieldtriptoolbox.org%2Fexample%2Fconnectivity_conditional_granger%2F&data=05%7C02%7Cfieldtrip%40science.ru.nl%7Cdfcfff4cecf443e0b05808dc8c9bbe5e%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638539847739380313%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=001av91qmwyS7Xfgokx2rno7%2FIG8H0EwfPjEgWu3Wmg%3D&reserved=0>
> ).
>
> In the tutorial, the code explicitly specifies whether to use a bivariate
> or multivariate model, whereas other sections of the documentation (e.g.,
> https://www.fieldtriptoolbox.org/tutorial/connectivityextended/#introduction
> <https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.fieldtriptoolbox.org%2Ftutorial%2Fconnectivityextended%2F%23introduction&data=05%7C02%7Cfieldtrip%40science.ru.nl%7Cdfcfff4cecf443e0b05808dc8c9bbe5e%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638539847739380313%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=wIJqzH6RAekNXNhpvtSnGR5X8B0yrAYCECpQ2T09I9Q%3D&reserved=0>)
> do not explicitly state this. My question is: when the model type
> (bivariate or multivariate) is not explicitly stated, which model does
> FieldTrip use by default?
>
> If you seek to do non-parametric spectral factorization on a data
> structure with N-channels -> multivariate, because the default of
> cfg.sfmethod = ‘multivariate'
> If you seek to do parametric estimation (where the data goes through
> ft_mvaranalysis, which by default does a multivariate AR model) -> also
> multivariate
>
> To investigate this, I conducted a small test. I found that specifying
> "cfg.granger.sfmethod = 'bivariate';" produced different directed
> connectivity results etween electrode pairs compared to when
> "cfg.granger.sfmethod" was not specified. This test was conducted using the
> non-parametric Granger method. However, when using the parametric method as
> the tutorial, I found that specifying "cfg.granger.sfmethod = 'bivariate';"
> in "ft_connectivityanalysis" resulted in an error.
>
>
> This is because sfmethod parameter’s value refers to how the spectral
> factorization is to be performed. This is essentially non-parametric, so
> the combination and cfg don’t make sense in your particular test case.
>
> A ‘bivariate’ decomposition of the spectral matrix, or equivalently a
> bivariate AR-model (i.e. estimated per pair of channels, if ft_mvaranalysis
> is supplied with a cfg.channelcmb field) will yield a different result for
> the estimated connectivity between channels A and B, as compared to a
> situation where the spectral matrix decomposition was done, with - say - a
> channel C present (or equivalently, when the AR-model is estimated in a
> multivariate sense). I think that there are ample explanations about this
> in the literature, but the gist of it is that any latent sources (e.g.
> source C not modelled along in the bivariate case) will influence the
> actual estimate between A and B.
>
> On the other hand, the higher the dimensionality of the decomposition
> challenge at hand (i.e. # of channels, often poorly unmixed) the less
> informative the spectral transfer functions / AR-model coefficients (or
> often the model does not converge meaningfully at all).
>
> Best wishes,
> Jan-Mathijs
>
> Here is the code for my non-parametric Granger causality test:
>
> cfg.toilim = [0.5 1.5];
> data = ft_redefinetrial(cfg, data);
>
> cfg = [];
> cfg.method ='mtmfft';
> cfg.taper = 'dpss';
> cfg.output = 'fourier';
> cfg.tapsmofrq = 2;
> mfreq = ft_freqanalysis(cfg, data);
>
> cfg = [];
> cfg.method = 'granger';
> cfg.granger.conditional = 'no';
> cfg.granger.sfmethod = 'bivariate';
> mgranger_bivar = ft_connectivityanalysis(cfg, mfreq);
>
> cfg = [];
> cfg.method = 'granger';
> mgranger_default = ft_connectivityanalysis(cfg, mfreq);
>
>
> Additionally, here is the code using the parametric Granger method, which
> does not allow explicit specification of "cfg.granger.sfmethod =
> 'bivariate';":
>
> cfg = [];
> cfg.order = 5;
> cfg.toolbox = 'bsmart';
> mdata = ft_mvaranalysis(cfg, data);
>
> cfg = [];
> cfg.method = 'mvar';
> mfreq_para = ft_freqanalysis(cfg, mdata);
>
> cfg = [];
> cfg.method = 'granger';
> % cfg.granger.conditional = 'no'; %error
> % cfg.granger.sfmethod = 'bivariate';%error
> granger_para = ft_connectivityanalysis(cfg, mfreq_para);
>
> I would greatly appreciate any guidance or recommendations on how to
> address this issue and understand the default model used by FieldTrip when
> not explicitly specified. My core question is whether the traditional
> bivariate model is used by default when the model type (bivariate or
> multivariate) is not specified, regardless of whether a parametric or
> non-parametric Granger method is employed. If the default is indeed the
> bivariate model, why do the results differ from those obtained when
> explicitly specifying the model as bivariate(with non-parametric Granger
> method)? Additionally, any insights on resolving the error encountered when
> specifying the bivariate model in the parametric method would be extremely
> helpful.
>
> Thank you very much for your time and assistance.
>
> Best regards,
>
> Chengyuan Wu
>
>
> _______________________________________________
> fieldtrip mailing list
> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
> https://doi.org/10.1371/journal.pcbi.1002202
>
>
> _______________________________________________
> fieldtrip mailing list
> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
> https://doi.org/10.1371/journal.pcbi.1002202
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20240707/635f5150/attachment.htm>
More information about the fieldtrip
mailing list