#### Transcript Detecting Connectivity between Images: the 'Bubbles' Task

Detecting sparse connectivity: the 'bubbles' task in the fMRI scanner Keith Worsley, Math + Stats, Arnaud Charil, Montreal Neurological Institute, McGill Philippe Schyns, Fraser Smith, Psychology, Glasgow Jonathan Taylor, Stanford and Université de Montréal What is ‘bubbles’? Nature (2005) Subject is shown one of 40 faces chosen at random … Happy Sad Fearful Neutral … but face is only revealed through random ‘bubbles’ First trial: “Sad” expression Sad 75 random Smoothed by a bubble centres Gaussian ‘bubble’ What the subject sees 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 Subject is asked the expression: Response: “Neutral” Incorrect Your turn … Trial 2 Subject response: “Fearful” CORRECT Your turn … Trial 3 Subject response: “Happy” INCORRECT (Fearful) Your turn … Trial 4 Subject response: “Happy” CORRECT Your turn … Trial 5 Subject response: “Fearful” CORRECT Your turn … Trial 6 Subject response: “Sad” CORRECT Your turn … Trial 7 Subject response: “Happy” CORRECT Your turn … Trial 8 Subject response: “Neutral” CORRECT Your turn … Trial 9 Subject response: “Happy” CORRECT Your turn … Trial 3000 Subject response: “Happy” INCORRECT (Fearful) Bubbles analysis E.g. Fearful (3000/4=750 trials): 1 Trial 4 + 5 + 2 + 3 + + 6 + 7 + … + 750 1 = Sum 300 0.5 200 0 100 250 200 150 100 50 Correct trials Proportion of correct bubbles =(sum correct bubbles) /(sum all bubbles) 0.75 Thresholded at proportion of 0.7 correct trials=0.68, 0.65 scaled to [0,1] 1 Use this as a 0.5 bubble mask 0 Results Mask average face Happy Sad Fearful Neutral But are these features real or just noise? Need statistics … Statistical analysis Correlate bubbles with response (correct = 1, incorrect = 0), separately for each expression Equivalent to 2-sample Z-statistic for correct vs. incorrect bubbles, e.g. Fearful: Z~N(0,1) Trial 1 2 3 4 5 6 7 … 750 1 0.5 0 1 1 Response 0 1 statistic 4 2 0 -2 0 1 1 … 1 0.75 Very similar to the proportion of correct bubbles: 0.7 0.65 Comparison Both depend on average correct bubbles, rest is ~ constant Z=(Average correct bubbles -average incorrect bubbles) / pooled sd 4 2 0 -2 Proportion correct bubbles = Average correct bubbles / (average all bubbles * 4) 0.75 0.7 0.65 Results Thresholded at Z=1.64 (P=0.05) Happy Average face Sad Fearful Neutral Z~N(0,1) statistic 4.58 4.09 3.6 3.11 2.62 2.13 1.64 Multiple comparisons correction? Need random field theory … Euler Characteristic = #blobs - #holes Excursion set {Z > threshold} for neutral face EC = 0 30 Euler Characteristic 20 0 -7 -11 13 14 9 1 0 Heuristic: At high thresholds t, the holes disappear, EC ~ 1 or 0, E(EC) ~ P(max Z > t). Observed Expected 10 0 -10 -20 -4 -3 -2 -1 0 Threshold 1 2 • Exact expression for E(EC) for all thresholds, • E(EC) ~ P(max Z > t) 3 4accurate. is extremely The details If ¡Z (s)¢ » N(0; 1) is an isot ropic Gaussian random ¯eld, s 2 < 2 , wit h ¸ I 2£ 2 = V @Z , µ ¶ @s P max Z (s) ¸ t ¼ E(E C(S \ f s : Z (s) ¸ tg)) s2 S Z 1 1 = E C(S) e¡ z 2 =2 dz (2¼) 1=2 t Intrinsic volumes or Minkowski functionals ¸ 1=2 + 1 Perimet er(S) e¡ t 2 =2 2 2¼ ¸ + Area(S) te¡ t 2 =2 : (2¼) 3=2 If Z (s) is whit e noise convolved wit h an isot ropic Gaussian ¯lt er of Full Widt h at Half Maximum FWHM t hen ¸ = 4 log 2 . FW HM 2 The brain mapping version If Z (s) is whit e noise smoot hed wit h an isot ropic Gaussian ¯lt er of Full Widt h at Half Maximum FWHM µ ¶ P max Z (s) ¸ t ¼ E(E C(S \ f s : Z (s) ¸ tg)) s2 S Z 1 1 = E C(S) e¡ z 2 =2 dz EC0(S) (2¼) 1=2 Resels0(S) t Resels1(S) Resels2(S) Perimet er(S) (4 log 2) 1=2 e¡ FWHM 2¼ Area(S) 4 log 2 + te¡ t 2 =2 : FWHM 2 (2¼) 3=2 + Resels (Resolution elements) 1 2 t 2 =2 EC1(S) EC2(S) EC densities Math version If T (s) = f (Z 1 (s); : : : ; Z n (s)) is a funct ion³ of i.i.d. Gaussian random ¯elds ´ Z i (s) » Z (s) » N(0; 1), s 2 < D , wit h V @Z ( s) = ¤ D £ D (s), t hen replace @s resels by Lipschitz-K illing curvature L d (S; ¤ ): µ ¶ P max T (s) ¸ t s2 S XD ¼ E(E C(S \ f s : T (s) ¸ tg)) = L d (S; ¤ )½d (t); d= 0 where ½d (t) is t he same EC density for t he isot ropic case wit h ¤ (s) = I D £ D (Taylor & Adler, 2003). Bot h Lipschit z-K illing curvat ure L d (S; ¤ ) and EC density ½d (t) are de¯ned implicit ly as coe± cient s of a power series expansion of t he volume of a t ube as a funct ion of it s radius. In t he case of Lipschit z-K illing curvat ure, t he t ube is about t he search region S in t he Riemannian met ric induced by ¤ (s); in t he case of t he EC densit ies, t he t ube is about t he reject ion region f ¡ 1 ([t; 1 )) and volume is replaced by Gaussian probability. Fort unat ely t here are simple ways of est imat ing Lipschit z-K illing curvat ure from sample dat a (Taylor & Worsley, 2007), and simple ways of calculat ing EC densit ies. The details … 2 Tube(S,r) r S 3 A B 6 Λ is big TubeΛ(S,r) S Λ is small r ν 2 U(s1) S S s1 Tube s2 s3 Tube U(s2) U(s3) Z2 R r Tube(R,r) Z1 N2(0,I) Tube(R,r) R t-r z t z1 RR Tube(R,r) r z2 z3 Summary Random field theory results For searching in D (=2) dimensions, P-value of max Z is (Adler, 1981; W, 1995): P(max Z > z) Resels (=Lipschitz-Killing curvature/c) is Image area / (bubble FWHM)2 = 146.2 (unitless) Euler characteristic density(×c) is ~ E( Euler characteristic of thresholded set ) = Resels × Euler characteristic density (+ boundary terms) (4 log(2))D/2 zD-1 exp(-z2/2) / (2π)(D+1)/2 See forthcoming book Adler, Taylor (2007) Results, corrected for search Thresholded at Z=3.92 (P=0.05) Happy Average face Sad Fearful Neutral Z~N(0,1) statistic 4.58 4.47 4.36 4.25 4.14 4.03 3.92 Bubbles task in fMRI scanner Correlate bubbles with BOLD at every voxel: Trial 1 2 3 4 5 6 7 … 3000 1 0.5 0 fMRI 10000 0 Calculate Z for each pair (bubble pixel, fMRI voxel) – a 5D “image” of Z statistics … Discussion: thresholding Thresholding in advance is vital, since we cannot store all the ~1 billion 5D Z values Resels=(image resels = 146.2) × (fMRI resels = 1057.2) for P=0.05, threshold is Z = 6.22 (approx) The threshold based on Gaussian RFT can be improved using new non-Gaussian RFT based on saddle-point approximations (Chamandy et al., 2006) Model the bubbles as a smoothed Poisson point process The improved thresholds are slightly lower, so more activation is detected Only keep 5D local maxima Z(pixel, voxel) > Z(pixel, 6 neighbours of voxel) > Z(4 neighbours of pixel, voxel) Discussion: modeling The random response is Y=1 (correct) or 0 (incorrect), or Y=fMRI The regressors are Xj=bubble mask at pixel j, j=1 … 240x380=91200 (!) Logistic regression or ordinary regression: logit(E(Y)) or E(Y) = b0+X1b1+…+X91200b91200 But there are only n=3000 observations (trials) … Instead, since regressors are independent, fit them one at a time: logit(E(Y)) or E(Y) = b0+Xjbj However the regressors (bubbles) are random with a simple known distribution, so turn the problem around and condition on Y: E(Xj) = c0+Ycj Equivalent to conditional logistic regression (Cox, 1962) which gives exact inference for b1 conditional on sufficient statistics for b0 Cox also suggested using saddle-point approximations to improve accuracy of inference … Interactions? logit(E(Y)) or E(Y)=b0+X1b1+…+X91200b91200+X1X2b1,2+ … MS lesions and cortical thickness Idea: MS lesions interrupt neuronal signals, causing thinning in down-stream cortex Data: n = 425 mild MS patients 5.5 Average cortical thickness (mm) 5 4.5 4 3.5 3 2.5 Correlation = -0.568, T = -14.20 (423 df) 2 1.5 0 10 20 30 40 50 Total lesion volume (cc) 60 70 80 Charil et al, NeuroImage (2007) MS lesions and cortical thickness at all pairs of points Dominated by total lesions and average cortical thickness, so remove these effects Cortical thickness CT, smoothed 20mm Lesion density LD, smoothed 10mm Find partial correlation(lesion density, cortical thickness) removing total lesion volume linear model: CT-av(CT) ~ 1 + TLV + LD, test for LD Repeat of all voxels in 3D, nodes in 2D Subtract average cortical thickness ~1 billion correlations, so thresholding essential! Look for high negative correlations … Thresholding? Cross correlation random field Correlation between 2 fields at 2 different locations, searched over all pairs of locations one in R (D dimensions), one in S (E dimensions) sample size n Cao & Worsley, Annals of Applied Probability (1999) MS lesion data: P=0.05, c=0.300, T=6.46 Cluster extent rather than peak height (Friston, 1994) fit a quadratic Y to the peak: L D (clust er) » c Â® k Peak height Choose a lower level, e.g. t=3.11 (P=0.001) Find clusters i.e. connected components of excursion set Z D=1 Measure cluster LD extent extent by t Distribution: Distribution of maximum cluster extent: Bonferroni on N = #clusters ~ E(EC). s Cao, Advances in Applied Probability (1999) How do you choose the threshold t for defining clusters? If signal is focal i.e. ~FWHM of noise If signal is broad i.e. >>FWHM of noise Choose a high threshold i.e. peak height is better Choose a low threshold i.e. cluster extent is better Conclusion: cluster extent is better for detecting broad signals Alternative: smooth data with filter that matches signal (Matched Filter Theorem)… try range of filter widths … scale space search … correct using random field theory … a lot of work … Cluster extent is easier! Thresholding? Cross correlation random field Correlation between 2 fields at 2 different locations, searched over all pairs of locations one in R (D dimensions), one in S (E dimensions) µ ¶ XD XE P max Correlat ion > c ¼ L d (R) L e (S) ½d;e (c) R ;S ½d;e (c) = d= 0 e= 0 (d ¡ 1)!e!2n ¡ ¼d +2 e + 1 d¡ e¡ 2 X X (¡ 1) k cd+ e¡ 1¡ 2k (1 ¡ c2 ) n ¡ d ¡ e¡ 1 2 +k k i ;j ¡ ( n ¡ d + i )¡ ( n ¡ e + j ) 2 2 i !j !(k ¡ i ¡ j )!(n ¡ 1 ¡ d ¡ e + i + j + k)!(d ¡ 1 ¡ k ¡ i + j )!(e ¡ k ¡ j + i )! Cao & Worsley, Annals of Applied Probability (1999) MS lesion data: P=0.05, c=0.300, T=6.46 0.1 correlation 0 5 x 10 2.5 2 -0.1 1.5 -0.2 1 -0.3 threshold -0.4 -0.5 0 50 100 150 Different hemisphere 0.1 5 x 10 2.5 0 correlation Histogram Same hemisphere -0.1 2 -0.2 1.5 -0.3 1 0.5 -0.4 0 -0.5 0 threshold 50 100 150 0.5 0 ‘Conditional’ histogram: scaled to same max at each distance 0.1 1 -0.1 0.6 -0.2 0.4 -0.3 -0.4 -0.5 0 threshold 50 100 150 distance (mm) 1 0 0.8 correlation correlation 0 0.1 0.8 -0.1 0.6 -0.2 0.4 -0.3 0.2 -0.4 0 -0.5 0 threshold 50 100 150 distance (mm) 0.2 0 Science (2004) fMRI activation detected by correlation between subjects at the same voxel The average nonselective time course across all activated regions obtained during the first 10 min of the movie for all five subjects. Red line represents the across subject average time course. There is a striking degree of synchronization among different individuals watching the same movie. Voxel-by-voxel intersubject correlation between the source subject (ZO) and the target subject (SN). Correlation maps are shown on unfolded left and right hemispheres (LH and RH, respectively). Color indicates the significance level of the intersubject correlation in each voxel. Black dotted lines denote borders of retinotopic visual areas V1, V2, V3, VP, V3A, V4/V8, and estimated border of auditory cortex (A1).The face-, object-, and building-related borders (red, blue, and green rings, respectively) are also superimposed on the map. Note the substantial extent of intersubject correlations and the extension of the correlations beyond visual and auditory cortices. What are the subjects watching during high activation? Faces … … buildings … … hands Thresholding? Homologous correlation random field Correlation between 2 equally smooth fields at the same location, searched over all locations in R (in D dimensions) Cao & Worsley, Annals of Applied Probability (1999) P-values are larger than for the usual correlation field (correlation between a field and a scalar) E.g. resels=1000, df=100, threshold=5, usual P=0.051, homologous P=0.139