0

関数に事前分布情報を入力したい。関数本体を変更することで個々の分布を手動で入力できますが、これを行う一般的な方法を探していますか? たとえば、事前分布が与えられた関数の事後分布をプロットしたいとします。

set.seed(1)
n <- 10
pars <- runif(n)
y <- NA
for (i in 1:n)
  y[i] <- rbinom(1,1, prob=pars[i])

plotPosterior <- function(pars,y,mean=0,vari=4)
{  
  x <- seq(-3,3,by = .1)
  logLik <- NA
  for (i in seq(along.with=x))
    logLik[i] <- sum(y*log(1 + exp(pars-x[i])) - (y-1)*log(1 + exp(x[i]-pars)))

  posterior <- logLik * dnorm(x,mean=mean,sd=sqrt(vari))
  plot(x,posterior,type="l")
}
plotPosterior(pars,y,0,4)

正規分布の平均分散パラメータを入力できます。しかし、たとえばベータ配布を使用したい場合は、関数を書き直す必要があります。dnorm(mean=xx,sd=yy)代わりに、「 」や「 」などの分布を入力する方法が必要dbeta(shape1=xx, shape2=yy)です...
私が見る唯一の実行可能な方法はdnorm(x,mean=mean,sd=sqrt(vari))、入力として機能に入る方法です。しかし、事前に事前に指定したくありませんx。これを行う他の方法はありますか？

4

1 に答える 1

1

明確にするために、コメントから抽出された実用的なソリューションを次に示します。

set.seed(1)
n <- 10
pars <- runif(n)
y <- NA
for (i in 1:n)
  y[i] <- rbinom(1,1, prob=pars[i])

plotPosterior <- function(pars,y, fun = dnorm, 
                          params.fun = list(mean=0, sd=2))
{  
  x <- seq(-3,3,by = .1)
  logLik <- NA
  for (i in seq(along.with=x))
    logLik[i] <- sum(y*log(1 + exp(pars-x[i])) - (y-1)*log(1 + exp(x[i]-pars)))

  posterior <- logLik * do.call(fun, c(list(x), params.fun))
  plot(x,posterior,type="l")
}

plotPosterior(pars, y) # default params and function
plotPosterior(pars, y, fun = dbeta, params.fun = list(shape1=2, shape2=3))
于 2013-06-18T21:32:19.237 に答える