如何在多级分析中用不同颜色显示不同级别
How to display different levels in a multilevel analysis with different colors
我是多级分析的初学者,我试图了解如何使用 base-R
中的绘图函数绘制图形。我理解下面 fit
的输出,但我在可视化方面苦苦挣扎。 df
只是一些简单的测试数据:
t <- seq(0, 10, 1)
df <- data.frame(t = t,
y = 1.5+0.5*(-1)^t + (1.5+0.5*(-1)^t) * t,
p1 = as.factor(rep(c("p1", "p2"), 10)[1:11]))
fit <- lm(y ~ t * p1, data = df)
# I am looking for an automated version of that:
plot(df$t, df$y)
lines(df$t[df$p1 == "p1"],
fit$coefficients[1] + fit$coefficients[2] * df$t[df$p1 == "p1"], col = "blue")
lines(df$t[df$p1 == "p2"],
fit$coefficients[1] + fit$coefficients[2] * df$t[df$p1 == "p2"] +
+ fit$coefficients[3] + fit$coefficients[4] * df$t[df$p1 == "p2"], col = "red")
它应该知道它必须包含p1
并且有两行。
结果应如下所示:
Edit: Predict est <- predict(fit, newx = t)
给出了与 fit 相同的结果,但我仍然不知道 "how to cluster".
编辑 2 @Keith:公式 y ~ t * p1
读作 y = (a + c * p1) + (b + d * p1) * t
。因为 "first blue line" c, d
都是零。
这就是我要做的。我也包括了一个 ggplot2
版本的情节,因为我发现它更适合我思考情节的方式。
此版本将占 p1
中的级别数。如果您想补偿模型参数的数量,您只需调整构建 xy
的方式以包含所有相关变量。我应该指出,如果您省略 newdata
参数,将对提供给 lm
.
的数据集进行拟合
t <- seq(0, 10, 1)
df <- data.frame(t = t,
y = 1.5+0.5*(-1)^t + (1.5+0.5*(-1)^t) * t,
p1 = as.factor(rep(c("p1", "p2"), 10)[1:11]))
fit <- lm(y ~ t * p1, data = df)
xy <- data.frame(t = t, p1 = rep(levels(df$p1), each = length(t)))
xy$fitted <- predict(fit, newdata = xy)
library(RColorBrewer) # for colors, you can define your own
cols <- brewer.pal(n = length(levels(df$p1)), name = "Set1") # feel free to ignore the warning
plot(x = df$t, y = df$y)
for (i in 1:length(levels(xy$p1))) {
tmp <- xy[xy$p1 == levels(xy$p1)[i], ]
lines(x = tmp$t, y = tmp$fitted, col = cols[i])
}
library(ggplot2)
ggplot(xy, aes(x = t, y = fitted, color = p1)) +
theme_bw() +
geom_point(data = df, aes(x = t, y = y)) +
geom_line()
我是多级分析的初学者,我试图了解如何使用 base-R
中的绘图函数绘制图形。我理解下面 fit
的输出,但我在可视化方面苦苦挣扎。 df
只是一些简单的测试数据:
t <- seq(0, 10, 1)
df <- data.frame(t = t,
y = 1.5+0.5*(-1)^t + (1.5+0.5*(-1)^t) * t,
p1 = as.factor(rep(c("p1", "p2"), 10)[1:11]))
fit <- lm(y ~ t * p1, data = df)
# I am looking for an automated version of that:
plot(df$t, df$y)
lines(df$t[df$p1 == "p1"],
fit$coefficients[1] + fit$coefficients[2] * df$t[df$p1 == "p1"], col = "blue")
lines(df$t[df$p1 == "p2"],
fit$coefficients[1] + fit$coefficients[2] * df$t[df$p1 == "p2"] +
+ fit$coefficients[3] + fit$coefficients[4] * df$t[df$p1 == "p2"], col = "red")
它应该知道它必须包含p1
并且有两行。
结果应如下所示:
Edit: Predict est <- predict(fit, newx = t)
给出了与 fit 相同的结果,但我仍然不知道 "how to cluster".
编辑 2 @Keith:公式 y ~ t * p1
读作 y = (a + c * p1) + (b + d * p1) * t
。因为 "first blue line" c, d
都是零。
这就是我要做的。我也包括了一个 ggplot2
版本的情节,因为我发现它更适合我思考情节的方式。
此版本将占 p1
中的级别数。如果您想补偿模型参数的数量,您只需调整构建 xy
的方式以包含所有相关变量。我应该指出,如果您省略 newdata
参数,将对提供给 lm
.
t <- seq(0, 10, 1)
df <- data.frame(t = t,
y = 1.5+0.5*(-1)^t + (1.5+0.5*(-1)^t) * t,
p1 = as.factor(rep(c("p1", "p2"), 10)[1:11]))
fit <- lm(y ~ t * p1, data = df)
xy <- data.frame(t = t, p1 = rep(levels(df$p1), each = length(t)))
xy$fitted <- predict(fit, newdata = xy)
library(RColorBrewer) # for colors, you can define your own
cols <- brewer.pal(n = length(levels(df$p1)), name = "Set1") # feel free to ignore the warning
plot(x = df$t, y = df$y)
for (i in 1:length(levels(xy$p1))) {
tmp <- xy[xy$p1 == levels(xy$p1)[i], ]
lines(x = tmp$t, y = tmp$fitted, col = cols[i])
}
library(ggplot2)
ggplot(xy, aes(x = t, y = fitted, color = p1)) +
theme_bw() +
geom_point(data = df, aes(x = t, y = y)) +
geom_line()