绘制具有随机和固定模型的面板数据混合效应模型
Plotting Panel data Mixed Effect model with Random and Fixed models
我正在研究面板数据模型,我现在正在使用 lme4
包中的混合模型,我还使用了基于随机、固定、LSDV、Fisrt_diff 等的模型...
我有一个绘制所有模型系数的函数。在 ggplot 中,但是从 lme4
绘制系数是一个问题,我可以让它工作:
有没有办法让下面的代码适用于所有型号,包括型号 mixed?
library(plm)
library(lme4)
library(ggplot2)
mixed <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
fixed = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "within")
random = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "random")
pool = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "pooling")
first_diff = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "fd")
# Function to extract point estimates
ce <- function(model.obj) {
extract <- summary(get(model.obj))$coefficients[2:nrow(summary(get(model.obj))$coefficients), 1:2]
return(data.frame(extract, vars = row.names(extract), model = model.obj))
}
# Run function on the three models and bind into single data frame
coefs <- do.call(rbind, sapply(paste0(list(
"fixed", "random", "pool", "first_diff"
)), ce, simplify = FALSE))
names(coefs)[2] <- "se"
gg_coef <- ggplot(coefs, aes(vars, Estimate)) +
geom_hline(yintercept = 0, lty = 1, lwd = 0.5, colour = "red") +
geom_errorbar(aes(ymin = Estimate - se, ymax = Estimate + se, colour = vars),
lwd = 1, width = 0
) +
geom_point(size = 3, aes(colour = vars)) +
facet_grid(model ~ ., scales="free") +
coord_flip() +
guides(colour = FALSE) +
labs(x = "Coefficient", y = "Value") +
ggtitle("Raw models coefficients")
gg_coef
您当前代码的错误是
data(sleepstudy)
mixed <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
coefficients(summary(mixed))
Estimate Std. Error t value
(Intercept) 251.40510 6.823773 36.842535
Days 10.46729 1.545958 6.770744
Days 在 sleepstudy 数据集中是数字,并使用了连续预测变量。使用你的 ce 函数,这 returns 一个错误,因为行名被删除, 2:nrow(..).
要获得与您的其他模型相似的估计值,请将天数设置为因数并将随机效应设置为(1|天)。我认为(天数 | 主题)没有意义。
sleepstudy$Days = factor(sleepstudy$Days)
mixed <- lmer(Reaction ~ Days + (1 | Subject), sleepstudy)
并且我们使用 drop=FALSE 稍微更改了您的 ce 代码,以防止空 row.names
ce <- function(model.obj) {
summ.model <- summary(get(model.obj))$coefficients
extract <- summ.model[2:nrow(summ.model),drop=FALSE, 1:2]
return(data.frame(extract, vars = row.names(extract), model = model.obj))
}
coefs <- do.call(rbind, sapply(paste0(list(
"fixed", "random", "pool", "first_diff","mixed"
)), ce, simplify = FALSE))
names(coefs)[2] <- "se"
运行你剩下的:
gg_coef <- ggplot(coefs, aes(vars, Estimate)) +
geom_hline(yintercept = 0, lty = 1, lwd = 0.5, colour = "red") +
geom_errorbar(aes(ymin = Estimate - se, ymax = Estimate + se, colour = vars),
lwd = 1, width = 0
) +
geom_point(size = 3, aes(colour = vars)) +
facet_grid(model ~ ., scales="free") +
coord_flip() +
guides(colour = FALSE) +
labs(x = "Coefficient", y = "Value") +
ggtitle("Raw models coefficients")
gg_coef
我正在研究面板数据模型,我现在正在使用 lme4
包中的混合模型,我还使用了基于随机、固定、LSDV、Fisrt_diff 等的模型...
我有一个绘制所有模型系数的函数。在 ggplot 中,但是从 lme4
绘制系数是一个问题,我可以让它工作:
有没有办法让下面的代码适用于所有型号,包括型号 mixed?
library(plm)
library(lme4)
library(ggplot2)
mixed <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
fixed = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "within")
random = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "random")
pool = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "pooling")
first_diff = plm(Reaction ~ Days, data = sleepstudy, index = c("Subject", "Days"), model = "fd")
# Function to extract point estimates
ce <- function(model.obj) {
extract <- summary(get(model.obj))$coefficients[2:nrow(summary(get(model.obj))$coefficients), 1:2]
return(data.frame(extract, vars = row.names(extract), model = model.obj))
}
# Run function on the three models and bind into single data frame
coefs <- do.call(rbind, sapply(paste0(list(
"fixed", "random", "pool", "first_diff"
)), ce, simplify = FALSE))
names(coefs)[2] <- "se"
gg_coef <- ggplot(coefs, aes(vars, Estimate)) +
geom_hline(yintercept = 0, lty = 1, lwd = 0.5, colour = "red") +
geom_errorbar(aes(ymin = Estimate - se, ymax = Estimate + se, colour = vars),
lwd = 1, width = 0
) +
geom_point(size = 3, aes(colour = vars)) +
facet_grid(model ~ ., scales="free") +
coord_flip() +
guides(colour = FALSE) +
labs(x = "Coefficient", y = "Value") +
ggtitle("Raw models coefficients")
gg_coef
您当前代码的错误是
data(sleepstudy)
mixed <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
coefficients(summary(mixed))
Estimate Std. Error t value
(Intercept) 251.40510 6.823773 36.842535
Days 10.46729 1.545958 6.770744
Days 在 sleepstudy 数据集中是数字,并使用了连续预测变量。使用你的 ce 函数,这 returns 一个错误,因为行名被删除, 2:nrow(..).
要获得与您的其他模型相似的估计值,请将天数设置为因数并将随机效应设置为(1|天)。我认为(天数 | 主题)没有意义。
sleepstudy$Days = factor(sleepstudy$Days)
mixed <- lmer(Reaction ~ Days + (1 | Subject), sleepstudy)
并且我们使用 drop=FALSE 稍微更改了您的 ce 代码,以防止空 row.names
ce <- function(model.obj) {
summ.model <- summary(get(model.obj))$coefficients
extract <- summ.model[2:nrow(summ.model),drop=FALSE, 1:2]
return(data.frame(extract, vars = row.names(extract), model = model.obj))
}
coefs <- do.call(rbind, sapply(paste0(list(
"fixed", "random", "pool", "first_diff","mixed"
)), ce, simplify = FALSE))
names(coefs)[2] <- "se"
运行你剩下的:
gg_coef <- ggplot(coefs, aes(vars, Estimate)) +
geom_hline(yintercept = 0, lty = 1, lwd = 0.5, colour = "red") +
geom_errorbar(aes(ymin = Estimate - se, ymax = Estimate + se, colour = vars),
lwd = 1, width = 0
) +
geom_point(size = 3, aes(colour = vars)) +
facet_grid(model ~ ., scales="free") +
coord_flip() +
guides(colour = FALSE) +
labs(x = "Coefficient", y = "Value") +
ggtitle("Raw models coefficients")
gg_coef