From: Dennis <dcflyer_at_gmail.com>

Date: Wed, 23 Jun 2010 15:27:46 -0700

model

{

}

}

}

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 23 Jun 2010 - 22:30:04 GMT

Date: Wed, 23 Jun 2010 15:27:46 -0700

The attached are codes of R and Winbugs.

library(R2WinBUGS) library(MCMCpack) library(coda)

## the input of the dataset

X=read.csv("X.csv",header=FALSE)

X=as.matrix(X) # transform the data into matrix

## parameter setting

N=nrow(X) # # of servers T=ncol(X) # Time m=sum(X)/(N*T) # mean of the training set

M=matrix(m,nrow=N,ncol=T)

s=sum((X-M)^2)/(N*T) # std of the training set

K=3 # # of clusters alpha=0.5 # parameter for Dirichlet distn sigmae=0.5 # var of cluster mean mu q1=rep(1/K,K) # prior for Z(n,1)

## MCMC sampling

data=list("X","m","s","N","T","K","alpha","sigmae","q1")
inits=function(){list(a0=rbeta(1,1,1),

qx=matrix(rgamma(K^2,alpha,1),nrow=K), sigma0.r=rbeta(K,1,1))}
model.sim=bugs(data,inits,model.file="model.txt",parameters=c("mu","sigma"),
n.chains=3,n.iter=3500,n.burnin=500,n.thin=1,bugs.directory="C:/Users/t-wec/Desktop/WinBUGS14",codaPkg=T,debug=T)
mcmcout=read.bugs(model.sim)

summary(mcmcout)

model

{

# cluster parameters mu and tau

tau1 <- (1-a*a)*taue

taue <- 1/sigmae

for (k in 1:K)

{

# cluster mean mu

mu[k,1] ~ dnorm(m,tau1)

for (t in 2:T) { mu[k,t] ~ dnorm(meanmu[k,t],taue) meanmu[k,t] <- m*(1-a)+a*mu[k,t-1] } # cluster varicance tau sigma0.r[k] ~ dbeta(1,1) sigma.r[k] <- s*sigma0.r[k] sigma[k] <- sigma.r[k]*sigma.r[k] tau[k] <- 1/sigma[k]

}

# cluster indicator Z and observation X

for (n in 1:N)

{

Z[n,1] ~ dcat(q1[1:K]) X[n,1] ~ dnorm(mu[Z[n,1],1],tau[Z[n,1]]) for (t in 2:T) { Z[n,t] ~ dcat(q[Z[n,t-1],1:K]) X[n,t] ~ dnorm(mu[Z[n,t],t],tau[Z[n,t]]) }

}

# prior on transition matrix Q

# each row of Q has a Dirichlet prior realized by Gamma

for (k in 1:K)

{

for (l in 1:K) { q[k,l] <- qx[k,l]/sum(qx[k,1:K]) qx[k,l] ~ dgamma(alpha,1) }

}

# prior on regression coefficient: uniform on [-1,1]

a0 ~ dbeta(1,1)

a <- a0*2-1

}

Wei Chen

[[alternative HTML version deleted]]

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 23 Jun 2010 - 22:30:04 GMT

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.2.0, at Thu 24 Jun 2010 - 12:50:35 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*