[FieldTrip] Inquiry About Default Model in Granger Causality Analysis
Schoffelen, J.M. (Jan Mathijs)
janmathijs.schoffelen at donders.ru.nl
Mon Jul 22 12:55:23 CEST 2024
Hi Jeff,
you wrote that you want to integrate across the frequency range because you intend to obtain time-domain Granger causality. Why taking the hassle to go through the frequency domain, and not just take the log ratio of the residuals of two appropriately defined MVAR models? At least, that’s what I would do.
Best wishes,
Jan-Mathijs
On 7 Jul 2024, at 08:34, Jeff via fieldtrip <fieldtrip at science.ru.nl> wrote:
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<mailto: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<mailto: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%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=sg71ZeIIafWd6E%2F5MY2LBus1ocVFSm05vDlXOK6zn1g%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%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=0sSeZIQt84BJChw0a6KmAX1didW%2F5DRKOMUB2oqDpLw%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://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fmailman.science.ru.nl%2Fmailman%2Flistinfo%2Ffieldtrip&data=05%7C02%7Cfieldtrip%40science.ru.nl%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=g4mtp2vClwZAjq6YXfF3W07rz%2FwTxTgs%2BSoiVEDlJDg%3D&reserved=0>
https://doi.org/10.1371/journal.pcbi.1002202<https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fdoi.org%2F10.1371%2Fjournal.pcbi.1002202&data=05%7C02%7Cfieldtrip%40science.ru.nl%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=GmYVIQ5TUOltas9TUlNd3yCvwZ9UGjnCFA7nK47zZ8k%3D&reserved=0>
_______________________________________________
fieldtrip mailing list
https://mailman.science.ru.nl/mailman/listinfo/fieldtrip<https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fmailman.science.ru.nl%2Fmailman%2Flistinfo%2Ffieldtrip&data=05%7C02%7Cfieldtrip%40science.ru.nl%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=g4mtp2vClwZAjq6YXfF3W07rz%2FwTxTgs%2BSoiVEDlJDg%3D&reserved=0>
https://doi.org/10.1371/journal.pcbi.1002202<https://eur01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fdoi.org%2F10.1371%2Fjournal.pcbi.1002202&data=05%7C02%7Cfieldtrip%40science.ru.nl%7C6c9f6323074d4275689208dcaa3cc9ab%7C084578d9400d4a5aa7c7e76ca47af400%7C1%7C0%7C638572425255159702%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C0%7C%7C%7C&sdata=GmYVIQ5TUOltas9TUlNd3yCvwZ9UGjnCFA7nK47zZ8k%3D&reserved=0>
_______________________________________________
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/20240722/0cc9365a/attachment.htm>
More information about the fieldtrip
mailing list