Mixed effect model for the spatiotemporal analysis of longitudinal - - PowerPoint PPT Presentation

Stéphanie Allassonnière

with J.B. Schiratti, J. Chevallier, I. Koval, V. Debavalaere and

• S. Durrleman

Université Paris Descartes & Ecole Polytechnique

Mixed effect model for the spatiotemporal analysis of longitudinal manifold valued data

• Represent and analyse geometrical elements upon which deformations

can act

• Describe the observed objects as geometrical variations of one or several

representative elements

• Quantify this variability inside a population

Deformable template model from Grenander

• How does the deformation act?
• What is a representative element?
• How to quantify the geometrical variability ?

2

Computational Anatomy

3

One solution:

• Quantify the distance between observations using deformations
• Provide a statistical model to approximate the generation of the observed

population from the atlas

• Propose a statistical learning algorithm
• Optimise the numerical estimation

Computational Anatomy

4

Bayesian Mixed Effect model

• First model :

– One observation per subject – Image or shape (viewed as currents) – Deformations either linearized or diffeomorphic – Homogeneous or heterogeneous populations (mixture models)

5

Bayesian Mixed Effect model

6

Bayesian Mixed Effect model

7

Bayesian Mixed Effect model

• First model :

– One observation per subject – Image or shape (viewed as currents) – Deformations either linearized or diffeomorphic – Homogeneous or heterogeneous populations (mixture models)

Ø Limitations

Ø One observation per subject Ø Corresponding acquistion time

8

Longitudinal Data Analysis

• Longitudinal model :

– Several observation per subject – Image, shape, etc – Atlas = representative trajectory and population variability

Longitudinal Data Analysis

9

Temporal marker of progression

(e.g. time since drug injection, seeding, birth, etc..)

subject #1 subject #2 subject #3 Time (age) Regression

(e.g. compare measurements at same time-point)

No temporal marker of progression

(e.g. in aging, neurodegenerative diseases, etc..)

How to learn representative trajectories of data changes from longitudinal data?

Learning spatiotemporal distribution of trajectories

• Find temporal correspondences
• Compare data at corresponding

stages of progression

Linear mixed-effects models

[Laird&Ware’82, Diggle et al., Fitzmaurice et al.]

Needs to disentangle differences in:

• Measurements
• Dynamics of measurement

changes vector-valued data manifold-valued data

(normalized data, positive matrices, shapes, etc;;)

Spatiotemporal Statistical Model

10

[Schiratti et al. IPMI’15, NIPS’15]

• Statistical model inclinding:
• a representative trajectory of

data changes

• spatiotemporal variations in:
• measurement values
• pace of measurement

changes

• Orthogonality condition ensures

identifiability (unique space/time decomposition)

• Time is not a covariate but a

random variable

t0 t

v0

p0

T0(t) = Expp0,t0(v0)(t)

vi

pij

vij

yij = Ti(ψi(t)) + εij

Ti(t) = ExpT0(t)(PT0

t0,t(vi))

ψi(t) = t0 + αi(t − t0 − τi)

Acceleration factor Time-shift Space-shift

αi ∼ log N(0, σ2

α)

τi ∼ N(0, σ2

τ)

vi = (A1| . . . |AK) si Ak⊥v0 Random effects: Fixed effects: (p0, t0, v0)

(σ2

α, σ2 τ, A1, ...AK)

and

11

Spatiotemporal Statistical Model

yij = Ti(ψi(t)) + εij

Ti(t) = ExpT0(t)(PT0

t0,t(vi))

T0(t) = Expp0,t0(v0)(t)

ψi(t) = t0 + αi(t − t0 − τi) αi ∼ log N(0, σ2

α)

τi ∼ N(0, σ2

τ)

vi = (A1| . . . |AK) si Ak⊥v0 (p0, t0, v0)

{

Submanifold value observations Parallel curve Representative trajectory Linear time reparametrization

{

Hidden random variables: Acceleration factor Time shift Space shift Parameters: Mean trajectory parametrization and prior parameter

(σ2

α, σ2 τ, A1, ...AK)

{

12

Spatiotemporal Statistical Model

vi

pij

vij

p0

Σ P0tT ΣP0t

Interest: Parallel transport keep invariant the structure of the distribution, but updated it in time Comparison with previous work:

v0

p0

vi vi

pij

vij v0

p0

Spatiotemporal Statistical Model

13

• The straight line model M = R

y time

b

t0 yij = (a + ai)(ti,j − t0) + b + bi | {z } +εi,j

Measurement of the ith subject at time t0 Laird & Ware (1982) Schiratti et al. (2015)

time

b

t0 yij = (a + ai)(ti,j −t0 − τi | {z }) + b + εi,j

Time at which measurement of the ith subject reaches ¯

b

x x

Spatiotemporal Statistical Model

14

• The logistic curve model
• Geodesic are logistic curves
• It is not equivalent to a linear model on the logit of the observations (i.e. the

Riemannian log at p0 = 0.5), since p0 is estimated

• If we fix p0 = 0.5 in our model à end up with our previous linear case (different

from Laird&Ware)

M =]0, 1[, g(p)(u, v) = uv p2(1 − p)2 γ0(t) = 1 + (1 − p0)/p0 exp ⇣ −

v0 p0(1−p0) (t − t0)

⌘ yij = γ0 ⇣ t0 + αi(t − t0 − τi) ⌘ + εij

Spatiotemporal Statistical Model

15

• The propagation model
• Geodesics are logistic curves in each coordinate
• Parametric family of geodesics seen as a model of propagation of an effect
• The parallel curve in the direction of the space-shift vi writes

à The parallel changes the relative timing of the effect onset across coordinates

M =]0, 1[N, g(p)(u, v) =

N

X

k=1

ukvk p2

k(1 − pk)2

γδ(t) = ⇣ γ0(t), γ0(t − δ1), . . . , γ0(t − δN−1) ⌘

✓ γ0 ✓ t + vi,1 v0 ◆ , γ0 ✓ t − δ1 + vi,2 v0 ◆ , ..., γ0 ✓ t − δN−1 + vi,N v0 ◆◆

Parameter Estimation

16

• Maximum Likelihood:
• EM:
• Distribution from the curved exponential family

maxθ p(y|θ) = Z p(y, z|θ)dz p(yi|zi, θ)p(zi|θ) θk+1 = argmaxθ

N

X

i=1

Z log ✓ p(yi, zi|θ) | {z } ◆ p(zi|yi, θk)dzi θk+1 = argmaxθ ( φ(θ)T

N

X

i=1

Z S(yi, zi)p(zi|yi, θk)dzi − N log(C(θ)) )

log p(yi, zi|θ) = φ(θ)T S(yi, zi) − log(C(θ))

y = (y1, ..., yN), z = (z1, ...zN), θ = (σ2

z, σ2 ε, A1, ..., AK, p0, t0, v0)

17

Parameter Estimation: stochastic algorithm

• SA-EM: replaces integration by one simulation of the hidden variable:

sample from , and a stochastic approximation of the sufficient statistics Maximization step (unchanged)

• MCMC-SAEM: replaces sampling by a single Markov Chain step
• For each subject, sample the random effect w.r.t a transition kernel of a

geometrically ergodic Markov chain targeting the conditional distribution

zi,k+1 p(zi|yi, θk) Sk+1 = (1 − ∆k)Sk + ∆k 1 N

N

X

i=1

S(yi, zi,k+1) ! θk+1 = argmaxθ

• φ(θ)T Sk+1 − log(C(θ))

[Delyon, Lavielle, Moulines.’99] [Allassonnière et al.’10]

q(zi|yi, θk)

18

Parameter Estimation: stochastic algorithm

• SA-EM: replaces integration by one simulation of the hidden variable:

sample from , and a stochastic approximation of the sufficient statistics Maximization step (unchanged)

• MCMC-SAEM: replaces sampling by a single Markov Chain step
• For each subject, sample the random effect w.r.t a transition kernel of a

geometrically ergodic Markov chain targeting the conditional distribution As long as “converges towards” as

zi,k+1 p(zi|yi, θk) Sk+1 = (1 − ∆k)Sk + ∆k 1 N

N

X

i=1

S(yi, zi,k+1) ! θk+1 = argmaxθ

• φ(θ)T Sk+1 − log(C(θ))

[Chevallier & Allassonnière, preprint ‘19]

˜ qk(zi, θk)

˜ qk

q(zi|yi, θk)

k → ∞

• Theoretical properties of the sampler:

Under mild conditions: – Drift property – Small set – Geometric ergodicity uniformly on any compact set of the parameters

• Theoretical properties of the estimation algorithm:

– a.s. convergence towards the MAP estimator – Normal asymptotic behaviour: speed – Normal asymptotoc behaviour with optimal speed with averaging sequences

19

Parameter Estimation: stochastic algorithm

1/ p ∆k

1/ √ k

Model of Alzheimer’s disease progression

20

The average trajectory of data changes

• Neuropsychological tests ADAS-

Gog from ADNI

• 248 subjects who converted

from MCI to AD

• 6 time-points per subjects on

average (min 3, max 11)

• Data points

with propagation logistic model yij ∈]0, 1[4

[Schiratti et al. IPMI’15, NIPS’15]

Model of Alzheimer’s disease progression

[Schiratti et al. IPMI’15, NIPS’15]

Acceleration factor 𝛽"

Distinguish fast vs. slow progressers Distinguish early vs. late onset individuals

Time-shift

+1σ −1σ

Model of Alzheimer’s disease progression

[Schiratti et al. IPMI’15, NIPS’15] Variability in the relative timing and ordering of the events

Decomposition vector A1

+1σ −1σ

Decomposition vector A2

23

Computational comparisons

• Comparison of: MCMC-SAEM – STAN – MONOLIX
• Number of iterations:
• MCMC-SAEM: 1 000 000 (6s / 1 000 iterations)
• STAN: 15 000 (25min / 1 000 iterations)
• MONOLIX: 20 000 (3,5 min / 1 000 iterations)

24

Computational comparisons

• Comparison of: MCMC-SAEM – STAN – MONOLIX

True values MCMC-SAEM STAN Monolix

Model of Alzheimer’s disease progression

[Schiratti et al. IPMI’15, NIPS’15]

26

Comparison AD vs Controls

27

Comparison MCI vs Controls

28

Each snapshot corresponds to the neuronal loss within 5 years (from 70-75 to 85-90) Neuronal loss within 10 years

ADNI Data à Alzheimer’s Disease cohort Measures of the cortical thickness for MCI converters

Model of propagation on a graph

[Koval et al, Frontiers in Neurosciences’17]

Spatiotemporal Statistical Model

29

[Schiratti et al. IPMI’15, NIPS’15]

• Geodesic representative trajectory
• > What for treated pathologies?

t0 t

v0

p0

T0(t) = Expp0,t0(v0)(t)

vi

pij

vij

yij = Ti(ψi(t)) + εij

Ti(t) = ExpT0(t)(PT0

t0,t(vi))

ψi(t) = t0 + αi(t − t0 − τi)

Acceleration factor Time-shift Space-shift

αi ∼ log N(0, σ2

α)

τi ∼ N(0, σ2

τ)

vi = (A1| . . . |AK) si Ak⊥v0 Random effects: Fixed effects: (p0, t0, v0)

(σ2

α, σ2 τ, A1, ...AK)

and

30

Spatiotemporal Statistical Model Piecewise geodesic trajectories

Spatiotemporal Statistical Model Piecewise geodesic trajectories

[Chevallier et al, NIPS’17]

Spatiotemporal Statistical Model Piecewise geodesic trajectories

[Chevallier et al, NIPS’17]

Spatiotemporal Statistical Model Piecewise geodesic trajectories

[Chevallier et al, NIPS’17]

Spatiotemporal Statistical Model Piecewise geodesic trajectories

35

The individual parameters are related to the real age of conversion of the individuals

Prediction tool

• Generic statistical model to learn spatiotemporal distribution of

trajectories on manifolds: – Calibrated on longitudinal data sets using MCMC-SAEM – Automatically finds temporal correspondences among similar events that may happen at different age/time – Estimates the variability of the data at the corresponding events

• It allows us to position disease progression within the life and history of

the patient

• Future work:

– Derive instances of the model for more complex manifold-valued data (e.g. textured shapes data), metamorphosis and mixtures of all of these!

36

Conclusion

37

Thank you!

vi

pij

vij

t0 t

v0

p0

T0(t) = Expp0,t0(v0)(t)