2016-12-14 96 views
3

我是多級分析的初學者,並嘗試瞭解如何使用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,並且有兩條線。
結果應該是這樣的: Fit with two levels (interactions needs to be included)

編輯:預測est <- predict(fit, newx = t)給出了相同的結果擬合但我仍然不知道「如何羣集」。

編輯2 @Keith:公式y ~ t * p1讀數爲y = (a + c * p1) + (b + d * p1) * t。對於「第一條藍線」c, d都是零。

+0

將't'和'y'永遠是一樣的嗎?你會總是有一個單一的因素列? –

+0

獲得擬合值的更好方法是通過'predict()'。傳遞一個data.frame作爲'newdata'來預測和驚訝。 –

+0

@KithithHughitt不,我真正的問題會更復雜,爲簡化的示例看[這裏](http://stackoverflow.com/questions/40765869/analysis-using-linear-regression-based-on-subgroups)。這個問題的主要目的是「如何用正確的行數獲得圖表」。即使在這種簡單的情況下('p1'只有2級!)我失敗了;-) – Christoph

回答

3

這就是我該怎麼做的。我還包括一個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]) 
} 

enter image description here

library(ggplot2) 
ggplot(xy, aes(x = t, y = fitted, color = p1)) + 
    theme_bw() + 
    geom_point(data = df, aes(x = t, y = y)) + 
    geom_line() 

enter image description here