如何将ddply函数用于线性模型?
x1 <- c(1:10, 1:10)
x2 <- c(1:5, 1:5, 1:5, 1:5)
x3 <- c(rep(1,5), rep(2,5), rep(1,5), rep(2,5))
set.seed(123)
y <- rnorm(20, 10, 3)
mydf <- data.frame(x1, x2, x3, y)
require(plyr)
ddply(mydf, mydf$x3, .fun = lm(mydf$y ~ mydf$X1 + mydf$x2)) 这将生成以下错误:
错误(model.frame.default= mydf$y ~ mydf$X1 + mydf$x2,drop.unused.levels = TRUE):变量mydf$X1的无效类型(NULL)
感谢你的帮助。
发布于 2011-09-23 02:02:07
这是你需要做的。
mods = dlply(mydf, .(x3), lm, formula = y ~ x1 + x2)mods是包含回归结果的两个对象的列表。你可以从mods中提取你需要的东西。例如,如果要提取系数,可以编写
coefs = ldply(mods, coef)这给了你
x3 (Intercept) x1 x2
1 1 11.71015 -0.3193146 NA
2 2 21.83969 -1.4677690 NA编辑。如果你想要ANOVA,那你就可以
ldply(mods, anova)
x3 Df Sum Sq Mean Sq F value Pr(>F)
1 1 1 2.039237 2.039237 0.4450663 0.52345980
2 1 8 36.654982 4.581873 NA NA
3 2 1 43.086916 43.086916 4.4273907 0.06849533
4 2 8 77.855187 9.731898 NA NA发布于 2011-09-23 02:31:37
Ramnath解释的是完全正确的。但我会再详细解释一下。
ddply需要一个数据帧,然后返回一个数据帧。lm()函数以数据帧作为输入,但作为回报返回一个线性模型对象。通过?lm查看lm的文档可以看出这一点。
值
lm返回类"lm“的对象或c类的多个响应(”mlm“、"lm")。
因此,您不能只是将lm对象推入数据帧。您的选择是要么将lm的输出强制放入数据帧,要么将lm对象推入列表而不是数据框架。
因此,为了说明这两种选择:
下面是如何将lm对象推入列表(与Ramnath演示的非常类似):
outlist <- dlply(mydf, "x3", function(df) lm(y ~ x1 + x2, data=df))另一方面,如果只想提取系数,则可以创建一个运行回归的函数,然后以如下数据帧的形式返回系数:
myLm <- function( formula, df ){
lmList <- lm(formula, data=df)
lmOut <- data.frame(t(lmList$coefficients))
names(lmOut) <- c("intercept","x1coef","x2coef")
return(lmOut)
}
outDf <- ddply(mydf, "x3", function(df) myLm(y ~ x1 + x2, df))发布于 2011-09-23 02:47:20
使用这个
mods <- dlply(mydf, .(x3), lm, formula = y ~ x1 + x2)
coefs <- llply(mods, coef)
$`1`
(Intercept) x1 x2
11.7101519 -0.3193146 NA
$`2`
(Intercept) x1 x2
21.839687 -1.467769 NA
anovas <- llply(mods, anova)
$`1`
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
x1 1 2.039 2.0392 0.4451 0.5235
Residuals 8 36.655 4.5819
$`2`
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
x1 1 43.087 43.087 4.4274 0.0685 .
Residuals 8 77.855 9.732
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 https://stackoverflow.com/questions/7523427
复制相似问题