2016-07-05 4 views
2

私は今後この10年間の回帰をサンプルします。良い線形回帰プロット(適合ライン、信頼/予測バンドなど)を作成する

date<-as.Date(c("2015-12-31", "2014-12-31", "2013-12-31", "2012-12-31")) 
value<-c(16348, 14136, 12733, 10737) 
#fit linear regression 
model<-lm(value~date) 
#build predict dataframe 
dfuture<-data.frame(date=seq(as.Date("2016-12-31"), by="1 year", length.out = 10)) 
#predict the futurne 
predict(model, dfuture, interval = "prediction") 

これに信頼バンドを追加するにはどうすればよいですか?

+0

可能な複製[スロープの95%信頼区間を計算する方法Rの線形回帰モデル](http://stackoverflow.com/questions/15180008/how-to-calculate-the-95-confidence-interval-for-the-slope-in​​-a-linear- regressio) – majom

+0

これは上記のものと重複していません。 OPは信頼区間ではなく信頼区間を求めていた。 – Alex

+0

@ Zheyuan Li私はそれがそれと重複しているとは思わないでしょう。最初に、用語は奇妙に使用されます。もし誰かがそれを探して、正しい言葉を使っているなら、彼はおそらくそれを見つけられません。第二に、ここでは、他の答えに言及されていない線形回帰の特別な場合があります。驚いたことに(私にとっても)、あなたは ''信頼バンド線形回帰 'を探すと良い結果を得られません。だから私はその質問がうまくいくと思う。しかし、私はまた、これにはたくさんの解決策があると確信しています。しかし、あなたがそれらを見つけることができない場合、彼らはどのように役立つことができますか? – Alex

答えて

3

次のコードは、見栄えの良い回帰プロットを生成します。コードに沿った私のコメントは、すべてを明確に説明するべきです。コードはvaluemodelをあなたの質問に使用します。

## all date you are interested in, 4 years with observations, 10 years for prediction 
all_date <- seq(as.Date("2012-12-31"), by="1 year", length.out = 14) 

## compute confidence bands (for all data) 
pred.c <- predict(model, data.frame(date=all_date), interval="confidence") 

## compute prediction bands (for new data only) 
pred.p <- predict(model, data.frame(date=all_date[5:14]), interval="prediction") 

## set up regression plot (plot nothing here; only set up range, axis) 
ylim <- range(range(pred.c[,-1]), range(pred.p[,-1])) 
plot(1:nrow(pred.c), numeric(nrow(pred.c)), col = "white", ylim = ylim, 
    xaxt = "n", xlab = "Date", ylab = "prediction", 
    main = "Regression Plot") 
axis(1, at = 1:nrow(pred.c), labels = all_date) 

## shade 95%-level confidence region 
polygon(c(1:nrow(pred.c),nrow(pred.c):1), c(pred.c[, 2], rev(pred.c[, 3])), 
     col = "grey", border = NA) 

## plot fitted values/lines 
lines(1:nrow(pred.c), pred.c[, 1], lwd = 2, col = 4) 

## add 95%-level confidence bands 
lines(1:nrow(pred.c), pred.c[, 2], col = 2, lty = 2, lwd = 2) 
lines(1:nrow(pred.c), pred.c[, 3], col = 2, lty = 2, lwd = 2) 

## add 95%-level prediction bands 
lines(4 + 1:nrow(pred.p), pred.p[, 2], col = 3, lty = 3, lwd = 2) 
lines(4 + 1:nrow(pred.p), pred.p[, 3], col = 3, lty = 3, lwd = 2) 

## add original observations on the plot 
points(1:4, rev(value), pch = 20) 

## finally, we add legend 
legend(x = "topleft", legend = c("Obs", "Fitted", "95%-CI", "95%-PI"), 
     pch = c(20, NA, NA, NA), lty = c(NA, 1, 2, 3), col = c(1, 4, 2, 3), 
     text.col = c(1, 4, 2, 3), bty = "n") 

regression plot

JPEGは、コードによって生成された:あなたはプロットから見ることができるようにあなたが予測をしようとして、

  • 自信バンドが広く育つ

    jpeg("regression.jpeg", height = 500, width = 600, quality = 100) 
    ## the above code 
    dev.off() 
    ## check your working directory for this JPEG 
    ## use code getwd() to see this director if you don't know 
    

    観察されたデータから遠ざかります。

  • 予測間隔は信頼区間よりも広い。

あなたはpredict.lm()が内部的信頼/予測区間を計算する方法についての詳細をお知りになりたい場合はは、How does predict.lm() compute confidence interval and prediction interval?を読み、そこに私の答え。

visregパッケージの単純な使用のアレックスのデモのおかげで、しかし、私はまだRベースを使用することを好む。あなたが値に興味があるなら

+0

私はこれらが信頼できるバンドではないと主張します。 – Alex

+0

あなたは私とは違う言葉を使うと思います。予測のために、wikipediaはhttps://en.wikipedia.org/wiki/Confidence_and_prediction_bands#Prediction_bandsのように予測バンドと呼ばれます。 「信頼区間」という用語は、ここで説明する回帰パラメータの間隔に使用します(https://en.wikipedia.org/wiki/Confidence_interval)。 – Alex

+0

@ ZheyuanLiはこれらの2つのプロットを1つの連続として表示できますか?私は、あなたのボトムプロットのように始めるが、あなたのトッププロットのように続けていることを意味する? – Oposum

2

あなたは、単にvisreg::visreg

library(visreg) 
visreg(model) 

enter image description here

を使用することができます。

> head(visreg(model)$fit) 
     date value visregFit visregLwr visregUpr 
1 2012-12-31 13434.5 10753.10 9909.073 11597.13 
2 2013-01-10 13434.5 10807.81 9974.593 11641.02 
3 2013-01-21 13434.5 10862.52 10040.033 11685.00 
4 2013-02-01 13434.5 10917.22 10105.389 11729.06 
5 2013-02-12 13434.5 10971.93 10170.658 11773.21 
6 2013-02-23 13434.5 11026.64 10235.837 11817.44 
+0

なぜ、Y軸に5桁の値があるのですか? – Oposum

+0

申し訳ありませんが、Y軸に以前のプロットと比較して異なる値が表示されるのはなぜですか? – Oposum

+0

@ Oposum他のプロットが将来のデータを使用したからだと思います。しかし、これをやっているのであれば、自信を持つバンドではなく、予測バンドを使うべきです。信頼帯は、回帰直線の推定に使用したデータに関連しています。一方、新しいデータには予測バンドを使用する必要があります。 – Alex

関連する問題