2016-10-04 14 views
1

信頼区間で予測をプロットするには少し助けが必要です。私は他の全ての変数は(/モードを意味する)一定に保持されているときvs(0と1)の効果をプロットしたい次の例線形モデル予測の95%信頼区間をプロットするときにプロットが正しくない

library(Hmisc) 
data("mtcars") 

mfit = lm(mpg ~ vs + disp + cyl, data = mtcars) 

#disp and cyl at their mean 
newcar = data.frame(vs = c(0,1), disp = 230, cyl = 6.188) 

pmodel <- predict(mfit, newcar, se.fit=TRUE) 

を考えてみましょう。私はここで間違ってやっている何

plot(1:2, pmodel$fit[1:2], ylim=c(0,1), pch=19, xlim=c(.5,2.5), xlab="X", 
    ylab = "Predicted values", xaxt = "n", main = "Figure1") 
arrows(1:2, (pmodel$fit[1:2] - 1.96 * pmodel$fit[1:2]), 
     1:2, (pmodel$fit[1,1] + 1.96 * pmodel$fit[1:2]), 
     length=0.05, angle=90, code=3) 
axis(1, at=c(1,2), labels=c("0","1")) 

は、私は以下のこのコードを実行するこれを行うには?ありがとう!

+0

: 'ライブラリ(visreg)。 visreg(mfit、xvar = "vs") ' – agenis

答えて

1

ylim = c(0, 1)がありますが、これは間違っています。信頼区間を描くときは、ylimがCIの下限と上限をカバーしていることを確認する必要があります。あなたのコードの

## lower and upper bound of CI 
lower <- with(pmodel, fit - 1.96 * se.fit) 
upper <- with(pmodel, fit + 1.96 * se.fit) 

## x-location to plot 
xx <- 0:1 

## set `xlim` and `ylim` 
xlim <- range(xx) + c(-0.5, 0.5) ## extends an addition 0.5 on both sides 
ylim <- range(c(lower, upper)) 

## produce figure 
plot(xx, pmodel$fit, pch = 19, xlim = xlim, ylim = ylim, xaxt = "n", 
    xlab = "X", ylab = "Predicted values", main = "Figure1") 
arrows(xx, lower, xx, upper, length = 0.05, angle = 90, code = 3) 
axis(1, at = xx) 

enter image description here

他のいくつかのコメント:それはfit - 1.96 * se.fitないfit - 1.96 * fitある

  • 1:2ではなく、0:1のx位置に直接プロットすることができます。
  • fit[1:2]ではなく、fit[1,1]です。 ggplotで
1

:あなたもvisregパッケージのように、まさにこれを行うために設計された機能を使用することができ

df <- data.frame(x=1:2, pred=pmodel$fit[1:2]) 
df$lower <- df$pred - 1.96 * pmodel$se.fit[1:2] 
df$upper <- df$pred + 1.96 * pmodel$se.fit[1:2] 
ggplot(df, aes(x, pred, group=x, col=as.factor(x))) + geom_point(size=2) + 
    geom_errorbar(aes(ymin=lower, ymax=upper), width=0.1) 

enter image description here

関連する問題