Surface-based Exploratory Group Analysis in FreeSurfer

Outline

• Processing Stages • Command-line Stream • Assemble Data • Design/Contrast (GLM Theory) • Analyze • Visualize • Interactive/Automated GUI (QDEC) • Correction for multiple comparisons 2

Aging Exploratory Analysis

In which areas does thickness change with age?

Cortical Thickness vs Aging Salat, et al, 2004, Cerebral Cortex 3

Aging Thickness Study

N=40 (all in fsaverage space)

p<.01

Positive Age Correlation Negative Age Correlation 4

Surface-based Measures

• Morphometric (eg, thickness) • Functional • PET • MEG/EEG • Diffusion (?) sampled just under the surface 5

Processing Stages

• Specify Subjects and Surface measures • Assemble Data: • Resample into Common Space • Smooth • Concatenate into one file • Model and Contrasts (GLM) • Fit Model (Estimate) • Correct for multiple comparisons • Visualize 6

The General Linear Model (GLM)

Is Thickness correlated with Age?

Thickness Dependent Variable, Measurement HRF Amplitude IQ, Height, Weight y1 y2 Subject 1 Subject 2 Age x1 x2 Independent Variable Of course, you’d need more then two subjects … 7

Thickness y1 y2

Linear Model

Intercept: b Slope: m System of Linear Equations y1 = 1 * b + x1 * m y2 = 1 * b + x2 * m Age Matrix Formulation y1 y2 = 1 x1 1 x2 * b m x1 Intercept = Offset x2 X = Design Matrix  = Regression Coefficients = Parameter estimates = “betas” = Intercepts and Slopes = beta.mgh (mri_glmfit) Y = X*   b m mri_glmfit output: beta.mgh

8

Hypotheses and Contrasts

Is Thickness correlated with Age?

Does m = 0?

Null Hypothesis: H0: m=0 Thickness y1 Intercept: b Slope: m  y1 y2 = 1 x1 1 x2 * b m b m m= [0 1]* b m



= C*



C=[0 1]: Contrast Matrix y2 Age mri_glmfit output: gamma.mgh

x1 x2 9

Thickness

More than Two Data Points

Intercept: b Slope: m y1 = 1* b + x1* m y2 = 1* b + x2* m y3 = 1* b + x3* m y4 = 1* b + x4* m Age y1 y2 y3 y4 = 1 x1 1 x2 1 x3 1 x4 * b m Y = X*  • Model Error • Noise • Uncertainty • rvar.mgh

10

t-Test and p-values 

Y = X*  = C* 

t = σ

2 p-value/significance • value between 0 and 1 • closer to 0 means more significant

C



C(X



T β X)

 1

C T

p-value FreeSurfer stores p-values as –log10(p): • 0.1=10 -1  sig=1, 0.01=10 -2  sig=2 • sig.mgh files • Signed by sign of  • p-value is for an unsigned test 0.1 =10 1.0

0.05 =10 1.3

0.01 =10 2.0

0.001=10 3.0

sig 1.0

1.3

2.0

3.0

11

Thickness

Two Groups

Intercept: b1 Slope: m1 Do groups differ in Intercept?

Do groups differ in Slope?

Is average slope different than 0?

… Slope: m2 Age Intercept: b2 12

Thickness

Two Groups

Intercept: b1 Slope: m1 Slope: m2 y11 y12 y21 y22 = 1 0 0 x11 0 1 0 0 x12 1 0 0 x21 1 0 x22 * b1 b2 m1 m2 Age Intercept: b2 y11 = 1* b1 + 0 * b2 + x11* m1 + 0 * m2 y12 = 1* b1 + 0 * b2 + x12* m1 + 0 * m2 y21 = 0 * b1 + 1* b2 + 0 * m1 + x21* m2 y22 = 0 * b1 + 1* b2 + 0 * m1 + x22* m2 Y = X*  13

Two Groups

Do groups differ in Intercept?

Does b1=b2?

Does b1-b2 = 0?

C = [ +1 -1 0 0 ],  = C*  y11 y12 y21 y22 = 1 0 x11 0 1 0 x12 0 0 1 0 x21 0 1 0 x22 * b1 b2 m1 m2 Do groups differ in Slope?

Does m1=m2?

Does m1-m2=0?

C = [ 0 0 +1 -1 ],  = C*  Y = X*   Thickness Intercept: b1 Slope: m1 b1 b2 m1 m2 Is average slope different than 0?

Does (m1+m2)/2 = 0?

C = [ 0 0 0.5

0.5

],  = C*  Slope: m2 Age Intercept: b2 14

Surface-based Group Analysis in FreeSurfer

• Create your own design matrix and contrast matrices • Create an FSGD File • FreeSurfer creates design matrix • You still have to specify contrasts • QDEC • Limited to 2 discrete variables, 2 levels max • Limited to 2 continuous variables 15

Command-line Processing Stages

• Assemble Data (mris_preproc) • Resample into Common Space • Smooth • Concatenate into one file • Fit Model (Estimate) (mri_glmfit) • Correct for multiple comparisons • Visualize (tksurfer) } • Model and Contrasts (GLM) (FSGD) recon-all -qcache 16

Subject ID

Specifying Subjects $SUBJECTS_DIR

bert fred jenny margaret … 17

FreeSurfer Directory Tree

bert Subject ID bem stats morph mri rgb scripts surf tiff label orig T1 brain wm aseg lh.aparc_annnot

rh.aparc_annnot

lh.white

rh.white

lh.thickness

rh.thickness

SUBJECTS_DIR environment variable lh.sphere.reg

rh.sphere.reg

18

Example: Thickness Study

1. $SUBJECTS_DIR/ bert /surf/lh.thickness

2. $SUBJECTS_DIR/ fred /surf/lh.thickness

3. $SUBJECTS_DIR/ jenny /surf/lh.thickness

4. $SUBJECTS_DIR/ margaret /surf/lh.thickness

5. … 19

FreeSurfer Group Descriptor (FSGD) File

• Simple text file • List of all subjects in the study • Accompanying demographics • Like a spreadsheet • Automatic design matrix creation • You must still specify the contrast matrices • Integrated with tksurfer Note: Can specify design matrix explicitly with --design 20

FSGD Format

GroupDescriptorFile 1 Class Male Class Female Variables Age Weight IQ Input bert Male 10 100 1000 Input fred Male 15 150 1500 Input jenny Input margaret Female 20 200 2000 Female 25 250 2500 • One Discrete Factor (Gender) with Two Levels (M&F) • Three Continuous Variables: Age, Weight, IQ Class = Group Note: Can specify design matrix explicitly with --design 21

Female Group

FSGDF



X (Automatic)

Male Age Female Age Male Group

X =

1 0 10 0 100 1 0 15 0 150 0 1 0 20 0 0 1500 0 200 1000 0 0 0 2000 0 1 0 25 Age 0 250 Weight 0 2500 IQ

C

= [-1 1 0 0 0 0 0 0] Tests for the difference in intercept/offset between groups

C

= [ 0 0 -1 1 0 0 0 0] Tests for the difference in age slope between groups DO D S – Different Offset, Different Slope 22

Another FSGD Example

• Two Discrete Factors – Gender: Two Levels (M&F) – Handedness: Two Levels (L&R) • One Continuous Variable: Age GroupDescriptorFile 1 Class MaleRight Class MaleLeft Class FemaleRight Class FemaleLeft Variables Age Input bert MaleLeft 10 Input fred MaleRight 15 Input jenny Input margaret FemaleRight 20 FemaleLeft 25 Class = Group 23

Interaction Contrast

• Two Discrete Factors (no continuous, for now) – Gender: Two Levels (M&F) – Handedness: Two Levels (L&R) L  2 • Four Regressors (Offsets) – MR (  1 ), ML (  2 ), FR (  3 ), FL (  4 ) R  1 GroupDescriptorFile 1 Class MaleRight Class MaleLeft Class FemaleRight Class FemaleLeft Input bert MaleLeft Input fred MaleRight Input jenny Input margaret FemaleRight FemaleLeft M       3  1 ) F  4  2 )  1 +  2 +  3  4 C = [-1 +1 +1 -1]    4  3   24

Number of Regressors

Each Group/Class: • Has its own Intercept • Has its own Slope for each continuous variable • DODS = Different offset, different slope NRegressors = NClasses*(NVariables+1) NRegressors

C

= [-1 1 0 0 0 0 0 0] Tests for the difference in intercept/offset between groups

C

= [ 0 0 -1 1 0 0 0 0] Tests for the difference in age slope between groups 25

Factors, Levels, Groups, Classes

Factors can be Discrete or Continuous: • Continuous Variables: Age, IQ, Volume, etc • Discrete Factors: Gender, Handedness, Diagnosis • Discrete Factors have Levels: • Gender: Male and Female • Handedness: Left and Right • Diagnosis: Normal, MCI, AD Group or Class: Specification of All Discrete Factors: • Left-handed Male MCI • Right-handed Female Normal 26

Assemble Data: mris_preproc

mris_preproc --help --fsgd FSGDFile --hemi lh --meas thickness : Specify subjects thru FSGD File : Process left hemisphere : $SUBJECTS_DIR/subjectid/surf/hemi.thickness

--target fsaverage --o lh.thickness.mgh

: common space is subject fsaverage : output “volume-encoded surface file” Lots of other options!

lh.thickness.mgh – file with thickness maps for all subjects  Input to Smoother or GLM 27

Surface Smoothing

• mri_surf2surf --help • Loads lh.thickness.mgh • 2D surface-based smoothing • Specify FWHM (eg, fwhm = 10 mm) • Saves lh.thickness.sm10.mgh

• Can be slow (~10-60min) • recon-all -qcache 28

mri_glmfit

• Reads in FSGD File and constructs X • Reads in your contrasts (C1, C2, etc) • Loads data (lh.thickness.sm10.mgh) • Fits GLM (ie, computes  ) • Computes contrasts (  =C*  ) • t or F ratios, significances • Significance -log10(p) (.01  2, .001  3) 29

mri_glmfit

mri_glmfit --y lh.thickness.sm10.mgh --fsgd gender_age.txt

--C age .mtx –C gender .mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir

mri_glmfit --help 30

mri_glmfit

mri_glmfit --y lh.thickness.sm10.mgh

--fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir Input file (output from smoothing).

Stack of subjects, one frame per subject 31

mri_glmfit

mri_glmfit --y lh.thickness.sm10.mgh

--fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir •FreeSurfer Group Descriptor File (FSGD) •Group membership •Covariates 32

mri_glmfit

mri_glmfit --y lh.thickness.sm10.mgh --fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir

• Contrast Matrices • Simple text/ASCII files • Test hypotheses 33

mri_glmfit

mri_glmfit --y lh.thickness.sm10.mgh --fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir • Perform analysis on left hemisphere of fsaverage subject • Masks by fsaverage cortex.label

• Computes FWHM in 2D 34

mri_glmfit --y lh.thickness.sm10.mgh --fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir

mri_glmfit

Output directory: lh.gender_age.glmdir

/ beta.mgh – parameter estimates rvar.mgh – residual error variance etc … age / sig.mgh – -log10(p), uncorrected gamma.mgh, F.mgh

gender / sig.mgh – -log10(p) gamma.mgh, F.mgh

35

File->LoadOverlay

Visualization with tksurfer

Threshold: -log10(p), Eg, 2=.01

uncorrected Saturation: -log10(p), Eg, 5=.00001

False Dicovery Rate Eg, .01

View->Configure->Overlay 36

Visualization with tksurfer

File-> Load Group Descriptor File … 37

Problem of Multiple Comparisons

p < 0.10

p < 0.01

p < 10 -7 p value is probability that a voxel is falsely activated • Threshold too liberal: many false positives • Threshold too restrictive: lose activation (false negatives) 38

p<.10

Clusters

p<.01

p<10 -7 • True signal tends to be clustered • False Positives tend to be randomly distributed in space • Cluster – set of spatially contiguous voxels that are above a given threshold.

Cluster-forming Threshold

p<.00001

sig<5 p<.0001

sig<4 Unthresholded p<.001

sig<3 As threshold lowers, clusters may expand or merge and new clusters can form. No way to say what the threshold should be.

p<.0001

sig<4

Cluster Table, Uncorrected

38 clusters ClusterNo Area(mm 2 ) X Y Z Structure Cluster 1 3738.82 -11.1 34.5 27.2 superiorfrontal Cluster 2 5194.19 -32.4 -23.3 15.7 insula Cluster 3 1271.30 -25.9 -75.0 19.0 superiorparietal Cluster 4 775.38 -44.4 -9.7 51.3 precentral Cluster 5 440.56 -33.0 -36.8 37.5 supramarginal … How likely is it to get a cluster of a certain size under the null hypothesis?

Correction for Multiple Comparisons

• Cluster-based – Monte Carlo simulation – Permutation Tests • False Discovery Rate (FDR) – built into tksurfer and QDEC. (Genovese, et al, NI 2002) 42

Cluster-based Corr. for Multiple Comparisons

1. Simulate data under Null Hypothesis: – – Synthesize Gaussian noise and then smooth (Monte Carlo) Permute rows of design matrix (Permutation, orthog) 2. Analyze, threshold, cluster, max cluster size 3. Repeat 10,000 times 4. Analyze real data, get cluster sizes 5. P(cluster) = #MaxClusterSize > ClusterSize/10000 mri_glmfit-sim 43

Cluster Table, Corrected

p<.0001

sig<4 22 clusters out of 38 have cluster p-value < .05

ClusterNo Area(mm 2 ) X Y Z Structure Cluster P Cluster 1 3738.82 -11.1 34.5 27.2 superiorfrontal .0001

Cluster 2 5194.19 -32.4 -23.3 15.7 insula .0003

Cluster 3 1271.30 -25.9 -75.0 19.0 superiorparietal .0050

Cluster 4 775.38 -44.4 -9.7 51.3 precentral .0100

Cluster 5 440.56 -33.0 -36.8 37.5 supramarginal .0400

… Note the difference between the Cluster Forming Threshold (p<.0001) and the Cluster p-value.

Surface-based Corr. for Multiple Comparisons

45 • 2D Cluster-based Correction at p < .05

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces

Surface-based Corr. for Multiple Comparisons

46 • 2D Cluster-based Correction at p < .05

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces Original mri_glmfit command: mri_glmfit --y lh.thickness.sm10.mgh --fsgd gender_age.txt

--C age.mtx –C gender.mtx

--surf fsaverage lh --cortex --glmdir lh.gender_age.glmdir lh.gender_age.glmdir/ beta.mgh – parameter estimates rvar.mgh – residual error variance etc … age/ sig.mgh – -log10(p), uncorrected gamma.mgh, F.mgh

gender/ sig.mgh – -log10(p) gamma.mgh, F.mgh

Surface-based Corr. for Multiple Comparisons

47 • 2D Cluster-based Correction at p < .05

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces • Use pre-cached simulation results • positive contrast • voxelwise threshold = 2 (p<.01) • Can use another simulation or permutation

Surface-based Corr. for Multiple Comparisons

48 • 2D Cluster-based Correction at p < .05

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces Cluster-wise threshold p<.05

Surface-based Corr. for Multiple Comparisons

49 • 2D Cluster-based Correction at p < .05

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces Bonferroni correct over two hemispheres

Correction for Multiple Comparisons Output

50 lh.gender_age.glmdir

mri_glmfit-sim --glmdir lh.gender_age.glmdir

--cache pos 2 --cwpvalthresh .05

--2spaces age gender sig.mgh – pre-existing uncorrected p-values cache.th20.pos.sig.

cache.th20.pos.sig.

cluster.mgh

ocn.annot

– map of significance of clusters – annotation of significant clusters cache.th20.pos.sig.

cluster.summary

– text file of cluster table (clusters, sizes, MNI305 XYZ, and their significances) • Only shows clusters p<.05

Tutorial 1. Command-line Stream

• Create an FSGD File for a thickness study • Age and Gender • Run • mris_preproc • mri_surf2surf • mri_glmfit • mri_glmfit-sim • tksurfer

2. QDEC – same data set

51

QDEC GUI

• Load QDEC Table File • List of Subjects • List of Factors (Discrete and Cont) • Choose Factors • Choose Input (cached): • Hemisphere • Measure (eg, thickness) • Smoothing Level • “Analyze” • Builds Design Matrix • Builds Contrast Matrices • Constructs Human-Readable Questions • Analyzes • Displays Results 52