densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · 2014. 6. 11. · e etti...

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


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

Contenuti

1 Un raffronto

2 Il dataset densitometry

3 Analisi step by stepinterpretare gli effettisimulare la risposta


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta


MassimoBorelli

Un raffronto


Analisi step bystep


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 un

esperimento aleatorio• γ effetto di region nel subject• δ effetto del subject• � residuo


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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

)


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

Analisi step by step

step by step

• normalità (no, rimedio)• omoschedasticità• correlazioni (s̀ı, rimedieremo)• interpretazione del modello• simulazione della risposta


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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%


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

Omoschedasticità intensity

> library(car)

> leveneTest(intensity region)

Df F value Pr(>F)

group 2 1.06 0.3493123


67

89

1011


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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.


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

da tenere a mente:selezione del modello

• brainstem• cerebellum• cortex

• brainstemcerebellum• cortex


MassimoBorelli

Un raffronto


Analisi step bystep


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)


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

un modello ad effetti misti

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


MassimoBorelli

Un raffronto


Analisi step bystep


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 + δ + �


67

89

1011


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta


> 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 + �


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta


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


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta

gli effetti casuali: shrinking

6 7 8 9 10

intensitiy N001 cortex


5 N001 7.94 cortex6 N001 9.38 cortex


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


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


MassimoBorelli

Un raffronto


Analisi step bystep


simulare larisposta brainstem cerebellum cortex

67

8910

densitometry

67

8910

simulazione


5 N001 7.94 cortex6 N001 9.38 cortex


MassimoBorelli

Un raffronto


Analisi step bystep


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?

Un raffrontoIl dataset densitometryAnalisi step by stepinterpretare gli effettisimulare la risposta

densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · 2014. 6. 11. · e etti...

Documents

Transcript of densitometry I modelli ad e etti mistiborelli/magnagraecia/29slides.pdf · 2014. 6. 11. · e etti...