set.seed(235)
n = 200
tau = 1
sig = 1
k = 2
phi = 0.9
pi = 0.9
mu = rep(0,n+1)
for (t in 2:(n+1))
  mu[t] = rnorm(1,phi*mu[t-1],tau)
lambda = sample(c(1,k),size=n,replace=TRUE,prob=c(pi,1-pi))
y = rnorm(n,mu,lambda*sig)	

plot(y,ylim=range(mu,y),col=lambda,pch=16)
lines(mu)




set.seed(235)
n = 200
tau = 1
sig = 1
k = 2
phi = 0.9
pi = 0.9
mu = rep(0,n+1)
lambda = rep(1,n+1)
for (t in 2:(n+1)){
  mu[t] = rnorm(1,phi*mu[t-1],tau)
  if (lambda[t-1]==1){
    lambda[t] = sample(c(1,k),size=1,prob=c(pi,1-pi))
  }else{
    lambda[t] = sample(c(1,k),size=1,prob=c(1-pi,pi))
  }  
}
y = rnorm(n,mu,lambda*sig)	

plot(y,ylim=range(mu,y),col=lambda,pch=16)
lines(mu)




set.seed(235)
n = 200
tau = 1
sig = 1
phi = 0.9
df  = 3
mu = rep(0,n+1)
for (t in 2:(n+1))
  mu[t] = rnorm(1,phi*mu[t-1],tau)
mu = mu[2:(n+1)]  
y = mu + sig*rt(n,df)

plot(y,ylim=range(mu,y),pch=16)
lines(mu,col=2)