Print

Print


Dear Peter:

 

Thank you so much for your help with the PEB model explanation. I still have a few questions:

 

1.      When I looked at the BMA structure, it has a field “Ep” and “Pp”. Based on a previous post, “Ep” should be the beta parameters of the GLM model, right? How about “Pp”? Is “Pp” the probability of “Ep”?

 

2.      In my case, BMA.Pp only has 3 values: 0, 0.5 and 1.0. If Pp is 0, does it mean the connection is switched off (removed)? Why Pp does not have any other values if it is a probability?

 

3.      When I review the BMA results using “spm_dcm_peb_review(BMA)” and choose group difference (covariate 2), what does the figure show? Does it show the BMA.Ep (beta values) for group difference?

 

4.      When I set the threshold to 0.5 (Pp>0.5), the figure (attached) shows that all (or almost all) parameters (connections) have a value. I am confused because it should remove the parameters with Pp = 0, right?

 

5.      When I select the threshold to be greater than 0.99 (Pp > .99), many parameters (connections) still survive (see attached figure). Is it reasonable? Does it mean Bayesian model reduction does not remove many parameters (connections)?

 

6.      How does Bayesian model reduction work? Should it remove the connections with very small BMA.Ep (because Ep represents the effects of group difference on the specific connection)? Or it remove the connections with very small BMA.Pp?

 

7.      How to see the group averaged connection strength (the connection strength that takes the covariance into account)?

 

8.      If I want to see the results of the full model (before reduction), should I use “spm_dcm_peb_review(PEB)”?

 

Thank you very much again!

 

Best regards,

Guoshi Li

 

 

From: Zeidman, Peter <[log in to unmask]>
Sent: Monday, January 7, 2019 5:50 AM
To: Li, Guoshi <[log in to unmask]>
Cc: [log in to unmask]
Subject: RE: [SPM] Running time issue of PEB model

 

Dear Guoshi

 

Thank you very much for your kind response and your suggestion works. After changing the PEB setting from M.Q='all' to M.Q='single', the PEB model estimation only took a few minutes and the nested PEB model search took about 2.5 hours to complete (in a local computer). Running in the cluster is not the issue because it also took forever to run if I run it in my local computer with M.Q='all'. So it makes a huge difference between M.Q='single' and M.Q='all'.

 

Great.

 

Something quite interesting is that it only estimate 725*2 parameters (instead of 729*2). After examining the field PEB.Pnames, I found four connections are missing: A(22, 22), A(27, 22), A(22, 27) and A(27, 27). Do you know why?

 

If you have a large number of regions in a DCM for fMRI, it will automatically use the functional connectivity as prior constraints for the effective connectivity. This effectively binds together connections which are highly correlated. I expect that after doing this, some of your parameters have sufficiently small prior variance that PEB treats them as having been switched off. For more detail on the functional connectivity parameters, see:

 

https://doi.org/10.1016/j.neuroimage.2012.12.005

https://doi.org/10.1162/NETN_a_00015 

 

Also, after I have the results (both PEB and BMA), how do I know which connections are significantly different between the two groups (sorry that I still don’t fully understand the PEB method).

 

The second level of the PEB model is simply a General Linear Model (GLM) with regressors corresponding to the group mean of the connections and the group difference. When you called spm_dcm_peb_bmc, it pruned the model, removing any parameters (weights or betas) which did not contribute. The results were stored in the structure named BMA. To review this, type:

 

spm_dcm_peb_review(BMA);

 

The drop-down menu lets you view the parameters relating to the group mean and the group difference. If you click on the bars for the group difference, it will tell you the probability for each connection exhibiting a group difference. There is no concept of ‘significance’ in Bayesian stats, but you might want to focus your discussion on parameters with probability greater than 0.9 or 0.95. Feel free to ask if you have further questions.

 

Best

Peter

 

 

-----Original Message-----
From: Zeidman, Peter <[log in to unmask]>
Sent: Saturday, January 5, 2019 4:23 PM
To: Li, Guoshi <[log in to unmask]>; [log in to unmask]
Subject: RE: [SPM] Running time issue of PEB model

 

PS on further thought - you have 2 covariates in your second level design matrix and a fully connected matrix of 27x27=729 connections. So that's 2*729=1,458 regressors / parameters in your second level GLM with just 200 subjects. I think you probably don't have enough degrees of freedom to estimate this model efficiently, and you may also be running out of RAM for the matrix inversions on your cluster. So as suggested, try running it locally. And ensure you set in the PEB settings:

 

M.Q='single' ;

 

Rather than M.Q='all' - otherwise you'll have a further 729 parameters quantifying the between-subject variability on each connection.

 

Best

Peter

 

-----Original Message-----

From: Zeidman, Peter

Sent: 05 January 2019 21:10

To: 'Guoshi Li' <[log in to unmask]>; [log in to unmask]

Subject: RE: [SPM] Running time issue of PEB model

 

Dear Guoshi

If your PEB analysis has taken 6 days, then something is wrong. It should take a few seconds! Minutes at most. The only exception is when you use the function spm_dcm_peb_fit instead of spm_dcm_fit, which iteratively re-estimates all subjects' DCMs using the group-average connectivity as priors.  If you're doing that, then just switch back to spm_dcm_fit.

 

Perhaps you could try the analysis on your local computer rather than the cluster?

 

Best

Peter

 

-----Original Message-----

From: SPM (Statistical Parametric Mapping) <[log in to unmask]> On Behalf Of Guoshi Li

Sent: 02 January 2019 16:46

To: [log in to unmask]

Subject: [SPM] Running time issue of PEB model

 

Dear Peter:

 

I used spectral DCM to estimate the effective connectivity (EC) of a pretty large model (27 nodes) of resting-state fMRI for two groups of subjects (each group has 100 subjects). I have completed the estimation of all 200 subjects (full connectivity) and now I need to compare the EC difference between the two groups. Previously I used traditional t-test in conjunction with network-based statistic (NBS; Zalesky et al., NeuroImage, 2010) to identify the significant EC links (edges) between the two groups. After reading your paper (Friston et al., NeuroImage, 2016) and some posts of this forum, I know there is an alternative approach (Bayesian Model Reduction with PEB) to compare the group difference. I am following the instructions in the Wikibooks (https://en.wikibooks.org/wiki/SPM/Parametric_Empirical_Bayes_(PEB)) and your previous responses related to the PEB method. Now the problem I have is the running time issue. I am running the PEB code in a cluster (SPM version 7219) and it has not been completed after 6 days. I can see a temporal data file (tmp.mat) generated and no error message was given (the PEB data file is not produced yet). I am not sure whether it is because I have too many nodes (27) and subjects (200) in my model or because there is something wrong in my handling of the PEB model. In your view, what is the reasonable time frame to perform the PEB analysis on such a large model? Below are my steps/codes.

 

Thank you very much for your help and best wishes for the new year!

 

Best regards,

Guoshi Li

 

___________________________________________________________________________________________

 

1)           Estimate the EC of 200 DCMs with full connectivity using spectral DCM approach

 

2)           Assemble all DCMs in a GCM cell array (I used DCM structure directly instead of mat file)

GCM = { DCM1

                DCM2

                 ⁞

                DCM200 }

 

3)           Estimate a second level PEB model

 

% Specify PEB model settings

M = struct();

M.alpha = 1;

M.beta  = 16;

M.hE    = 0;

M.hC    = 1/16;

M.Q     = 'all';

 

% Specify design matrix for N subjects.

M.X(:,1) = ones(200,1);

M.X(:,2) = [ones(100, 1); -1*ones(100, 1)];

 

% Choose field

field = {'A'};

 

% Estimate model

PEB  = spm_dcm_peb(GCM,M,field);

save('PEB_Model.mat','PEB');

 

4)           Search over nested PEB models

BMA = spm_dcm_peb_bmc(PEB(1));

save('BMA_Model.mat','BMA');