Ist es möglich, ein lmer-Modell ohne feste Effekte anzugeben?

9

Wie spezifiziere ich in R ein früheres Modell ohne globalen festen Effekt? Zum Beispiel, wenn ich so etwas sage

lmer(y ~ (1 | group) + (0 + x | group), data = my_df)

das angepasste Modell wird sein

yij=a+αi+βixij

Wie passe ich Modell

yij=αi+βixij ?

Löwe
quelle
Meine knifflige Antwort war, dass lmer(y~0+(1|group)+(0+x|group))das funktionieren würde, aber dies ergibt einen Fehler.
Mike Lawrence
4
Ich denke nicht, dass es möglich ist, ein Modell ohne festen Effekt mit lmer anzugeben, da das lme4-Paket nur gemischten Modellen gewidmet ist (mit mindestens einem festen Effekt und einem zufälligen Effekt). Ich kann mich jedenfalls nicht erinnern, Modelle mit zufälligen Effekten in der Dokumentation gesehen zu haben. Es könnte nützlich sein, in der Liste der R-Sig-Mixed-Models ( stat.ethz.ch/mailman/listinfo/r-sig-mixed-models )
danach zu
1
Nur aus Neugier, warum sollten Sie das tun wollen? Vielleicht gibt es eine andere Möglichkeit, sich Ihrem zugrunde liegenden Ziel zu nähern.
Jbowman
2
Mitte yzuerst :)
Makro

Antworten:

8

Wie @Mike Lawrence erwähnte, ist das Offensichtliche, was zu tun ist, wenn ein Modell ohne feste Effekte definiert wird, Folgendes:

lmer(y ~ -1 + (1|GroupIndicator))

das ist eigentlich ganz einfach; man definiert keinen Achsenabschnitt oder eine X-Matrix. Der Hauptgrund, warum dies nicht funktioniert, ist, dass, wie @maxTC betonte, "das lme4-Paket nur gemischten Modellen gewidmet ist ".

Insbesondere berechnet lmer () die profilierte Abweichung, indem es die bestrafte Regression der kleinsten Quadrate zwischen und sowie die sphärischen Zufallseffekte und löst (Gleichung (11), Lit. (2)). Dieses Optimierungsverfahren berechnet rechnerisch die Cholesky-Zerlegung des entsprechenden Systems unter Ausnutzung der Blockstruktur des Systems (Gleichung (5), Lit. (1)). Wenn Sie keine globalen festen Effekte festlegen, wird diese Blockstruktur praktisch so verzerrt, dass der Code von lmer () nicht damit umgehen kann. Unter anderem basiert der bedingte Erwartungswert von auf , löst aber nach yu0u & bgr; & bgr; RXXLXy^yu0uβ^β^fragt nach der Lösung eines nie existierenden Matrixsystems (die Matrix in Lit. (1) oder in Lit. (2)). Sie erhalten also eine Fehlermeldung wie:RXXLX

Error in mer_finalize(ans) : 
  Cholmod error 'invalid xtype' at file:../Cholesky/cholmod_solve.c, line 970

Denn schließlich gab es überhaupt nichts zu lösen.

Angenommen, Sie möchten die lmer () - Profilabweichungskostenfunktion nicht neu schreiben, basiert die einfachste Lösung auf dem CS-101-Axiom: Müll rein , Müll raus .

 N = length(y); Garbage <- rnorm(N); 
 lmer(y ~ -1 + Garbage + (1|GroupIndicator));

Wir definieren also einen variablen , der nur Rauschen ist. Nach wie vor wird lmer () angewiesen, keinen festen Achsenabschnitt zu verwenden, sondern nur die von uns definierte X-Matrix (in diesem Fall die einspaltige Matrix Garbage). Diese zusätzliche Gaußsche Rauschvariable wird erwartungsgemäß nicht mit unseren Probenmessfehlern sowie mit Ihrer zufälligen Effektvarianz korreliert. Es ist unnötig zu erwähnen, dass die Wahrscheinlichkeit unerwünschter, aber statistisch signifikanter zufälliger Korrelationen umso geringer ist, je strukturierter Ihr Modell ist.Garbage

Lmer () hat also eine Placebo- Variable (Matrix) zum Spielen. Sie erwarten, dass das zugehörige Null ist, und Sie mussten Ihre Daten in keiner Weise normalisieren (zentrieren, aufhellen usw.). . Wahrscheinlich wird es auch nicht schaden , eine zufällige Initialisierung der Placebo- Matrix zu versuchen . Ein letzter Hinweis zum "Müll": Die Verwendung von Gaußschem Rauschen war nicht "zufällig"; Es hat die größte Entropie unter allen Zufallsvariablen gleicher Varianz, so dass die geringste Chance besteht, einen Informationsgewinn zu erzielen.β X.XβX

Dies ist natürlich eher ein Rechen-Trick als eine Lösung, aber es ermöglicht dem Benutzer , ein früheres Modell ohne globalen festen Effekt effektiv anzugeben. Entschuldigung für die Hoffnung auf die beiden Referenzen. Im Allgemeinen denke ich, dass Ref. (1) die beste Wahl für jeden ist, um zu erkennen, was lmer () tut, aber Ref. (2) ist näher am Geist des tatsächlichen Codes.

Hier ist ein bisschen Code, der die obige Idee zeigt:

library(lme4)
N= 500;                 #Number of Samples
nlevA = 25;             #Number of levels in the random effect 
set.seed(0)             #Set the seed
e = rnorm(N); e = 1*(e - mean(e) )/sd(e); #Some errors

GroupIndicator =  sample(nlevA, N, replace=T)   #Random Nvel Classes 

Q = lmer( rnorm(N) ~ (1| GroupIndicator ));      #Dummy regression to get the matrix Zt easily
Z = t(Q@Zt);            #Z matrix

RA <-  rnorm(nlevA )                        #Random Normal Matrix
gammas =c(3*RA/sd(RA))                      #Colour this a bit

y = as.vector(  Z %*% gammas +  e )         #Our measurements are the sum of measurement error (e) and Group specific variance

lmer_native <- lmer(y ~ -1 +(1| GroupIndicator)) #No luck here.
Garbage <- rnorm(N)                              #Prepare the garbage

lmer_fooled <- lmer(y ~ -1 + Garbage+(1| GroupIndicator)) #OK...
summary(lmer_fooled)                             #Hey, it sort of works!

Verweise:

  1. Lineare gemischte Modelle und bestrafte kleinste Quadrate von DM Bates und S. DebRoy, Journal of Multivariate Analysis, Band 91, Ausgabe 1, Oktober 2004 ( Link zum Preprint )
  2. Berechnungsmethoden für gemischte Modelle von Douglas Bates, Juni 2012. ( Link zur Quelle )
usεr11852
quelle
Gibt es einen Grund, den von Ihnen vorgeschlagenen Ansatz dem von Marco in den Kommentaren erwähnten vorzuziehen?
Russellpierce
2
Ja; Sie ändern die Daten nicht, sodass keine Rücktransformationen erforderlich sind. Als Ergebnis werden alle Standarddiagnosen und Goodies von lmer () z. Angepasste Variablen, Residuen, zufällige Effektstufen usw. usw. sind direkt interpretierbar, da sie Ihrem "wahren" Datensatz entsprechen und nicht einem "geänderten".
usεr11852