Hello Martin,
Thanks for your confirmation. I downloaded the last version of FreeSurfer: the dev5-20120624, from the "..fswiki/LMEModels" page, but there is the old mris_preproc, so I will use the fsgd file.
Thanks for pointing out about the long.base - files!! That was a mistake of mine. And thanks for the shell command line ! So the fsgd file was created /only/ with long subjects, then I created the lh.thickness_sm10.mgh file and went further ! OK. =======================
Hi Jorge,
Thanks for showing clearly the steps ! Everything went fine. The X contains those columns as you said. our stats command is:
/stats = lme_mass_fit_vw(X,1,Y,144,cortex);
/where "144" is the 2*ones. We have 36 subjects and 2 time points, so we have 72 ones * 2 -> 144. (but I also tried with 72, 36, 2)
Unfortunately, there is again an error:
/File: lme_mass_fit.m Line: 162 Column: 14// //The variable stats in a parfor cannot be classified. Error in ==> lme_mass_fit_vw at 73 [stats1,st1] = lme_mass_fit(X,[],Xrows,Zcols,Y,ni,prs,e); / Also when applying: /[M,Y,144] = sortData(M,1,Y,sID); /there is an error:/An array for multiple LHS assignment cannot contain numeric value.
/Our M contains: <72x2 double> value The X: <72x4 double> value but the sID contains: <73x1 cell>
so we thought that might be the problem and we removed the cell "sID" , and we got and sID: <72x1 cell>. The error persisted.
did we miss something ?
Thank you !
Sincerely, Alex.
Le 06/12/2012 11:15 AM, Martin Reuter a écrit :
Yes,
in 5.2 there should be the two new options --qdec and --qdec-long
to allow passing qdec tables both in the regular cross sectional and in the longitudinal format in addition to the existing other ways of passing subjects lists (--s or --fsgd or --f ).
Now care needs to be taken in the longitudinal setting. What you want is to get the files (e.g. thickness) from the *.long.<base> directories into the target (e.g fsaverage space). Not from the cross sectional directories!
For this you will be able to use the longitudinal qdec file (with a column of 'fsid' and 'fsid-base' to group time points into subjects, see wiki).
In 5.1 you should do the following: create a text file with each longitudinal id (containing the .long.<base> ) per line and pass it using the --f flag. So in the example below you'd write: 002.long.002_base 003.long.003_base 003_m6.long.003_base
You can use a simple shell command to get such a list from a longitudinal qdec file: cat long.qdec.table.txt | awk '{if ($1 != "fsid" && substr($1,0,1) != "#") printf("%s.long.%s\n", $1, $2)}' > outsubjectsfile.txt
Then preproc will do its job. The order of subjects needs to be identical to the qdec file that you use to run the LME matlab tools to ensure thickenss maps are stacked in the same order as you pass the covariates.
Best, Martin
On Thu, 2012-12-06 at 04:27 +0000, jorge luis wrote:
Hi Alex
I think that Martin added the flag --qdec to mris_preproc for version 5.2. Martin, could you confirm this please?
The code for building the study design matrix X that is in the wiki is just an example of how lme tools can be used for that purpose but the actual instance of your code will depend on your particular Qdec table. For instance, usually, Freesurfer's longitudinal Qdec tables are of the form:
fsid fsid-base sID time group... 002 002_base 002 0 0 003 003_base 003 0 1 003_m6 003_base 003 0.56 1 003_m12 003_base 003 1.12 1 003_m24 003_base 003 2.2 1 005 005_base 005 0 0 005_m6 005_base 005 0.45 0 005_m12 005_base 005 1.2 0
were "fsid" is the Freesurfer's ID, "fsid-base" is the name of the subject's specific template, sID is the subject-specific ID uniquely identifying each subject, "time" is the time covariate and group is the group membership covariate (there can be more covariates of course). I suppose that you processed your longitudinal MRI scans using a Qdec table similar to that.
You can read that Qdec table into Matlab:
Qdec = fReadQdec('qdec.table.dat');
Now you need to build your numeric design matrix X from the cell string array Qdec to represent a specific longitudinal design. Here you don't need the columns fsid and fsid-base in Qdec so you simply remove them:
Qdec = rmQdecCol(Qdec,1); Qdec = rmQdecCol(Qdec,1);
The you can grab the subject-specific ID:
sID = Qdec(:,1);
and then remove that column:
Qdec = rmQdecCol(Qdec,1);
You can now convert your cell string array Qdec to a numeric matrix:
M = Qdec2num(Qdec);
Here, the data in M is already ordered according to time for each subject, otherwise, you can order the data using:
[M,Y,ni] = sortData(M,1,Y,sID);
where Y is the cortical thickness data matrix that you read with fs_read_Y. Please take a look at the help of sortData. Note that matrix M only have two columns and you need to build your design matrix from it. In your case it is quite simple:
X = [ones(length(M),1) M M(:,1).*M(:,2)];
That is to say X contains a column of 1s (the intercept), a column with the time covariate, a column with the group membership and a column with the interaction term.
Now, you are ready for fitting the lme model:
stats = lme_mass_fit_vw(X,1,Y,ni,cortex);
The 1 here (second input) means that you are using a single random effect for the intercept term which is located in column 1 of X. Then you can test the interaction term using:
CM.C = [0 0 0 1];
F_stats = lme_mass_F(stats,CM);
and write the significance map for visualization and FDR correction in tksurfer
fs_write_fstats(F_stats,mri,'sig.mgh','sig');
It must be recognized that some basic Matlab knowledge is required to apply the lme tools. But the gain in flexibility for building your design matrix and analysis is significant. As you see, the relatively difficult part is how to get from the Qdec variable to your design matrix X. It requires some study-specific Matlab code. Sorry, but this is what he have at this point.
Best -Jorge
______________________________________________________________ De: Alex Hanganu <al.hanganu@yahoo.ca> Para: FS Mailing List <Freesurfer@nmr.mgh.harvard.edu>; Jorge Luis Bernal-Rusiel <jbernal0019@yahoo.es> Enviado: Miércoles 5 de diciembre de 2012 18:20 Asunto: Re: [Freesurfer] Longitudinal analysis - contrast Can you please help with this ? Thanks in advance ! best, Alex. Le 05/12/2012 6:12 PM, Alex Hanganu a écrit : > Hi Jorge, > > On the page presenting the > "../fswiki/LinearMixedEffectsModels", I think there is a > misunderstanding - > You show the command: > > "mris_preproc --qdec ...." > > but "--qdec" is not recognized, and in the help menu of > "mris_preproc" I don't see this flag. So I guess you meant > the "--fsgd" flag. > > so, in the fsgd file, in our example, with 2 groups, and 2 > time points, would it be more correct to put the classes for > each subject (exmpl.1) or for 2 groups (exmpl.2), or this is > not important at this stage ? > > For the previous error, we found the answer, it was because > in the qdec file, we used the "tab" between the lines. We > changed that to "space", and it was ok. > > Further, after performing the step: > [M,Y,ni] = sortData(M,2,Y,sids); > > we get the Error: > File: sortData.m Line: 59 Column: 6 > Expression or statement is incorrect--possibly unbalanced (, > {, or [. > > Can you please help with this? > > Also, in the "lme_mass_fit_vw.m" help, we saw that we need > to create the X - ordered design matrix. Is this matrix the > qdec file ? > > As we see it now, the steps that we should perform in > Matlab, are these: > > [Y,mri] ... > lhsphere ... > lhcortex ... > Qdec = fRead... > Qdec = rmQdecCol... > sids... > Qdec = rmQdecCol... > M = ... > [M,Y,ni] ... > > lme_mass_fit_vw(M,months,Y,Group,[],OutputFileName,[],[],[]) > > where M - should be something like '-1' '1' '1' '-1' > months - is the column from qdec file, where is shown the > difference between time points > Y - is the lh.thickness_sm10.mgh file > Group - should be the Group column from the qdec file > [] - are showing that the default values are to be taken. > > But this approach doesn't seem to be entirely correct, cause > the M (matrix) includes either the group, either the months > (depending how we play with the command). > > Thank you !! > > Best regards, > Alex. > > > Le 05/12/2012 2:30 PM, jorge luis a écrit : > > > Hi Alex > > > > This error is likely due to the Qdec variable being empty. > > So, nothing was read into this variable when you applied > > > > Qdec = fReadQdec('qdec.table.dat'); > > > > Please check that. If you don't find a solution to this > > then send me your Qdec table data file and I will check it > > out. > > > > Best > > -Jorge > > > > > >