2011-07-05 12 views
45

Mam model regresji dla niektórych szeregów czasowych danych badających wykorzystanie narkotyków. Celem jest, aby dopasować splajn do szeregów czasowych i wypracować 95% CI itp model jest następujący:Wartości współczynników regresji wyodrębnienia

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) 

Wyjście podsumowanie mg jest:

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 

Używam Wartość a2 dla sprawdzenia, czy dane będące przedmiotem badania są powiązane.

Czy można wyodrębnić tę wartość z Pr(>|t|) (w tym modelu 0.33329) i zapisać w skalarnie, aby wykonać test logiczny?

Alternatywnie, czy można to wyjaśnić inną metodą? Obiekt

+0

. @ John - Dlaczego użyłeś 'Pr (> | t |)' wartości 'a2', a nie jednej z pierwszych trzech kolumn? –

Odpowiedz

53

przechowuje te wartości w postaci matrix o nazwie 'coefficients'. Więc wartość jesteś po można uzyskać z:

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

lub, bardziej ogólnie/czytelnie, coef(summary(mg))["a2","Pr(>|t|)"]. Zobacz, dlaczego ta metoda jest preferowana, patrz here.

18

Przydaje się tutaj pakiet broom (używa formatu "schludny").

da ładnie sformatowane dane.frame ze współczynnikami, t statystyki itp. Działa również dla innych modeli (np. Plm, ...).

przykład z broom jest github repo:

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

Aby odebrać OP z tego: '' td [1, "estimate"] '' lub '' td [td $ term == "(Przechwytywanie)", "estimate"] '' lub nawet '' tdt <- as .data.table (td); setkey (tdt); tdt ["(Intercept)", "estimate"] '' – PatrickT

0

prostu przechodzić modelu regresji do następującej funkcji:

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) 
    } 

stosowania w sposób następujący:

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

plot_coeffs(model) 

enter image description here

Powiązane problemy