3dLME to evaluate a group effect on seed base analysis

Dear Gang,
I have produced several correlation maps (SBA) for each of my subjects,
Sometime, but not always I have different sessions per subject as well as runs
I would like to measure the overall effect as well as the inter individual variability

I have tried

3dLME \
-prefix 3dLME_glt.nii.gz \
-jobs 20 \
-mask mask_mean_func_overlapp.nii.gz \
-model 'Sess*run' \
-ranEff '~1' \
-num_glt 2 \
-gltLabel 1 'Group_Effect' \
-gltCode 1 'Sess : 1' \
-gltLabel 2 'Interindividual_Variability' \
-gltCode 2 'Subj : 1' \
-dataTable @design_matrix.txt \
-resid resid.nii.gz

but i am really not sure of what I am doing
sor far I got:


Read 100 items
Error in seq.default(2, length(sepTerms), 2) : 
  wrong sign in 'by' argument
Calls: process.LME.opts -> gl_Constr -> glfConstr -> seq -> seq.default
Execution halted

not sure why?
design_matrix.txt look like

Subj Sess run InputFile
A 1 0 pathtodata.nii.gz

Thank you so much for your help
Clément

Clément,

I assume that both variables, Sess and run, are within-individual factors. In this case, the 3dLMEr program is more appropriate. I recommend following this blog post for guidance on model construction.

A few specific points regarding your effect specification:

  1. The following specification:

    -gltCode 1 'Sess : 1' \  
    

    should be corrected to:

    -gltCode 1 'Sess : 1*1' \  
    
  2. The intent behind this specification is unclear, but as written, it does not make sense:

    -gltCode 2 'Subj : 1' \  
    

Let us know if you need further clarification.

Gang Chen

Hey Gang,
Thank you so much for you reply!
It seems to allow the task to run now but I ran into another problem =(

anyway,
just few question before we found a way to work around this other problem,

if gltCode 2 is not meaningful for 3dLME,
I am not sure which glt would test for:

Inter-individual variability?
Inter-session variability?
Inter-run variability?

I also read the doc you send to me, however, I can't test 3dLMEr for the same reason mention previously...

However, would

3dLMEr \
-prefix 3dLME_glt.nii.gz \
-jobs 20 -mask mask_mean_func_overlapp.nii.gz \
-model "Sess*run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)" \
-dataTable @/design_matrix.txt \
-resid resid.nii.gz

be the way to go for testing my 3 hypotheses?
Thank you again!

Clément

Hi Clément,

I'm not very familiar with installation issues involving a docker container. Hopefully, someone more knowledgeable can assist.

Regarding your model for 3dLMEr:

-model "Sess*run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)"  

this specification includes four variance components:

  1. Cross-individual variance: (1|Subj)
  2. Cross-session-individual variance: (1|Sess:Subj)
  3. Cross-run-individual variance: (1|run:Subj)
  4. Residual variance (cross-session-run-individual)

If you're looking to extract these variance components, 3dLMEr currently does not provide an option for such extraction.

Gang Chen

Hey Gang,
Thank you, now that I solved the docker issue I manage to run it !
However, I seems that I am missing the group effect, that takes into account run and session as within subject variable
It is like the output of 3dttest against zero, to see where does the activation occur
Thank you again!
Clément

Clément,

I seems that I am missing the group effect, that takes into account run and session as within subject variable

Are you looking to extract the overall average from the model? If you have two runs coded as r1 and r2, try adding the following line:

-gltCode group.mean 'run : 0.5*r1 +0.5*r2' \  

Let us know if this works for you.

Gang Chen

thank you,
and if I add:

-gltCode group.mean 'Sess: 0.5*1 +0.5*2' \

This will be the Session effect taking into account the different runs?
Also, If I have way more session 01 than 02, should I write something like

-gltCode group.mean 'Sess: 0.8*1 +0.2*2' \

Thank you again

if I add:
-gltCode group.mean 'Sess : 0.5*1 +0.5*2' \
This will be the Session effect taking into account the different runs?

Yes, this would give the average effect across both sessions and runs. In other words, it would be equivalent to:

-gltCode group.mean2 'run : 0.5*1 +0.5*2' \

If I have way more session 01 than 02, should I write something like
-gltCode group.mean 'Sess : 0.8*1 +0.2*2' \

No, the unequal number of samples is already accounted for in the model. The average session/run effect should still be specified as:

-gltCode group.mean 'Sess : 0.5*1 +0.5*2'  \

Gang Chen

Thank you very much Gang, it helped a lot!
Sadly I had a result for the glt that was completely empty

what does it means Interpret type III hypotheses with care.
++ Smallest FDR q [0 Sess Chi-sq] = 0.00440489
*+ WARNING: Smallest FDR q [1 run Chi-sq] = 0.999865 ==> few true single voxel detections
++ Smallest FDR q [2 Sess:run Chi-sq] = 1.2213e-11
*+ WARNING: mri_fdrize: will not process only 0 values (min=20)
++ Smallest FDR q [4 groupeffect Z] = 0

Talking with chat gpt, I was suggested:

2️⃣ Type II Hypothesis (Main Effects Without Interaction)

To test main effects without interaction bias, remove Sess:run from the model:

3dLMEr -model "Sess+run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)"

Why?

    This forces Sess and run to be tested independently, avoiding distortions from unbalanced interactions.
    You can still examine interactions separately by running a different model just for interactions:

    3dLMEr -model "Sess*run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)" -gltCode "Sess:run" '1'

3️⃣ Testing the Group Effect

To test if the group mean is different from zero, modify gltCode:

-gltCode group.mean 'Sess : 1*01 + 1*02'

OR, if you want to correct for the imbalance in Sess:

-gltCode group.mean 'Sess : X*01 + Y*02'

where X and Y are weights adjusted to reflect the number of observations.

I indeed, tested 3dLMEr with Sess+run instate of Sess*run and it did gave me the expected group results.

Does it make sense? Do the other proposals of chat make sense to you ?

Thank you again,
Clément

Clément,

The chatbot's explanations/suggestions seem off target.

Could you describe your data structure in more detail? How are sessions and runs coded in the data table? Does each session have two runs, or are runs nested within sessions?

Could you also provide a contingency table for sessions and runs in terms of the number of participants as below?

contingency table
           session1 session2
run1          ?        ?
run2          ?        ?

Gang chen

Dear Gang,
thank you again

here is my datatable

Subj Sess Run InputFile
147BCBB 01 0 sub-147BCBB_ses-01_task-rest_bold_correlations_fish.nii.gz
169BAB 01 0 sub-169BAB_ses-01_task-rest_bold_correlations_fish_resampled.nii.gz
184CB 01 0 sub-184CB_ses-01_task-rest_bold_correlations_fish_resampled.nii.gz
208CBF 01 0 sub-208CBF_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
208CBF 01 1 sub-208CBF_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
208CBF 01 2 sub-208CBF_ses-01_task-rest_run-03_bold_correlations_fish_resampled.nii.gz
263BCE 02 0 sub-263BCE_ses-02_task-rest_bold_correlations_fish_resampled.nii.gz
263BCE 01 0 sub-263BCE_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
263BCE 01 1 sub-263BCE_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
276BC 01 0 sub-276BC_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
276BC 01 1 sub-276BC_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
283CA 01 0 sub-283CA_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
283CA 01 1 sub-283CA_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
283CCA 02 0 sub-283CCA_ses-02_task-rest_bold_correlations_fish_resampled.nii.gz
283CCA 01 0 sub-283CCA_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
283CCA 01 1 sub-283CCA_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
283EA 02 0 sub-283EA_ses-02_task-rest_bold_correlations_fish_resampled.nii.gz
283EA 01 0 sub-283EA_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
283EA 01 1 sub-283EA_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
285AAA 01 0 sub-285AAA_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
285AAA 01 1 sub-285AAA_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
285AB 01 0 sub-285AB_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
285AB 01 1 sub-285AB_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
285D 02 0 sub-285D_ses-02_task-rest_bold_correlations_fish_resampled.nii.gz
285D 01 0 sub-285D_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
285D 01 1 sub-285D_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz
285D 01 2 sub-285D_ses-01_task-rest_run-03_bold_correlations_fish_resampled.nii.gz
285D 01 3 sub-285D_ses-01_task-rest_run-04_bold_correlations_fish_resampled.nii.gz
285E 01 0 sub-285E_ses-01_task-rest_run-01_bold_correlations_fish_resampled.nii.gz
285E 01 1 sub-285E_ses-01_task-rest_run-02_bold_correlations_fish_resampled.nii.gz

Clément,

Here is the contingency table for sessions and runs based on your data table:

...            Run
            0  1  2  3
Session  1 13 10  2  1
         2  4  0  0  0

The modeling issue here is not about the type of hypothesis test (e.g., Type I, II, III), as the chatbot suggested. Instead, the missing cells in the table indicate that there is no way to estimate interaction effects between Session and Run. Given this limitation, the most complex model you can fit is:

-model "Sess+Run+(1|Subj)+(1|Sess:Subj)+(1|Run:Subj)" \

Regarding the group mean, equal weights should be used, contrary to the chatbot's suggestion.

Gang Chen

that make sense,
one last question:
Looking at the output I just want to make that I understand sure who is who

  -- At sub-brick #0 'Sess Chi-sq' datum type is float:            0 to        29.398
     statcode = fict;  statpar = 2
  -- At sub-brick #1 'run Chi-sq' datum type is float:            0 to        15.154
     statcode = fict;  statpar = 2
  -- At sub-brick #2 'Sess:run Chi-sq' datum type is float:            0 to        363.37
     statcode = fict;  statpar = 2

If i used Sess+ run, I loose image 3;
which make sense and made me think I could probably simplify?

-model "Sess+run+(1|Subj)+(1|Sess:Subj)"

anyway,

I would have guess based on what you said previously the meaning of the 3 images is:

  • sub-brick #0 Cross-individual variance: (1|Subj)
  • sub-brick #1 Cross-session-individual variance: (1|Sess:Subj)
  • sub-brick #2 Cross-run-individual variance: (1|run:Subj)
  • residual image
    Residual variance (cross-session-run-individual)

I am just confused because image 1 says:
Sess and not Subj

Sorry for all my questions,
Clément

...
  -- At sub-brick #0 'Sess Chi-sq' datum type is float:            0 to        29.398
     statcode = fict;  statpar = 2
  -- At sub-brick #1 'run Chi-sq' datum type is float:            0 to        15.154
     statcode = fict;  statpar = 2
  -- At sub-brick #2 'Sess:run Chi-sq' datum type is float:            0 to        363.37
     statcode = fict;  statpar = 2

Does this output above come from the full model below, without any missing combinations between Session and run?

-model "Sess*run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)"

The first three sub-bricks represent population-level effects: the main effect of Session, the main effect of run, and their interaction Sess:run.

As mentioned earlier, the full model accounts for four variance components:

  1. Inter-individual variance (1|Subj)
  2. Inter-session variance within individuals (1|Sess:Subj)
  3. Inter-run variance within individuals (1|run:Subj)
  4. Residual variance (within-session, within-run, within-individual)

As I mentioned earlier, 3dLMEr does not provide an option to extract these variance components. Therefore, do not conflate the first three sub-bricks of population-level effects with these variance components.

Gang Chen

Does this output above come from the full model below, without any missing combinations between `Session` and `run`?

yes, it is

I understand now, thank you for the explanations,
If I am looking for the main effect of subject, should I add something like

if their is 10 subject for example:

-gltCode group.mean 'Subj : 0.1*147BCBB +0.1*169BAB ....'  \

or should I create another column with labels 1 for example for a subject and 0 for the other?
then do that for all subjects and average it?

Clément

If you want to make inferences about a specific set of 10 participants, create a separate 3dLMEr script that includes only those individuals.

Gang Chen

I see, but otherwise with this script is their way to get the subject effect in addition to session and run?

or is it the interaction interaction Sess:run?

If yes, it would means that with the model
-model "Sess+run+(1|Subj)+(1|Sess:Subj)+(1|run:Subj)"

Subject effect is not availalble?

Subject effect is not availalble?

Could you clarify what you mean by "subject effect"? A precise definition would help ensure we're on the same page.

Gang Chen

of course,
I am referring to a difference in signal that can be attributed to subject differences.
As my understanding of Sess effect and run effect are difference in signal that can be attributed to Sess or run differences (e.g difference between Sess 1 and 2, which in this case is specific to the difference between this two sessions). My aim is to extract difference in signal that can be attributed to subject differences at the population level.

The intra-class correlation (ICC) measures the proportion of total variability that can be attributed to differences between individuals. If that aligns with your goal, you might try using 3dICC with the same model specification as your 3dLMEr setup.

Gang Chen