% cobbmcmc.m
% from page 5 of 11/21/09 Workshop handout
dp=[0 3 0 3 3;3 0 3 0 3;0 3 0 3 3;3 0 3 0 3;4 4 4 4 0]
np=ones(size(dp))
p=zeros(size(dp))
i=find(dp~=0)
p(i)=np(i)./dp(i)
Twostepp=zeros(5,5);
for i=1:5
    for j=1:5
        Twostepp(i,j)=sum(p(i,:)*p(:,j))
    end
end
% note that Twostepp=p^2
% Chapman-Kolmogorov equations
k=1;n=4;
p1=p^k*p^(n-k)
p2=p^4
% metropolizing a Markov chain
metropolp=p;
pt=p';
pdiff=zeros(size(p));
pdiff=p-pt;
i=find(pdiff>0);
metropolp(i)=pt(i)
sump=sum(metropolp')
metropolp=metropolp+diag(1-sump)


