So stellen Sie benutzerdefinierte Kontraste mit lmer in R ein

9

Ich verwende lmer in R, um die Auswirkung von condition ( cond) auf ein Ergebnis zu überprüfen . Hier sind einige aus Daten, wobei s der Objektidentifizierer und a, bund csind Bedingungen.

library("tidyr")
library("dplyr")
set.seed(123)
temp <- data.frame(s = paste0("S", 1:30), 
                   a = rnorm(30, -2, 1), 
                   b = rnorm(30, -3, 1), 
                   c = rnorm(30, -4, 1)) 

Ich würde gerne vergleichen

  1. Ebene azum Mittelwert der Ebenen bund cund
  2. Level bzu Level c.

Meine Frage ist, wie ich die Kontraste so einstelle, dass der Achsenabschnitt den Mittelwert der drei Bedingungen widerspiegelt und die beiden berechneten Schätzungen die Unterschiede gemäß 1. und 2. direkt widerspiegeln.

Ich habe es mit versucht

c1 <- cbind(c(-0.5, 0.25, 0.25), c(0, -0.5, 0.5))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c1))

wo cond2scheint in Ordnung zu sein, ist es aber cond1nicht.

Folgen Wie werden diese benutzerdefinierten Kontraste interpretiert? Ich habe versucht, stattdessen die verallgemeinerte Inverse zu verwenden, aber diese Schätzungen sind auch nicht sinnvoll.

c2 <- t(ginv(c1))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c2))

Ich habe auch Helmert-Kontraste ausprobiert, aber die Mittel stimmen immer noch nicht überein.

gather(temp, cond, result, a, b, c) %>%
  mutate(cond = factor(cond, levels = c("c", "b", "a"))) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = contr.helmert))

Was ist der richtige Weg, um dies zu tun?

M4RT1NK4
quelle
Das klingt nach einem Helmert-Kontrast (c ist die erste Ebene, dann b, dann a).
Michael M
Ich habe auch Helmert ausprobiert, aber die Zahlen sind nicht das Mittel, nach dem ich suche. Ich habe die Frage so bearbeitet, dass sie Helmert-Kontraste enthält, danke.
M4RT1NK4

Antworten:

13

Für die folgenden Schritte benötigen wir den Datenrahmen im Langformat. Der Datenrahmen datenthält die abhängige Variable result, die kategorische Prädiktor cond(Ebene: a, b, und c), und der Zufallsfaktor s.

library(tidyr)
dat <- gather(temp, cond, result, a, b, c)

Im Folgenden werde ich zwei Ansätze veranschaulichen, um eine Kontrastmatrix zu erstellen, die den Bedingungen entspricht, die Sie vergleichen möchten:

  1. ab+c2
  2. bc

Benutzerdefinierte Kontraste

Die Matrix matentspricht den Pegeldifferenzen.

mat <- rbind(c(1, -0.5, -0.5),     # a vs. (b + c) / 2
             c(0, 1, -1))          # b vs. c

Um die tatsächliche Kontrastmatrix zu erstellen, berechnen wir die verallgemeinerte Inverse mit ginv(von MASS).

library(MASS)
cMat <- ginv(mat)
#            [,1]          [,2]
# [1,]  0.6666667 -7.130169e-17
# [2,] -0.3333333  5.000000e-01
# [3,] -0.3333333 -5.000000e-01

Diese Kontrastmatrix cMatkann in verwendet werden lmer.

library(lme4)
res <- lmer(result ~ cond + (1|s), data = dat, 
            contrasts = list(cond = cMat))
coef(summary(res))    
#              Estimate Std. Error    t value
# (Intercept) -2.948115  0.0946025 -31.163182
# cond1        1.351517  0.2006822   6.734612
# cond2        1.153918  0.2317279   4.979625

Wie Sie sehen können, entsprechen die Schätzungen mit festem Effekt den oben angegebenen Unterschieden. Darüber hinaus repräsentiert der Achsenabschnitt den Gesamtmittelwert.

Helmert kontrastiert mit contr.helmert

Sie können auch die integrierte contr.helmertFunktion verwenden, um die Kontrastmatrix zu erstellen.

cHelmert <- contr.helmert(3)
#   [,1] [,2]
# 1   -1   -1
# 2    1   -1
# 3    0    2

Die Reihenfolge entspricht jedoch nicht der in der Frage angegebenen. Daher müssen wir die Reihenfolge der Spalten und Zeilen umkehren. Die erste Spalte entspricht bvs. aund die zweite entspricht cvs. dem Mittelwert von bund a.

cHelmert2 <- cHelmert[c(3:1), 2:1]
#   [,1] [,2]
# 3    2    0
# 2   -1    1
# 1   -1   -1

Vergleichen Sie die Kontrastmatrix cHelmert2mit cMat. Sie werden feststellen, dass die Spalten skalierte Versionen der anderen Matrix sind.

Das Ergebnis von lmerist:

library(lme4)
res2 <- lmer(result ~ cond + (1|s), data = dat, 
             contrasts = list(cond = cHelmert2))
coef(summary(res2))    
#               Estimate Std. Error    t value
# (Intercept) -2.9481150 0.09460250 -31.163182
# cond1        0.4505056 0.06689407   6.734612
# cond2        0.5769590 0.11586393   4.979625

Diese Kontrastmatrix ermöglicht dieselben Vergleiche wie die benutzerdefinierte Kontrastmatrix. Da jedoch die Werte in der Matrix unterschiedlich sind, sind auch die Koeffizienten mit festen Effekten unterschiedlich. Es überrascht nicht, dass die Werte gleich sind.t

Sven Hohenstein
quelle
Vielen Dank! Nur um sicherzugehen, dass ich das jetzt verstehe - wenn ich die erste Ebene mit den restlichen Ebenen in einer 4-Ebenen-Variablen vergleichen matmöchte , wäre das c(1, -1/3, -1/3, -1/3)? Deshalb setze ich die Zahlen immer so, wie sie in der Formel (a + (b + c + d) / 3) stehen, und ginvskaliere sie dann entsprechend, sodass die Koeffizienten die Differenz direkt widerspiegeln. Und als Sie die Reihenfolge im Helmert-Beispiel geändert haben, war das nur passend zur Frage? Ansonsten sollten die Ergebnisse unabhängig von der Reihenfolge der Kontraste gleich sein, oder?
M4RT1NK4
@ M4RT1NK4 Ihre Formel und der entsprechende Kontrast sind korrekt. Die Reihenfolge der Spalten wurde nur geändert, um der Reihenfolge der Spalten in der Frage zu entsprechen. Die Reihenfolge der Zeilen ist jedoch wichtig, da die erste Ebene die Referenzebene ist. In Ihrem Beispiel ist die Referenzstufe die dritte Stufe.
Sven Hohenstein
@SvenHohenstein Ich hatte eine verwandte Frage basierend auf dieser Antwort. stats.stackexchange.com/questions/357781/…
mat