model
{
#distribution of Ys
###################
for (i in 1:N) {
ysigmadet<-exp(th[i,1]+th[i,2])*(1-rhoep*rhoep);
Yisigma2[i,1,1] <- exp(th[i,2])/ysigmadet;
Yisigma2[i,2,2] <- exp(th[i,1])/ysigmadet;
Yisigma2[i,1,2] <- -rhoep*exp(0.5*th[i,1]+0.5*th[i,2])/ysigmadet; Yisigma2[i,2,1] <- Yisigma2[i,1,2];
Y[i,1:2]~ dmnorm(muy[],Yisigma2[i,,]);
}
muy[1]<-0;
muy[2]<-0;
thmean[1,1] <- mu1;
thmean[1,2] <- mu2;
th[1,1]~dnorm(thmean[1,1],itaua2);
th[1,2]~dnorm(thmean[1,2],itaub2);
sig1[1]<-exp(0.5*th[1,1]);
sig2[1]<-exp(0.5*th[1,2]);
q[1]~dnorm(psi0,itau2);
rhoep[1]<-(exp(q[1])-1)/(exp(q[1])+1);
for (i in 2:N) {
thmean[i,1] <- mu1 + phi1*(th[i-1,1]-mu1);
thmean[i,2] <- mu2 + phi2*(th[i-1,2]-mu2);
th[i,1]~dnorm(thmean[i,1],itaua2);
th[i,2]~dnorm(thmean[i,2],itaub2);
sig1<-exp(0.5*th[i,1]);
sig2<-exp(0.5*th[i,2]);
qmean<-psi0+psi*(q[i-1]-psi0);
q[i]~dnorm(qmean[i],itau2);
rhoep[i]<-(exp(q[i])-1)/(exp(q[i])+1);
}
#distribution of phi, mu, rhoep
###########################
phi1star ~ dbeta(20,1.5);
phi1 <- 2*phi1star -1;
phi2star ~ dbeta(20,1.5);
phi2 <- 2*phi2star -1;
psistar ~ dbeta(20,1.5);
psi <- 2*psistar -1;
itaua2 ~ dgamma(2.5,0.025);
taua <- sqrt(1/itaua2);
itaub2 ~ dgamma(2.5,0.025);
taub <- sqrt(1/itaub2);
itau2 ~ dgamma(2.5,0.025);
tau <- sqrt(1/itau2);
mu1 ~ dnorm(0,0.04);
mu2 ~ dnorm(0,0.04);
psi0~dnorm(0.7,0.1);
}
list(phi1star=0.99, phi2star=0.99,mu1=0,mu2=0,itaua2=100,itaub2=100,psistar=0.99,psi0=1.9,itau2=100)