wake-up-neo.net

Regressionskoeffizientenwerte extrahieren

Ich habe ein Regressionsmodell für einige Zeitreihendaten, die die Drogennutzung untersuchen. Der Zweck ist, einen Spline an eine Zeitserie anzupassen und 95% CI usw. auszuarbeiten. Das Modell lautet wie folgt:

id <- ts(1:length(drug$Date))
a1 <- ts(drug$Rate)
a2 <- lag(a1-1)
tg <- ts.union(a1,id,a2)
mg <-lm (a1~a2+bs(id,df=df1),data=tg) 

Die zusammenfassende Ausgabe von mg lautet:

Call:
lm(formula = a1 ~ a2 + bs(id, df = df1), data = tg)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.31617 -0.11711 -0.02897  0.12330  0.40442 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)        0.77443    0.09011   8.594 1.10e-11 ***
a2                 0.13270    0.13593   0.976  0.33329    
bs(id, df = df1)1 -0.16349    0.23431  -0.698  0.48832    
bs(id, df = df1)2  0.63013    0.19362   3.254  0.00196 ** 
bs(id, df = df1)3  0.33859    0.14399   2.351  0.02238 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Ich verwende den Pr(>|t|)-Wert von a2, um zu testen, ob die untersuchten Daten automatisch korreliert werden.

Kann man diesen Wert von Pr(>|t|) (in diesem Modell 0.33329) extrahieren und in einem Skalar speichern, um einen logischen Test durchzuführen?

Kann alternativ mit einer anderen Methode gearbeitet werden?

55
John

Ein summary.lm-Objekt speichert diese Werte in einer matrix mit dem Namen 'coefficients'. Auf den Wert, nach dem Sie suchen, können Sie also mit: 

a2Pval <- summary(mg)$coefficients[2, 4]

Oder allgemeiner/lesbarer coef(summary(mg))["a2","Pr(>|t|)"]. hier , warum diese Methode bevorzugt wird.

61
wkmor1

Das Paket broom ist hier praktisch (es verwendet das "Tidy" -Format).

tidy(mg) liefert ein schön formatiertes data.frame mit Koeffizienten, t-Statistiken usw. Funktioniert auch für andere Modelle (z. B. plm, ...).

Beispiel aus dem Github-Repo von broom:

lmfit <- lm(mpg ~ wt, mtcars)
require(broom)    
tidy(lmfit)

      term estimate std.error statistic   p.value
1 (Intercept)   37.285   1.8776    19.858 8.242e-19
2          wt   -5.344   0.5591    -9.559 1.294e-10

is.data.frame(tidy(lmfit))
[1] TRUE
24
Helix123

Übergeben Sie einfach Ihr Regressionsmodell in die folgende Funktion:

    plot_coeffs <- function(mlr_model) {
      coeffs <- coefficients(mlr_model)
      mp <- barplot(coeffs, col="#3F97D0", xaxt='n', main="Regression Coefficients")
      lablist <- names(coeffs)
      text(mp, par("usr")[3], labels = lablist, srt = 45, adj = c(1.1,1.1), xpd = TRUE, cex=0.6)
    }

Verwenden Sie wie folgt:

model <- lm(Petal.Width ~ ., data = iris)

plot_coeffs(model)

 enter image description here

1
Cybernetic

Um Ihre Frage zu beantworten, können Sie den Inhalt der Modellausgabe untersuchen, indem Sie das Modell als Variable speichern und im Umgebungsfenster darauf klicken. Sie können dann klicken, um zu sehen, was darin enthalten ist und was wo gespeichert ist.

Eine andere Möglichkeit ist, yourmodelname$ einzugeben und die Komponenten des Modells nacheinander auszuwählen, um zu sehen, was die einzelnen Komponenten enthalten. Wenn Sie zu yourmodelname$coefficients gelangen, werden alle von Ihnen gewünschten beta-, p- und t-Werte angezeigt.

0