マルコフ連鎖モンテカルロ法 (統計ライブラリー)の1.2.3を参考にしつつ。
rejection samplingでは、target distributionからの乱数を間接的に取ってきているわけですが、今度はxの関数h(x)の期待値に興味があるとしましょう。で、積分が難しいので、サンプリングしてどうこうしまようというのがこのimportance sampling。ここではh(x)=xとして、普通に期待値を計算させてきましょう。g(x)はh(x)p(x)になるべく近く、サンプリング可能であるとする、ということなので、またガンマ関数から取ってこれないとして、コーシー分布から取ってくるというのをやります。定数Mで覆えるというとかいうのはここではありません。
plot(function(x){dgamma(x,shape=10)},0,30) plot(function(x){dcauchy(x,scale=2.5,location=10)},0,30,add=TRUE)
> t <- 1000 > x <- rcauchy(t,scale=2.5,location=10) > w <- dgamma(x,shape=10) / (dcauchy(x,scale=2.5,location=10)) > sum(w*x) / t [1] 10.20139
