| rwupdate_move(y , yt , x , cdfperio_y , w , tp_wk , D , DtD , M , MDt , ytPy , delta , t , hyper , k , IN , out);
|
function [K , tp_wk , D , DtD , M , MDt , ytPy , out] = rwupdate_move(y , yt , x , cdfperio_y , w , tp_wk , D , DtD , M , MDt , ytPy , delta , t , hyper , k , IN , out);
%
% [K , tp_wk , D , DtD , M , MDt , ytPy , out] = rwupdate_move3(y , yt , x , cdfperio_y , w , tp_wk , D , DtD , M , MDt , ytPy , delta , t , hyper , k , IN , out);
%
N = length(y);
k1 = k(t) + 0;
k2 = k(t) + 0;
inv_delta = 1/delta(t);
std_sRW = sqrt(hyper.sRW);
cte1 = pi + hyper.walk;
cte2 = 0 - hyper.walk;
cte_hyper = (hyper.upsilon + N)*0.5;
Ik = inv_delta*eye(k2);
u = rand(k1 , 1);
v = rand(k1 , 1);
%delta_w = w(2) - w(1);
for j = 1 : k1
Dt = D';
DtD = Dt*D;
invM = DtD + Ik;
detinvM = det(invM);
M = inv(invM);
MDt = M*Dt;
ytPy = (yt*(IN - D*MDt))*y;
if (u(j) < hyper.lambda)
proposal = w(sum(rand > cdfperio_y)); % + delta_w*rand ; %M-H step 1 : q1(w'_{j,k}|w_{j,k})
%<-- bin index -------->%
%wk' drawn uniformely inside the frequency
%frontiers associated with the bin%
% proposal = min([proposal ; cte1]);
%
% proposal = max([proposal ; cte2]);
else
proposal = tp_wk(j) + std_sRW*randn; %M-H step 2 : q2(w'_{j,k}|w_{j,k})
proposal = min([proposal ; cte1]);
proposal = max([proposal ; cte2]);
end
% if( ~all(tp_wk - proposal) )
%
% tp_wk(j) = tp_wk(j);
%
% else
Dproposal = D;
Dproposal(: , j) = bsinu(proposal , x);
Dproposalt = Dproposal';
DtDproposal = Dproposalt*Dproposal;
invMproposal = DtDproposal + Ik;
detinvMproposal = det(invMproposal);
Mproposal = inv(invMproposal);
MproposalDt = Mproposal*Dproposalt;
ytPproposaly = (yt*(IN - Dproposal*MproposalDt))*y;
ratio = sqrt(abs(detinvM/detinvMproposal))*((hyper.gamma + ytPy)/(hyper.gamma + ytPproposaly))^(cte_hyper);
acceptance = min(1 , ratio);
if (v(j) < acceptance)
tp_wk(j) = proposal;
D = Dproposal;
DtD = DtDproposal;
M = Mproposal;
MDt = MproposalDt;
ytPy = ytPproposaly;
out.arwupdate = out.arwupdate + 1;
else
out.rrwupdate = out.rrwupdate + 1;
end
%tp_wk = sort(tp_wk);
% end
end
K = k1;
|
|