読者です 読者をやめる 読者になる 読者になる

Take a Risk:林岳彦の研究メモ

自らの研究に関連するエトセトラについてのメモ的ブログです。主にリスク学と統計学を扱っています。

Rで推定の幅を美しく描く

研究Hacks

前回の以下の記事で描いたグラフをもっと美しく描く方法についてのメモです。

変動性と不確実性について説明してみた - Take a Risk: 林岳彦の研究メモ


前回描いたのはこんなグラフでした。

ここで左図は各パラメータの事後分布の中央値のみを用いて描いた推定分布、右図は中央値だけでなくパラメータの事後分布からのサンプル50個も加えて描いた推定分布です。


右上図と同じ意図の図をものすごく細いライン(lwd=0.02)を用いてサンプル2000個で描いてやるとこう描けます*1

#plotting distribution with posterior samples
plot(x,dnorm(x,mean=median(theta.post.dist),sd=median(sqrt(sigma2.post.dist))),type="l",xlim=c(0,400),ylim=c(0,0.013),xlab="height",ylab="Probability",main="",lwd=0.5,col=2)
for(i in 2:2000){
lines(x,dnorm(x,theta.post.dist[i],sqrt(sigma2.post.dist[i])),type="l",xlim=c(0,400),ylim=c(0,0.013),lwd=0.02,col=4)
}

赤線がパラメータの中央値で描いたグラフ(variabilityのみ)で、青線(によって示される密度)が2000個のパラメータサンプルで描いた推定の幅(variability+uncertaintyを考慮)を表しています。

Rの図ウィンドウ内で推定分布が青く染上がっていく様もなかなかに美しいです。

追記:透過色を使う方法

コメント欄にて久保さんに透過色を使う方法について教えてもらいましたので早速やってみました。同じ図を500個のパラメータサンプルでlwd=1、色をcol="#0000ff10"で描いてみました。

#plotting distribution with posterior samples
plot(x,dnorm(x,mean=median(theta.post.dist),sd=median(sqrt(sigma2.post.dist))),type="l",xlim=c(0,400),ylim=c(0,0.013),xlab="height",ylab="Probability",main="",lwd=0.5,col=2)
for(i in 2:500){
lines(x,dnorm(x,theta.post.dist[i],sqrt(sigma2.post.dist[i])),type="l",xlim=c(0,400),ylim=c(0,0.013),lwd=1,col="#0000ff10")
}

自分には透過色を使うという発想はなかったですね。

実はlwd=0.02みたいな極端な値を使うと図をPDF化したときにうまく反映されないという問題が発生して悩まされたので、lwd=1という普通の値でできるこちらの方法のほうが良いかもしれません。


ありがとうございました!>久保さま



関連文献

Rグラフィックス ―Rで思いどおりのグラフを作図するために―

Rグラフィックス ―Rで思いどおりのグラフを作図するために―

*1:Rのオブジェクトは前回の記事のものを引き継いでいます