Converting rDCM estimation output in DCM structure for PEB analysis

[PasqualeDellaRosa]

Dear Stephan,
I have written you and Peter an email some time ago, as the topic concerned rDCM in a PEB modeling framework, but I’m not sure you have received it so I will post our issue as a rDCM topic on the TAPAS GitHub discussion forum as well, hoping that you can reach out back to us.

We are currently working on a rs-fMRI dataset and our idea is to take rDCM estimated models at the 1st level to a 2nd group-level PEB approach.
As you have mentioned previously on the forum, rDCM output does not store the information that is required for PEB, so we wrapped up some code to dress the rDCM estimated output with a DCM structure resembling the DCM model estimation output in order to get it to 2nd level PEB. I will summarize the steps we took

DCM Analysis: First Level Model Specification   

A template DCM model was specified for 3 Models 

•  Model 1: 53 ROIs
•  Model 2: 48 ROIs
•  Model 3: 54 ROIs

DCM specification parameters were the following: 

•  i) bilinear modulatory effects, 
•  ii) one state per region, 
•  iii) no stochastic effects, 
•  iv) no input centering
•  v) CSD fitting.
  The Template model was then input into the ‘specify group’ batch, to perform model specification for all subjects.       

rDCM 1st subject-level model estimation (i.e. tapas_rdcm_estimate.m).

rDCM model estimation output conversion to DCM structure for 2nd group-level PEB

 see “convert_rDCM_to_DCM.m code attached 

 PEB 2nd group-level (i.e. specify/estimate batch from SPM’s GUI)

PEB specify/estimate parameters

• DCM Index : 1          
• Covariates : none    
• Fields : 1x1 cell        
• Precision Components : Single      
• Within:between ratio : 0     
• Uncertainty : 0.0625 
• Group ratio : 1         
• Max Iterations : 8     
• Review PEB parameters : no          
  The max iterations parameter was set to 8 in order to overcome a convergence issue.

1. we included some comments in the code for each step in order for you to tell us if we are going in the right direction and if our assumptions for converting the rDCM to DCM structure are plausible or if we are missing something that should necessarily be taken in account when going from rDCM to DCM for PEB with respect to the theoretical assumptions grounding rDCM on one hand, and DCM/PEB on the other. 
2. you have pointed out in a thread to ” store the information that is required for PEB (i.e., posterior means and covariances, negative free energy) in the format that the PEB routines would expect. Note that rDCM by default only stores the region-wise covariance matrices and not the full covariance matrix to save storage space. To obtain the full covariance matrix (which you would need for PEB), you need to set options.evalCp = 1; where options is one of the input arguments to tapas_rdcm_estimate.m ”, however there is still other information missing so we wrote some lines to feed rDCM estimate_output Ep and Cp variables  + pE and pC priors. Could you eventually tell us  

•  If we have included all the variables in the DCM structure for PEB and if the artificial DCM structure we create resembles expected PEB input with all necessary variables for PEB
•  If the pE and pC calculation step sounds fit
•  If we should expect something different in the output from the classical DCM/PEB pipeline or if something may stand out in a different way.

Thanking you very much in advance for any kind of help and input you may give us. I believe that it would be great to combine rDCM at the subject level with PEB for group analyses, however I would really appreciate your advice before moving on to be completely sure that our approach and steps are fully plausible and hold or if we should instead move outside of the SPM-PEB framework and approach group-level stats outside SPM anyway on values extracted from rDCM single-subject estimated matrices.

If you should need the DCM input and DCM output matrices going in and coming out of our code for PEB just let me know where I could upload them for you, as the files are too heavy to attach to an email. I’ll attach the convert_rDCM_to_DCM.m code and rDCM estimate_output.mat here. 

Thanking you very much in advance for your help and support! 
All the best
Pasquale

rDCM_PEB_PADR.zip

[StefanFraessle]

Dear Pasquale,

Thank you very much for your interest in rDCM and my apologies for not having responded earlier to your message. The pre-Christmas time is always extremely busy and I simply did not have the time to take a closer look at your issue earlier.

Regarding your question: I think it is a fantastic idea to combine rDCM with SPM-PEB. In fact, I have already implemented a function that will transform the rDCM output structure into an SPM-compatible format that can then also be used with the standard SPM routines. This will also included the PEB routines. The function that does the transformation from rDCM to SPM output will be included in the rDCM toolbox in the next release of TAPAS.

I have quickly looked at your code. I did not have the time to go through it thoroughly, but I think it looks fine in principle. The only problem is that there was a small bug in the code for computing Cp out of the region-wise posterior covariance matrices. Specifically, the cross-covariance terms were not assigned correctly. This is not an issue for classical rDCM analyses (as you typically work with the region-wise covariance matrices), but it will lead to false estimates when using PEB (which requires the full Cp). My sincere apologies for this - it will be corrected in the new version of the rDCM toolbox.

I think the best way forward would be if you use the update of the code that I have written to convert rDCM to SPM format. I will double-check when we will have our next planed TAPAS release. If this is too far in the future and/or you need the script more urgently, I can also share the function with you privately prior to the release.

Apart from that, there is a couple of things to consider:

1. The creation of Cp from the region-wise covariance matrices takes quite a bit of time and can become cumbersome for whole-brain networks. I still have to figure out a way to speed this process up. In your case, with 50ish regions, this should be fine.
2. If I'm not mistaken, PEB will essentially create one big GLM comprising small GLMs for every model parameter that you enter into the PEB routine. This is already a relatively large model for normal DCMs. However, for the type of whole-brain models that you typically investigate with rDCM, the number of parameters is orders of magnitude higher. It is not clear to me yet, if this works at all with the standard PEB routines (for computational reasons) or whether one would need to (strongly) narrow down the number of parameters that one would like to take to the PEB level. This is certainly something that necessitates thorough testing.

I hope this is helpful. And once again, my apologies for the slow response.

With very best wishes,
Stefan