Roys unionintersection principle, random fields, and brain imaging

1
Roys union-intersection principle, random
fields, and brain imaging
ISI Platinum Jubilee, Jan 1-4, 2008

• Keith Worsley,
• McGill
• Jonathan Taylor,
• Stanford and Université de Montréal

2
Sankhya (1937) 365-72
3
Deformation Based Morphometry (DBM) (Tomaiuolo
et al., 2004)

• n1 19 non-missile brain trauma patients, 3-14
  days in coma,
• n2 17 age and gender matched controls
• Data non-linear vector deformations needed to
  warp each MRI to an atlas standard
• Y1i(s) deformation vector of trauma patient i
  at point s
• Y2i(s) deformation vector of control patient i
  at point s
• Locate damage find regions where
• deformations are different
• Test statistic Hotellings T2 (1931)
• T(s)2 (Y1.(s) Y2.(s))S-1 (Y1.(s) Y2.(s))
• P.C. Mahalanobis with Harold Hotelling,
• at Gupta Nivas, 1939

4
Multivariate linear models for random field data

regressors p1
pgt1 q1
T F variables
qgt1 Hotellings T2 Wilks ?

Pillais trace
Roys max root
etc.
5
EC heuristic
6
(No Transcript)
7
After lots of messy algebra
8
Last case

regressors p1
pgt1 q1
T? F? variables
qgt1 Hotellings T2? Wilks ?

Pillais trace
Roys max root
etc.
?
9
Roys union-intersection principle
Almost the EC density of the Roys maximum root
field
?
Why ½? F(s,u)F(s,-u)
?
10
Non-negative least squares random field theory
11
fMRI data 120 scans, 3 scans hot, rest, warm,
rest,
First scan of fMRI data
1000
Highly significant effect, T6.59
hot
890
rest
880
500
870
warm
0
100
200
300
No significant effect, T-0.74
0
820
hot
T statistic for hot - warm effect
rest
800
5
warm
0
100
200
300
Drift
0
810
800
-5
790
T (hot warm effect) / S.d. t110 if no
effect
0
100
200
300
Time, seconds
12
Linear model regressors
Alternating hot and warm stimuli separated by
rest (9 seconds each).
2
warm
warm
hot
hot
1
0
-1
0
50
100
150
200
250
300
350
Hemodynamic response function (HRF)
0.4
Delays and disperses the stimuli by 6s
0.2
0
-0.2
0
50
Regressors x(t) stimuli HRF, sampled every 3
seconds
2
1
0
-1
0
50
100
150
200
250
300
350
Time, seconds
13
Linear model for fMRI time series with AR(p)
errors
? ?

unknown

parameters ? ?
?
14
Allowing for unknown delay of the HRF
x(t)
2
x2(t)
x1(t)
1
stimulus
0
-1
10
20
30
40
50
60
70
t
Time,
(seconds)
15
Non-negative coefficients is a cone alternative
x1(t)
z2
ß1 0, ß2 0
Cone angle ? angle between x1 and x2
Cone alternative
z1
Null 0
ß1 0 ß2 0
ß1 0, ß2 0
x2(t)
16
Example of three extremes
x3(t) spread 4 seconds
4
3D cone
0 25
4
x1(t) standard
x2(t) delayed 4 seconds
z3
z2
z1
17
General non-negative least squares problem
Footnote Woolrich et al. (2004) replace hard
constraints by soft constraints through a prior
distribution on ß, taking a Bayesian approach.
18
(No Transcript)
19
(No Transcript)
20
Pick a range of say p 150 plausible values of
the non-linear parameter ? ?1,,?p
21
Fitting the NNLS model

• Simple
• Do all subsets regression
• Throw out any model that does not satisfy the
  non-negativity constraints
• Among those left, pick the model with smallest
  error sum of squares
• For larger models there are more efficient
  methods e.g. Lawson Hanson (1974).
• The non-negativity constraints tend to enforce
  sparsity, even if regressors are highly
  correlated (e.g. PET).
• Why? Highly correlated regressors have huge
  positive and negative unconstrained coefficients
  non-negativity suppresses the negative ones.

22
Example n20, p150, but surprisingly it does
not overfit
30

Y
Yhat
25
Yhat component 1
Yhat component 2
20
15
Tracer
10
5
0

0
1000
2000
3000
4000
5000
6000
Time
Tend to get sparse pairs of adjacent regressors,
suggesting best regressor is somewhere inbetween.
23
P-values for testing the cone alternative?
1 when j?
24
P-values for PET data at a single voxel
p 2 Cone weights w1 ½, w2 ?/2p
j
j
25
The Beta-bar random field
1 when j?
From well-known EC densities of F field
From simulations at a single voxel
Same linear combination!
26
Roys union-intersection principle
x1(t)
z2
Cone
u
0
1
z1
U
x2(t)
27
Roys union-intersection principle
?
?
28
Functions of Gaussian fields Example
Z1N(0,1)
Z2N(0,1)
s2

3
2
1
0
-1
-2

-3
s1
Rejection regions,
Excursion sets,
Z2

4
Cone alternative
Search Region, S
2
3
0
2
Z1
Null
1
-2

Threshold t
-2
2
0
29
Euler characteristic heuristic again
Excursion sets, Xt
Search Region, S
EC blobs - holes 1 7
6 5 2 1
1 0
Observed
Expected
Euler characteristic, EC
Threshold, t
EXACT!
30
Beautiful symmetry
Steiner-Weyl Tube Formula (1930)
Taylor Gaussian Tube Formula (2003)

• Put a tube of radius r about the search region
  ?S and rejection region Rt

Z2N(0,1)
Rt
r
Tube(Rt,r)
Tube(?S,r)
r
?S
Z1N(0,1)
t
t-r

• Find volume or probability, expand as a power
  series in r, pull off coefficients

31
Tube(?S,r)
r
?S
Steiner-Weyl Volume of Tubes Formula (1930)
Lipschitz-Killing curvatures are just intrinisic
volumes or Minkowski functionals in the
(Riemannian) metric of the variance of the
derivative of the process
32
S
S
S
Edge length ?
Lipschitz-Killing curvature of triangles
Lipschitz-Killing curvature of union of triangles
33
Non-isotropic data?
ZN(0,1)
s2
s1

• Can we warp the data to isotropy? i.e. multiply
  edge lengths by ??
• Globally no, but locally yes, but we may need
  extra dimensions.
• Nash Embedding Theorem dimensions D
  D(D1)/2 D2 dimensions 5.
• Better idea replace Euclidean distance by the
  variogram
• d(s1, s2)2 Var(Z(s1)
  - Z(s2)).

34
(No Transcript)
35
Non-isotropic data? Use Riemannian metric of
Var(?Z)
ZN(0,1)
ZN(0,1)
s2
s1
Edge length ?
Lipschitz-Killing curvature of triangles
Lipschitz-Killing curvature of union of triangles
36
We need independent identically distributed
random fields e.g. residuals from a linear model

Lipschitz-Killing curvature of triangles
Lipschitz-Killing curvature of union of triangles
Taylor Worsley, JASA (2007)
37
Beautiful symmetry
Steiner-Weyl Tube Formula (1930)
Taylor Gaussian Tube Formula (2003)

• Put a tube of radius r about the search region
  ?S and rejection region Rt

Z2N(0,1)
Rt
r
Tube(Rt,r)
Tube(?S,r)
r
?S
Z1N(0,1)
t
t-r

• Find volume or probability, expand as a power
  series in r, pull off coefficients

38
Z2N(0,1)
Rejection region Rt
Tube(Rt,r)
r
Z1N(0,1)
t
t-r
Taylors Gaussian Tube Formula (2003)
39

From well-known EC densities of F field
From simulations at a single voxel
Same linear combination!
40
Power? S 1000cc brain, FWHM 10mm, P
0.05
Event
Block (20 seconds)
1
1
Cone angle ? 78.4o
Cone angle ? 38.1o
0.9
0.9
T-test on ß1
0.8
0.8
Beta-bar test
0.7
0.7
F-test on (ß1, ß2)
Cone weights w1 ½ w2 ?/2p
0.6
0.6
0.5
0.5
Power of test
Power of test
0.4
0.4
0.5
2
x1(t)
0.3
0.3
0
Response
0
Response
0.2
0.2
x2(t)
-0.5
-2
0.1
0.1
0
20
40
0
20
40
Time t (seconds)
Time t (seconds)
0
0
0
1
2
3
0
1
2
3
d
d
Shift
of HRF (seconds)
Shift
of HRF (seconds)
41
Cross correlation random field
42
Maximum canonical cross correlation random field
43
Mann-Whitney random field