MCMCを用いたパラメータ推定[Matlab]

Question

今、MCMCを用いてパラメータ推定のためのプログラムを考えています。Matlabを用いてコードを書いている最中なのですが、IfとEndが合わないという忠告が何個も出てしまいます。下記にコードを載せましたので、何が間違っているか教えていただけませんでしょうか。
求めたいパラメータは、ポアッソン分布に従うと仮定し,尤度を使ってデータの当てはまりを見ています。
Matlab初心者なので間違いも多いと思いますが、どうぞよろしくお願いいたします。

logpt =0; %空値
theta =0; %空値
theta0 = ceil(100.*rand()); % i = initial value of MCMC 100
lambda = 5; % true mean value of dispersal
updatewidth = 0.1; %sigma update width of MCMC
nsample = 100; %number of iterations for MCMC nsample 100000
nthin = 1; %thinning number of MCMC kannbiki, sukashi
individual = 100; %individual numbers 1000
xsize = 10; %space size of matrix
ysize = 10;
dat = zeros(xsize, ysize); %space
xinit = 5; %initial position of population
yinit = 5;

%making dammy data
for j = 1:individual
lam(1) = ceil(100.*rand()); % decide a initial number
P = poissrnd(lambda);
tx = xinit;
ty = yinit;
i = 1;
while i <= P
d = ceil(4.*rand());
if d==1 && ty -1 >0
ty = ty-1;
i = i+1;
elseif 1==2 && tx +1 <= xsize
tx = tx+1;
i = i+1;
elseif d==3 && ty+1 <= ysize
ty = ty+1;
i = i+1;
elseif d==4 && tx-1 > 0
tx = tx-1;
i = i+1;
end;
end;
dat(tx,ty) = dat(tx,ty) + (1);
end;

%MCMC start
for n = 1:nsample
for o = 1:nthin
%(ii) simulation with the previous step dispersal value
sp0 = zeros(xsize, ysize);
ind0 = zeros(2,individual);
for m = 1:individual
r = poissrnd(theta0);
tx = xinit;
ty = yinit;
i = 1;
while i <= r
d = ceil(4.*rand());

if d==1 && ty-1 > 0
ty = ty-1;
i = i+1;
elseif d==2 && tx+1 <= xsize
tx = tx+1;
i = i+1;
elseif d==3 && ty+1 <= ysize
ty = ty+1;
i = i+1;
elseif d==4 && tx-1 > 0
tx = tx-1;
i = i+1;
end;

end;
sp0(tx,ty) = sp0(tx,ty) + 1;
ind0(1,m) = tx;
ind0(2,m) = ty;
end;

%calculating likelihood
LOGP = 0;
for q = 1:xsize
for s = 1:ysize
if sp0(q,s) == 0
LOGP = LOGP + log(poisspdf(0.001,dat(q,s)));
else
LOGP = LOGP + log(poisspdf(sp0(q,s),dat(q,s)));
end;
end;
end;

%(iii) making a proposal value
tmp = updatewidth * randn(size(theta0));
%b= updatewidth = 0.1; sigma update width of MCMC

while tmp < 0
tmp = updatewidth * randn(size(theta0));
end;

%(iv) making proposal distances of each individuals
theta1 = poissrnd(tmp);

%(v) simulation with the proposal value
for p = 1:individual
tx = xinit;
ty = yinit;
i = 1;
while i <= theta1
d = ceil(4.*rand());

if d==1 && ty-1>0
ty = ty-1;
i = i+1;
elseif d==2 && tx+1 <= xsize
tx = tx+1;
i = i+1;
elseif d==3 && ty+1 <= ysize
ty = ty+1;
i = i+1;
elseif d==4 && tx-1 > 0
tx = tx-1;
i = i+1;
end;
end;

%calculating likelihood
prev = sp0(ind0(1,p),ind0(2,p));
dat_prev = dat(ind0(1,p),ind0(2,p));
prop = sp0(tx,ty);
dat_prop = dat(tx,ty);

if prop == 0

LOGP0 = log(poisspdf(prev,dat_prev)) + log(poisspdf(0.001,dat_prop));
else
LOGP0 = log(poisspdf(prev,dat_prev)) + log(poisspdf(prop,dat_prop));
end;
if prev == 1
LOGP1 = log(poisspdf(0.001,dat_prev)) + log(poisspdf(prop+1,dat_prop));
else
LOGP1 = log(poisspdf(prev-1,dat_prev)) + log(poisspdf((prop+1),dat_prop));
end;
if rand < exp(LOGP1 - LOGP0)
r = theta1;
LOGP0 = LOGP1;
end;
end;
%(vi) calculating the mean of theta
theta0 = mean(r);
end;
theta = theta0;
LOGPt = LOGP;
end;

Kules · Accepted Answer

確認ですが、
>IfとEndが合わないという忠告が何個も出てしまいます
というのは「警告」であって「エラー」ではないんですよね？
まあ警告ならほっといてもいいんですが…

もし警告だとすれば、それはifとendのインデントがあっていないことによると思われます。

スペースを_（アンダーバー）で表します。
if 1;
______if 0;
________A=5;
____end;
end;
とすると、おそらく2行目に警告が出ると思われます。
どうやら、ifとendがスペース2つ分以上離れていると出るみたいです。

で、解決方法ですが、Matlabのエディタを使っているのであれば
ctrl+A（全て選択）→ctrl+I（スマートインデント）で解決すると思われます。

参考になれば幸いです。

Kules · Answer

A No.1のKulesです。

警告がなくなったようで安心しました。

>肝心のプログラムが間違っているようで、上手く答えが出てくれません。。。
というのは、エラーが出るということではなく、何か「このような答えが出るはず」というパラメータを入れてもそのような答えが出ないということでしょうか？

残念ながら私の使用しているMatlabにはStatistics Toolboxが入っていないので、
poissrndやpoisspdfが動きません。
したがって、あなたの書いたコードが動くかどうかを確かめる術がありません。

また、MCMCがどのようなものなのかを知らないため、どのように動けば正解なのかもわかりません。

少なくともこのプログラム自体にはエラーはなさそう（poissrnd以降が私の環境で動かないので断言はできませんが）なので、あとはアルゴリズム上の問題かと。

お役に立てず申し訳ないです。

MCMCを用いたパラメータ推定[Matlab]

確認ですが、

A No.1のKulesです。

似たような質問が見つかりました

関連するカテゴリからQ&Aを探す

デイリーランキングこのカテゴリの人気デイリーQ&Aランキング

マンスリーランキングこのカテゴリの人気マンスリーQ&Aランキング