我正在运行一个简单的 lm() 回归，它返回 4 个系数。我想将它们两两绘制（彼此相邻，具有不同的颜色），就好像我会运行两个模型，每个模型给出两个系数。
最小的例子：
Y <- runif(100, 0, 1)
G <- round(runif(100, 0, 1))
C <- round(runif(100, 0, 1))
df <- data.frame(Y, G, C)
out <- lm(Y ~ G + factor(C) + G*factor(C), data=df)
summary(out)
coefplot(out) # from library(fixest)
但是，我希望 (Intercept) 和 G 在 X 轴上彼此相邻的同一位置，以及 factor(C)1 和 G:factor(C)1。一个可能的解决方案可能是将输出与 lm() 分开并将其绘制为：coefplot(list(out1, out1)). 这将如何工作？还有什么办法可以奏效？
1）使用 nlme 对每个 C 级别分别执行回归，然后发出两个 coefplot 调用。
library(nlme)
library(fixest)
L <- lmList(Y ~ G | C, df)
opar <- par(mfrow = 1:2)
ylim <- extendrange(unlist(confint(L)))
coefplot(L[[1]], col = "red", lwd = 2, ylim = ylim)
title(names(L)[1])
coefplot(L[[2]], col = "green", lwd = 2, ylim = ylim)
title(names(L)[2])
par(opar)
（图表后续）
2）或者这就是你想要的。L 和 ylim 来自上方
coefplot(L[[1]], col = "green", x.shift = -0.2, ylim = ylim)
coefplot(L[[2]], add = TRUE, col = "red", x.shift = 0.2)
legend("topright", legend = names(L), col = c("green", "red"), lty = 1)
（图表后续）
3)另一种方法是使用此处显示的模型。
fm <- lm(Y ~ C / G + 0, transform(df, C = factor(C)))
nms <- variable.names(fm)
ix <- grep("G", nms)
g <- list(intercept = nms[-ix], slope = nms[ix])
coefplot(fm, group = g)
以可重现形式使用的输入是
df <- structure(list(Y = c(0.677570635452867, 0.572633401956409, 0.102924682665616,
0.899824970401824, 0.24608773435466, 0.0420595335308462, 0.327920719282702,
0.954503649147227, 0.889539316063747, 0.6928034061566, 0.640506813768297,
0.994269776623696, 0.655705799115822, 0.708530468167737, 0.544066024711356,
0.59414202044718, 0.28915973729454, 0.147113647311926, 0.963024232536554,
0.902299045119435, 0.690705278422683, 0.795467417687178, 0.0246136845089495,
0.477795971091837, 0.758459537522867, 0.216407935833558, 0.318181007634848,
0.231625785352662, 0.142800022382289, 0.414546335814521, 0.413724326295778,
0.368845450924709, 0.152444747742265, 0.13880606344901, 0.233034099452198,
0.465962450252846, 0.265972640365362, 0.857827715342864, 0.0458311666734517,
0.442200074205175, 0.798924845643342, 0.12189925997518, 0.560947983758524,
0.20653138961643, 0.127531650243327, 0.753307864302769, 0.895045359153301,
0.374462775886059, 0.665115194628015, 0.0948406609240919, 0.383969637798145,
0.27438364457339, 0.814640038879588, 0.448516341391951, 0.810064353048801,
0.812389509519562, 0.794342321110889, 0.439831687603146, 0.754475158639252,
0.629221131559461, 0.710182401351631, 0.000624773325398564, 0.475316574098542,
0.220118885161355, 0.379816537722945, 0.612771003274247, 0.351797909243032,
0.111135424347594, 0.243619472719729, 0.66805558744818, 0.417646779678762,
0.788195834029466, 0.102864644257352, 0.434892741497606, 0.984956979984418,
0.893051114398986, 0.886469060787931, 0.175052650272846, 0.130695691565052,
0.653101925039664, 0.343516472261399, 0.656758127966896, 0.320373242488131,
0.187691119266674, 0.782294301316142, 0.0935949867125601, 0.46677904156968,
0.511505459900945, 0.59998895926401, 0.332823540316895, 0.488613033667207,
0.954473827499896, 0.4829023971688, 0.890350222121924, 0.914438186911866,
0.608734982321039, 0.410689776530489, 0.147094690939412, 0.935299803270027,
0.301228899974376), G = c(0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0,
0, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1,
0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0, 1,
1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0,
1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 0,
0, 0, 0, 0)), row.names = c(NA, -100L), class = "data.frame")
你可以推出你的coefand confint。
pdat <- cbind(coef(out), confint(out))
xs <- c(0.6, 1.1, 2.9, 3.4)
plot(xs, pdat[, 1], pch=20, xlim=c(0, 4), ylim=range(pdat), col=2:3, xaxt='n',
main='My Model', xlab='Coefficients', ylab='Beta')
axis(1, xs, F)
mtext(rownames(pdat), 1, .5, at=xs, cex=.7)
arrows(xs, pdat[, 2], xs, pdat[, 3], code=3, angle=90, length=.05, col=2:3)
abline(h=0)
grid()
数据：
set.seed(264211)
df <- data.frame(Y=runif(100, 0, 1), G=round(runif(100, 0, 1)), C=round(runif(100, 0, 1)))
out <- lm(Y ~ G + factor(C) + G*factor(C), data=df)