densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci...

29
I modelli ad effetti misti Massimo Borelli Un raffronto Il dataset densitometry Analisi step by step interpretare gli effetti simulare la risposta I modelli ad effetti misti Massimo Borelli May 30, 2014

Transcript of densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci...

Page 1: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

I modelli ad effetti misti

Massimo Borelli

May 30, 2014

Page 2: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Contenuti

1 Un raffronto

2 Il dataset densitometry

3 Analisi step by stepinterpretare gli effettisimulare la risposta

Page 3: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Page 4: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Page 5: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

i modelli ad effetti misti

• nei modelli (lineari) ad effetti misti ci sono piu componentidi natura aleatoria. Ad esempio:

y = (mx + q) + (γx + δ) + ε

• m, q effetti fissi, le stime dei parametri del modello

• δ, ε calcolati a posteriori come se fossero eventi di unesperimento aleatorio

• γ effetto di region nel subject• δ effetto del subject• ε residuo

Page 6: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

i modelli ad effetti misti

differenze

• effetti fissi: il software fornisce le stime puntuali (e glistandard error)

• effetti casuali: il software fornisce i parametri dellavariabile aleatoria che li genera

• deviazioni standard e correlazioni• analogia: approccio bayesiano

Page 7: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

un dataset non bilanciato

subject measure region

1 N001 91 brainstem2 N001 70 brainstem3 N001 92 cerebellum4 N001 81 cerebellum5 N001 63 cortex6 N001 88 cortex7 N003 77 brainstem.. .. .. ..

125 N130 56 cortex126 N130 33 cortex

subject

N001 6N003 18N004 12N005 6N007 18N010 6N013 6N101 6N104 6N107 6N109 6N111 6N120 12N130 12

Page 8: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Il dataset densitometry

densitometry pixel intensity

brai

n re

gion

cortex

brainstem

cerebellum

40 60 80 100 120

N001 N003

40 60 80 100 120

N004 N005

40 60 80 100 120

N007

cortex

brainstem

cerebellum

N010 N013 N101 N104 N107

cortex

brainstem

cerebellum

N109

40 60 80 100 120

N111 N120

40 60 80 100 120

N130

Page 9: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

quali saranno i risultati

√y = (mx + q) + (γx + δ) + ε

• q = 8.17

• m = {0, 0.66,−0.066} (brainstem, cerebellum, cortex)

• ε ∈ N(0, 0.84)

• δ, γ ∈ N = (0,

1.00 −0.64 −0.79−0.64 1.00 0.97−0.79 0.97 1.00

)

Page 10: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Analisi step by step

step by step

• normalita (no, rimedio)

• omoschedasticita

• correlazioni (sı, rimedieremo)

• interpretazione del modello

• simulazione della risposta

Page 11: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Non c’e normalita

qnorm

measure

40

60

80

100

120

-2 -1 0 1 2

brainstem

-2 -1 0 1 2

cerebellum

-2 -1 0 1 2

cortex

Page 12: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Trasformazione di Box e Cox

> library(MASS)

> boxcox (measure ∼ region)

> boxcox (measure ∼ region, lambda = seq(0, 1, 0.05) )

> intensity = sqrt(measure)

0.0 0.2 0.4 0.6 0.8 1.0

-161.5

-160.5

-159.5

λ

log-Likelihood

95%

Page 13: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Risposta trasformata: intensity

> shapiro.test(intensity[region ==”brainstem”])

W = 0.9183, p-value = 0.005314

> shapiro.test(intensity[region ==”cerebellum”])

W = 0.9601, p-value = 0.1484

> shapiro.test(intensity[region ==”cortex”])

W = 0.9727, p-value = 0.4049

brainstem cerebellum cortex

67

89

1011

Page 14: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Omoschedasticita intensity

> library(car)

> leveneTest(intensity region)

Df F value Pr(>F)

group 2 1.06 0.3493123

brainstem cerebellum cortex

67

89

1011

Page 15: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Tentativo: cercare le correlazioni

> modello = aov(intensity ∼region)

> TukeyHSD(modello)

> plot(TukeyHSD(modello))

cortex-cerebellum

In cortex e cerebellum la intensityappare scorrelata. Brainstempotrebbe, in senso statistico, esserecorrelata a cortex e cerebellum

-1.5 -1.0 -0.5 0.0 0.5 1.0

cortex-cerebellum

cortex-brainstem

cerebellum-brainstem 95% family-wise confidence level

Differences in mean levels of region

Page 16: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

Tentativo: cercare le correlazioni

> tentativo = lm( intensity[region == ”brainstem”] ∼intensity[region ==”cerebellum”] + intensity[region==”cortex”])

Estimate Std. Error t value Pr(>|t|)(Intercept) 0.69 1.40 0.49 0.62

intensitycerebellum 0.60 0.15 4.03 0.00intensitycortex 0.28 0.15 1.83 0.07

cerebellum

cerebellum e un predittore di brainstem per la intensity.

Page 17: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

da tenere a mente:selezione del modello

• brainstem

• cerebellum

• cortex

• brainstemcerebellum

• cortex

Page 18: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

un modello misto massimale

• intensity dipende da region

• fixed effect

• eterogenita dovuta a subject + eterogenita dovuta aregion ’within’ subject

• variabilita residua• random effect

y = q + mx + δ + ε

> library(lme4)

> mixed1 = lmer(intensity ∼ region + (region | subject))

> summary(mixed1)

Page 19: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

un modello ad effetti misti

> mixed1 = lmer(intensity ∼ region + (region | subject))

> summary(mixed1)

Page 20: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

interpretare gli effetti fissi

> fixef(mixed1)

Estimate

(Intercept) 8.17regioncerebellum 0.66

regioncortex -0.06

• brainstem: y = 8.17 + 0 + δ + ε

• cerebellum: y = 8.17 + 0.66 + δ + ε

• cortex: y = 8.17− 0.06 + δ + ε

brainstem cerebellum cortex

67

89

1011

Page 21: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

interpretare gli effetti casuali

> ranef(mixed1)

(Intercept) regioncerebellum regioncortex

N001 0.69 -0.19 -0.27N003 1.43 -0.29 -0.47N004 1.59 -0.91 -0.97N005 1.72 -0.14 -0.40N007 0.93 -0.62 -0.64

.. .. .. ..N120 -1.19 1.63 1.47N130 -1.57 0.11 0.35

Page 22: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

interpretare gli effetti casuali

> ranef(mixed1)

(Intercept) regioncerebellum regioncortex

N001 0.69 -0.19 -0.27.. .. .. ..

Per il subject N001:

• brainstem: y = 8.17 + 0.69 + 0 + ε

• cerebellum: y = 8.17 + 0.69 + 0.66− 0.19 + ε

• cortex: y = 8.17 + 0.69− 0.06− 0.27 + ε

Page 23: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

interpretare gli effetti casuali

Per il subject N001:

• brainstem: y = 8.86 + ε

• cerebellum: y = 9.33 + ε

• cortex: y = 8.53 + ε

subject intensity region

1 N001 9.54 brainstem2 N001 8.37 brainstem3 N001 9.59 cerebellum4 N001 9.00 cerebellum5 N001 7.94 cortex6 N001 9.38 cortex

Page 24: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

gli effetti casuali: shrinking

6 7 8 9 10

intensitiy N001 cortex

subject intensity region

5 N001 7.94 cortex6 N001 9.38 cortex

Page 25: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

simulazione

Groups Name Variance Std.Dev. Corr1 subject (Intercept) 1.72655 1.313982 regioncerebellum 0.62401 0.78994 -0.6353 regioncortex 0.60098 0.77523 -0.790 0.9754 Residual 0.69917 0.83616

La matrice di correlazione: 1.00 −0.64 −0.79−0.64 1.00 0.97−0.79 0.97 1.00

Page 26: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

simulazione

σ =

1.00 −0.64 −0.79−0.64 1.00 0.97−0.79 0.97 1.00

> sigma = matrix(c(1, -0.635, -0.790, -0.635, 1, 0.975,

-0.790, 0.975, 1), 3, 3)

> library(MBESS)

> SIGMA = cor2cov( sigma , sd = c(1.31398, 0.78994,0.77523) )

Σ =

1.73 −0.66 −0.80−0.66 0.62 0.60−0.80 0.60 0.60

Page 27: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

simulazione

> fissi = c( rep(8.17, 42), rep(8.83, 42),rep(8.12, 42) )

> effissi = matrix(fissi, 42, 3)

> pert = mvrnorm(n = 42, rep(0, 3), SIGMA)

> simul = effissi + pert

Page 28: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta brainstem cerebellum cortex

67

8910

densitometry

67

8910

simulazione

subject intensity region

5 N001 7.94 cortex6 N001 9.38 cortex

Page 29: densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · ... ad e etti misti ci sono piu componenti di natura aleatoria. Ad esempio: ... random e ect y = q + mx

I modelli adeffetti misti

MassimoBorelli

Un raffronto

Il datasetdensitometry

Analisi step bystep

interpretare glieffetti

simulare larisposta

parziale conclusione

il modello massimale• calcolato gli effetti fissi

• calcolato gli effetti casuali

cosa manca ancora?• trovare il modello minimale adeguato

• semplificare gli effetti fissi• stime REML o stime ML?

• semplificare gli effetti casuali• parametric bootstrap?