Lineares Modell: Vergleich der Vorhersagekraft zweier verschiedener Messmethoden

Ich interessiere mich für Vorhersagen Yund studiere verschiedene zwei Messtechniken X1und X2. Es könnte zum Beispiel sein, dass ich den Geschmack einer Banane vorhersagen möchte, indem ich entweder messe, wie lange sie auf dem Tisch liegt, oder indem ich die Anzahl der braunen Flecken auf der Banane messe.

Ich möchte wissen, welche der Messtechniken besser ist, sollte ich nur eine durchführen.

Ich kann ein lineares Modell in R erstellen:

m1 = lm(Y ~ X1)
m2 = lm(Y ~ X2)

X1Nehmen wir an, es ist ein überlegener Prädiktor für den Bananengeschmack als X2. Bei der Berechnung des der beiden Modelle ist das des Modells deutlich höher als das Modell$R^2$$R^2$m1m2 . Bevor ich eine Arbeit darüber schreibe, wie Methode X1besser ist als X2, möchte ich einen Hinweis darauf geben, dass der Unterschied nicht zufällig ist, möglicherweise in Form eines p-Werts.

Wie würde man das machen? Wie mache ich das, wenn ich verschiedene Bananenmarken verwende und zu einem Modell mit linearem Mischeffekt wechsle, bei dem die Bananenmarke als zufälliger Effekt verwendet wird?

Später

Eine Sache möchte ich hinzufügen, nachdem ich gehört habe, dass Sie lineare Modelle mit gemischten Effekten haben: Die Modelle und B I C können weiterhin zum Vergleichen der Modelle verwendet werden. Siehe dieses Papier zum Beispiel. Aus anderen ähnlichen Fragen auf der Website geht hervor, dass dieses Papier von entscheidender Bedeutung ist.$AIC, AIC_{c}$$BIC$

Ursprüngliche Antwort

$AIC$$AIC_{c}$$BIC$$AIC, AIC_{c}, BIC$ bei der Berichterstattung über die Ergebnisse.

Aufgrund dieser und dieser Antworten empfehle ich folgende Ansätze:

1. Erstellen Sie eine Streudiagramm-Matrix (SPLOM) Ihres Datensatzes, einschließlich Glättungselementen : pairs(Y~X1+X2, panel = panel.smooth, lwd = 2, cex = 1.5, col = "steelblue", pch=16). Überprüfen Sie, ob die Linien (die Glätter) mit einer linearen Beziehung kompatibel sind. Verfeinern Sie das Modell bei Bedarf.
2. Berechnen Sie die Modelle m1und m2. Führen Sie einige Modellprüfungen (Residuen usw.) durch: plot(m1)und plot(m2).
3. $AIC_{c}$$AIC$$AIC_{c}$R psclAICcabs(AICc(m1)-AICc(m2))$AIC_{c}$
4. Berechnen Sie Likelihood-Ratio-Tests für nicht verschachtelte Modelle. Das R Paketlmtest enthält die Funktionen coxtest(Cox-Test), jtest(Davidson-MacKinnon J-Test) und encomptest(einschließlich Test von Davidson & MacKinnon).

Einige Gedanken: Wenn die beiden Bananenmaße wirklich dasselbe messen, sind beide möglicherweise gleichermaßen für die Vorhersage geeignet , und es gibt möglicherweise kein "bestes" Modell.

Dieses Papier könnte auch hilfreich sein.

Hier ist ein Beispiel in R:

#==============================================================================
# Generate correlated variables
#==============================================================================

set.seed(123)

R <- matrix(cbind(
  1   , 0.8 , 0.2,
  0.8 , 1   , 0.4,
  0.2 , 0.4 , 1),nrow=3) # correlation matrix
U <- t(chol(R))
nvars <- dim(U)[1]
numobs <- 500
set.seed(1)
random.normal <- matrix(rnorm(nvars*numobs,0,1), nrow=nvars, ncol=numobs);
X <- U %*% random.normal
newX <- t(X)
raw <- as.data.frame(newX)
names(raw) <- c("response","predictor1","predictor2")

#==============================================================================
# Check the graphic
#==============================================================================

par(bg="white", cex=1.2)
pairs(response~predictor1+predictor2, data=raw, panel = panel.smooth,
      lwd = 2, cex = 1.5, col = "steelblue", pch=16, las=1)

Die Glätter bestätigen die linearen Beziehungen. Das war natürlich beabsichtigt.

#==============================================================================
# Calculate the regression models and AICcs
#==============================================================================

library(pscl)

m1 <- lm(response~predictor1, data=raw)
m2 <- lm(response~predictor2, data=raw)

summary(m1)

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.004332   0.027292  -0.159    0.874    
predictor1   0.820150   0.026677  30.743   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.6102 on 498 degrees of freedom
Multiple R-squared:  0.6549,    Adjusted R-squared:  0.6542 
F-statistic: 945.2 on 1 and 498 DF,  p-value: < 2.2e-16

summary(m2)

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.01650    0.04567  -0.361    0.718    
predictor2   0.18282    0.04406   4.150 3.91e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.021 on 498 degrees of freedom
Multiple R-squared:  0.03342,   Adjusted R-squared:  0.03148 
F-statistic: 17.22 on 1 and 498 DF,  p-value: 3.913e-05

AICc(m1)
[1] 928.9961

AICc(m2)
[1] 1443.994

abs(AICc(m1)-AICc(m2))
[1] 514.9977

#==============================================================================
# Calculate the Cox test and Davidson-MacKinnon J test
#==============================================================================

library(lmtest)

coxtest(m1, m2)

Cox test

Model 1: response ~ predictor1
Model 2: response ~ predictor2
                Estimate Std. Error   z value  Pr(>|z|)    
fitted(M1) ~ M2   17.102     4.1890    4.0826 4.454e-05 ***
fitted(M2) ~ M1 -264.753     1.4368 -184.2652 < 2.2e-16 ***

jtest(m1, m2)

J test

Model 1: response ~ predictor1
Model 2: response ~ predictor2
                Estimate Std. Error t value  Pr(>|t|)    
M1 + fitted(M2)  -0.8298   0.151702  -5.470 7.143e-08 ***
M2 + fitted(M1)   1.0723   0.034271  31.288 < 2.2e-16 ***

$AIC_{c}$m1$R^{2}$

$R^2, AIC$$BIC$$R^2$

$R^2$

Das Bananenbeispiel ist hier vermutlich scherzhaft, aber ich würde nicht erwarten, dass geradlinige Passungen überhaupt gut funktionieren ...

Die Inferenzmaschinerie, die andere als Antwort herausfahren, ist eine Sache von intellektueller Schönheit, aber manchmal braucht man keinen hochmodernen Vorschlaghammer, um eine Nuss zu knacken. Manchmal scheint es, dass jeder, der diese Nacht veröffentlicht, dunkler als der Tag ist, immer jemanden fragt: "Haben Sie das formal getestet? Was ist Ihr P-Wert?".

Führen Sie einen Cox-Test für nicht verschachtelte Modelle durch.

y <- rnorm( 10 )
x1 <- y + rnorm( 10 ) / 2
x2 <- y + rnorm( 10 )

lm1 <- lm( y ~ x1 )
lm2 <- lm( y ~ x2 )

library( lmtest )

coxtest( lm1, lm2 )
?coxtest

(Sie finden Verweise auf andere Tests).

Siehe auch diesen Kommentar und diese Frage . Erwägen Sie insbesondere die Verwendung von AIC / BIC.