0

f= theta ^(z_f+n+alpha-1)*(1-theta)^(n+1-z_f-k+ beta-1)theta以外のすべてのパラメーターがわかっているパラメーター theta の分布をシミュレートしようとしています。MCMC シミュレーションを実行するために、メトロ ポリッシュ ヘイスティング アルゴリズムを使用しています。私の提案密度は、パラメーター alpha と beta を持つベータ分布です。シミュレーション用の私のコードは次のとおりです。この目的のために、mhsample() という名前の組み込みの Matlab コードを使用しています。

clear
clc
alpha=2;
beta=2;
z_f=1;
n=6;
k=5;

nsamples = 3000;
pdf= @(x) x^(z_f+n+alpha-1)*(1-x)^(n+1-z_f-k+beta-1); % here x acts as theta
proppdf= @(x,y) betapdf(x, alpha, beta);
proprnd =@(x) betarnd(alpha,beta,1);

smpl = mhsample(0.1,nsamples,'pdf',pdf,'proprnd',proprnd,'proppdf',proppdf);
4

1 に答える 1

0

「コードが正常に動作しているかどうかを確認するにはどうすればよいですか」という質問の意味がよくわかりません。コードが実行されると仮定していますか? ただし、関数とシミュレーションを視覚的に比較するには、PDF と mhsample から取得したデータの両方を次のようにプロットできます。

% i'm assuming you ran the code above so that smpl and @pdf are both defined...

fplot(pdf,[0 1]); % fplot takes your function and plots it between x-limit [0,1]
figure % new figure
hist(smpl,30); % 30 here is bin size, change it to your preference

下の図:

  • smpl左の の出力、つまりシミュレーションのヒストグラム
  • pdfシミュレーションと比較するための右側の [0,1] に制限された関数

左が smpl の出力、右が [0,1] で囲まれた pdf

これらの 2 つの図は互いに似ており、ベータ版配布のようなものでもあるため、これは単なる推測にすぎません。

それよりも複雑な分析が必要な場合は、残念ながら私はまだ MCMC に習熟していません :)

于 2014-08-26T23:07:36.723 に答える