3dISC mixed effects modeling

Hello,

I am using 3dISC to model data from a mixed-effects design. I have 3 between-subject levels (patients with drug; patients with placebo; controls with placebo) and 4 within-subject levels (4 types of videos all subjects watched).

I am not interested in the between-group ISC and started by building a model disregarding the within-subject effects (averaging the pairwise correlation matrices from all videos using 3dCalc). Based on example 2 I came up with the following model:

3dISC -prefix ISC2b -jobs 12                  \
      -model  '0+grp+(1|Subj1)+(1|Subj2)'     \
      -gltCode ave     '0.333 0 0 0.333 0 0.333'            \
      -gltCode G11     '1 0 0 0 0 0'               \
      -gltCode G22     '0 0 0 1 0 0 '               \
      -gltCode G33     '0 0 0 0 0 1'               \
      -gltCode G11vG22 '1 0 0 -1 0 0'              \
      -gltCode G11vG33 '1 0 0 0 0 -1'              \
      -gltCode G22vG33 '0 0 0 1 0 -1'              \
      -dataTable                              \
      Subj1 Subj2    grp      InputFile       \
      s1     s2      G11      s1_2+tlrc       \
      ...          
      s1     s3      G12      s1_3+tlr       \
      ... 
      s1     s4      G13      s1_4+tlr       \
      ...
      s3     s6      G22      s3_6+tlr       \
      ...  
      s3     s4      G23      s3_4+tlr       \         
      ...
      s4     s7      G33     s4_s7+tlr       \

First I wanted to check if this looks correct or if the reverse group comparison constrasts need to be added:

-gltCode G11vG22 ‘-1 0 0 1 0 0’
-gltCode G11vG33 ‘-1 0 0 0 0 1’
-gltCode G22vG33 ‘0 0 0 -1 0 1’ \

Secondly, I wanted to asked if the best way to model for the within subject conditions is as described in Example 5.
For example, from those 4 types of videos, 3 types are social and 1 type is non-social. I wanted to look at the interaction effect of “group by socialness” and at the effect of group within social videos only.
Is there a way to run a single model with all of these effects or do I have to run separate ISC models (i.e. a model with all videos, as above; a model only with social videos; and a model for the socialvs.non-social).

Thanks for any input in advance!

Cheers,
Vasco

Vasco,

First I wanted to check if this looks correct or if the reverse group comparison constrasts need to be added

The original script looks fine and you don’t need to add those reversed versions because they are simply flipping the sign.

Is there a way to run a single model with all of these effects or do I have to run separate ISC models (i.e. a model with all
videos, as above; a model only with social videos; and a model for the socialvs.non-social).

You would have to create a separate model for each case.

Hi Gang,

Thank you for the clarifications. I understand that you similarly recommend that if using BML for this type of analysis, each condition contrast should be run separately for a computational time. Is there a recommended way to correct for the fact that I would be running 3 models (e.g. by adjusting quantiles when plotting and interpreting results)?

Also, I wanted to check if to adapt the MBA program for this analysis all I would have to do to adapt it is ensure that region-region (which I understand are ok if running the MBA analysis, but would become subject-subject under ISC and should not be included given I have between-subject variables?) and region-subject interactions are not included in the model.

Best regards,
Vasco

Hello again,

I have some additional questions after going through the MBA script from AFNI. I suppose this question might be more adequate for the stan forums, but if you could help that would be much appreciated. I came up with the following model for my analysis:

Y ~ 1 + group_dummy1 + group_dummy2 + (1 | ROI) + (1 + group_dummy1 + group_dummy2 | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Would this be correct for my situation? My understanding is that I could not run the model with a single grouping variable, but I don’t understand why not (example below).

Y ~ 1 + group + (1 | ROI) + (1 + group | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Similarly, if I were to include a variable like age and wanted to look at the age by group interaction (let’s assume two groups for simplicity here), from the MBA documentation I see that I should code a variable for that (group_age, corresponding to the first model below), but I don’t understand why it could not be coded as the second model below:

First model:
Y ~ 1 + group + age + group_age + (1 | ROI) +
(1 + group + age+ group_age| mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Second model:
Y ~ 1 + groupage + (1 | ROI) +
(1 + group
age | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Finally, if I wanted to include a within-subject factor of socialness (social and non-social videos), where I want to look at the ISC within each condition, the difference between conditions within each group, and the interaction of group by socialness, would the model below be the correct one (once again assuming 2 groups only)? Or should I include socialness by itself as well (i.e. group*socialness)? Or should I really stick with running separate models for each of the within-subject contrasts?

Y ~ 1 + group + group:socialness + (1 | ROI) +
(1 + group + group:socialness | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Thank you so much for your help again!

Best,
Vasco

Is there a recommended way to correct for the fact that I would be running 3 models (e.g. by
adjusting quantiles when plotting and interpreting results)?

Exploratory studies are usually hard to reach a decisive conclusion. Thus, I think that one should focus on effect estimation, not decision making. In light of this consideration, I would not worry about the issue as to whether to make adjustment for the 3 separate models.

For MBA-type analysis, do you mean something like the following?

https://www.sciencedirect.com/science/article/pii/S1053811919310651

a within-subject factor of socialness (social and non-social videos)

If you could obtain the ISC between the two levels of socialness, it would be much easier to handle from the perspectives of modeling and computational cost.

Hello Gang,

Thanks for the clarification regarding the three models.
Yes I can easily obtain the ISC between the two levels, so I would just do that. I’m still curious what a model with such a factor incorporated would look like.

And yes, the paper you attached is precisely what I am trying to accomplish, so I was wondering if the adaptation i mentioned to MBA would the correct one (i.e. remove the subject:subject interaction).

Cheers,
Vasco

Yes I can easily obtain the ISC between the two levels, so I would just do that. I’m still curious what a model
with such a factor incorporated would look like.

Maybe something like:

Y ~ groupsocialness + (groupsocialness | ROI) + (socialness |Subj1:Subj2) + (socialness | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

the paper you attached is precisely what I am trying to accomplish, so I was wondering if the adaptation i
mentioned to MBA would the correct one (i.e. remove the subject:subject interaction).

Try this:

Y ~ groupage + (groupage | ROI) + (1 | Subj1:Subj2) + (1 | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))

Let me know if this makes sense to you.

Hello Gang,

Thank you so much it all make sense. I still have two small questions.

1 - If I am including age, is the age the sum of ages of the 2 participants a given r corresponds to? Am I also correct in assuming that this summed age needs to be mean-centered across all pairs?
2- In my case I have 3 groups, and not interested in cross-group ISC, am I correct in assuming I can model this with just 2 dummy variables for this?

Best,
Vasco

If I am including age, is the age the sum of ages of the 2 participants a given r corresponds to? Am I also
correct in assuming that this summed age needs to be mean-centered across all pairs?

Center the age first and then sum each pair.

In my case I have 3 groups, and not interested in cross-group ISC, am I correct in assuming I can model
this with just 2 dummy variables for this?

Yes, you need 2 variables to code the three groups. Use the contrast coding of 1, -1 and 0.

Hi Gang,

Thank you so much. I coded the model as you specified with the 3 groups and was able to a pillot test on a subset of my sample.
Now I’m working on extracting effect estimates, summary info, plots, etc.
As I was going through the MBA code, I noticed that when summarizing the information for the posterior samples of Subject pairs, you apply some sort of standardization to all values based on the mean of all samples for that subject and the square root of 2, per below:

for(ii in 1:nSubjs) for(jj in 1:nSubjs) ps[,ii,jj] ← sqrt(2)*(ps[,ii,jj] - mm[ii,jj]) + mm[ii,jj]

I am wondering why there was a need to do this, as I have never come across this formula.

Best,
Vasco

I noticed that when summarizing the information for the posterior samples of Subject pairs, you apply some sort of
standardization to all values based on the mean of all samples for that subject and the square root of 2

Vasco, are you stealing my secret codes? :slight_smile: The rationale for that adjustment is due to the fact the data are used twice in the whole process.

Thank you for the clarification! I am, as you suggested, running the within-subject effects as separate models. However, I still have a question regarding this. Let assume for simplicity thay i have 3 r values per subject (rA, rB, and rC).

Firstly, if i wanted the overall effect regardless of video, should i input all 3 values or take the mean after r-to-z transformation and use that as input (i.e. rInput = (rA+rB+rC)/3) ?

Secondly, if I wanted a contrast of condition X > condition Y where rA and rB belong to conditon X and rC belondgs to condition Y, would it be correct to do (rA+rB)/2) - rC and use the resulting value as input.

Thank you once again for all your help, it has been invaluable!

Cheers,
Vasco

if i wanted the overall effect regardless of video, should i input all 3 values or take the mean after r-to-z
transformation and use that as input (i.e. rInput = (rA+rB+rC)/3) ?

That seems reasonable to me.

if I wanted a contrast of condition X > condition Y where rA and rB belong to conditon X and rC belondgs
to condition Y, would it be correct to do (rA+rB)/2) - rC and use the resulting value as input.

I’m not familiar with the context. If you believe that (rA+rB)/2) represents condition X just as rC does condition Y, then your suggestion may be OK.

Hi,

Thanks for confirming.

I am still wondering if it would be incorrect to use all 3 values as input in the first example. I ask because my intuition is that 3 values would be better than the mean of those 3 values. Is there any advantage to using the mean besides shorter computation time?

Vasco

I am still wondering if it would be incorrect to use all 3 values as input in the first example.

It is not about correctness. Rather, it’s more about the convenience of model specification and computational cost. For instance, you can run

-model ‘grp*video+(0+video|Subj1)+(0+video|Subj2)’ \

The gltCode specification would become more complex.

Hi Gang,

I am a bit confused with the group contrasts. Let’s say I have 3 groups (Control, Treatment1 and Treatment2).

Let’s assume my model is Y ~ group + (group | ROI) + (1 |Subj1:Subj2) + (1 | mm(Subj1, Subj2, weights = cbind(1, 1), scale = FALSE))
Group is coded as Control = 0, Treatment1 = 1, and Treatment2 = -1.

Now if I am understanding correctly, in my output I will have posterior estimates for Intercept (the mean of the 3 groups), Group1 (Overall mean of the 3 groups minus the mean of Group2 and Group0) , and Group2 (Overall mean of the 3 groups minus the mean of Group0 and Group1). From these I can then calculate the mean of each group individually and other differences. However, I am unclear on how I match the Group1 and Group2 to the -1,0,1 coding.

Is this interpretation correct or do the Groups in the output mean something else with this coding?

Cheers,
Vasco

Vasco,

I assume that 1) you’re writing your code; and 2) you’re quantitatively coding the three groups with -1, 0, and 1.

To properly handle three groups, you need two variables (plus an intercept). The following example demonstrates the coding strategy (with group3 as the base/reference level):

                 G1        G2

group1 1 0
group2 0 1
group3 -1 -1

With the above coding, the intercept is the average of the three groups, G1 shows the difference between group1 and the group average, and G2 indicates the difference between group2 and the group average.