2016-12-07 60 views
0

我期待從lmer模型中預測'條件',特別是ns樣條。我已經用mtcars數據集複製了這個問題(技術上很差的例子,但是爲了解決這個問題)。R - 從lmer模型中提取ns樣條對象並預測新數據

這裏是我嘗試用線性模型做:

data(mtcars) 
mtcarsmodel <- lm(wt ~ ns(drat,2) + hp + as.factor(gear), data= mtcars) 
summary(mtcarsmodel) 
coef(mtcarsmodel) 
test <- predict(mtcarsmodel, type = "terms") 

完美。但是,對於lmer預測,沒有相應的「條款」選項(unresolved issue here)。

mtcarsmodellmer <- lmer(wt ~ ns(drat,2) + (hp|as.factor(gear)), data= mtcars) 
summary(mtcarsmodellmer) 
coef(mtcarsmodellmer) 
ranef(mtcarsmodellmer) 

由於沒有相當於「預測,術語」功能,我要提取上述固定和隨機係數及係數適用於mtcars數據,但對如何提取的NS花不知道對象來自模型,並將其預測爲一些新數據。對於「多元」轉換變量也是如此。聚(drat,2) - 額外的榮譽,如果你能得到這個。

回答

1

這並不難。

library(lme4) 
library(splines) 
X <- with(mtcars, ns(drat, 2)) ## design matrix for splines (without intercept) 
## head(X) 
#    1   2 
#[1,] 0.5778474 -0.1560021 
#[2,] 0.5778474 -0.1560021 
#[3,] 0.5738625 -0.1792162 
#[4,] 0.2334329 -0.1440232 
#[5,] 0.2808520 -0.1704002 
#[6,] 0.0000000 0.0000000 

## str(X) 
# ns [1:32, 1:2] 0.578 0.578 0.574 0.233 0.281 ... 
# - attr(*, "dimnames")=List of 2 
# ..$ : NULL 
# ..$ : chr [1:2] "1" "2" 
# - attr(*, "degree")= int 3 
# - attr(*, "knots")= Named num 3.7 
# ..- attr(*, "names")= chr "50%" 
# - attr(*, "Boundary.knots")= num [1:2] 2.76 4.93 
# - attr(*, "intercept")= logi FALSE 
# - attr(*, "class")= chr [1:3] "ns" "basis" "matrix" 

fit <- lmer(wt ~ X + (hp|gear), data= mtcars) 

beta <- coef(fit) 
#$gear 
#   hp (Intercept)  X1   X2 
#3 0.010614406 2.455403 -2.167337 -0.9246454 
#4 0.014601363 2.455403 -2.167337 -0.9246454 
#5 0.006342761 2.455403 -2.167337 -0.9246454 
# 
#attr(,"class") 
#[1] "coef.mer" 

如果我們要預測ns來看,只是做

## use `predict.ns`; read `?predict.ns` 
x0 <- seq(1, 5, by = 0.2) ## example `newx` 
Xp <- predict(X, newx = x0) ## prediction matrix 
b <- with(beta$gear, c(X1[1], X2[1])) ## coefficients for spline 
y <- Xp %*% b ## predicted mean 

plot(x0, y, type = "l") 

enter image description here