# WinBUGS在统计分析中的应用（第一部分）

## 第一节 什么是WinBUGS?

WinBUGS对于研究Bayesian统计分析的人来说，应该不会陌生。至少对于MCMC方法是不陌生的。WinBUGS (Bayesian inference Using Gibbs Sampling）就是一款通过MCMC方法来分析复杂统计模型的软件。其基本原理就是通过Gibbs sampling和Metropolis算法，从完全条件概率分布中抽样，从而生成马尔科夫链，通过迭代，最终估计出模型参数。引入Gibbs抽样与MCMC的好处是不言而喻的，就是想避免计算一个具有高维积分形式的完全联合后验概率公布，而代之以计算每个估计参数的单变量条件概率分布。具体的算法思想，在讲到具体问题的时候再加以叙述，在此不过多论述。就不拿公式出来吓人了（毕竟打公式也挺费劲啊）。

## 第三节 如何得到WinBUGS?

WinBUGS目前是一款免费的软件，去http://www.mrc-bsu.cam.ac.uk/bugs/下载就好了。不过要用高级功能（如GeoBUGS）的话，还是去http://www.mrc-bsu.cam.ac.uk/bugs/winbugs/contents.shtml注册一下好了，挺方便的。系统会立即把注册码发到你邮箱（真是好人啊）。不过只可以用一个月。这倒无妨，到时在注册一下就好了。

## 第四节 初试WinBUGS

#MODEL
model
{
for (i in 1:N) {
O[i] ~ dpois(mu[i])
log(mu[i]) <- log(E[i]) + alpha0 + alpha1 * X[i]/10 +
b[i]
# Area-specific relative risk (for maps)
RR[i] <- exp(alpha0 + alpha1 * X[i]/10 + b[i])
}
# CAR prior distribution for random effects:
b[1:N] ~ car.normal(adj[], weights[], num[], tau)
for (k in 1:sumNumNeigh) {
weights[k] <- 1
}
# Other priors:
alpha0 ~ dflat()
alpha1 ~ dnorm(0, 1e-05)
tau ~ dgamma(0.5, 5e-04)
# prior on precision
sigma <- sqrt(1/tau)
# standard deviation
}
#DATA

list(N = 56, O = c(9, 39, 11, 9, 15, 8, 26, 7, 6,
20, 13, 5, 3, 8, 17, 9, 2, 7, 9, 7, 16, 31, 11, 7, 19, 15,
7, 10, 16, 11, 5, 3, 7, 8, 11, 9, 11, 8, 6, 4, 10, 8, 2,
6, 19, 3, 2, 3, 28, 6, 1, 1, 1, 1, 0, 0), E = c(1.4, 8.7,
3, 2.5, 4.3, 2.4, 8.1, 2.3, 2, 6.6, 4.4, 1.8, 1.1, 3.3, 7.8,
4.6, 1.1, 4.2, 5.5, 4.4, 10.5, 22.7, 8.8, 5.6, 15.5, 12.5,
6, 9, 14.4, 10.2, 4.8, 2.9, 7, 8.5, 12.3, 10.1, 12.7, 9.4,
7.2, 5.3, 18.8, 15.8, 4.3, 14.6, 50.7, 8.2, 5.6, 9.3, 88.7,
19.6, 3.4, 3.6, 5.7, 7, 4.2, 1.8), X = c(16, 16, 10, 24,
10, 24, 10, 7, 7, 16, 7, 16, 10, 24, 7, 16, 10, 7, 7, 10,
7, 16, 10, 7, 1, 1, 7, 7, 10, 10, 7, 24, 10, 7, 7, 0, 10,
1, 16, 0, 1, 16, 16, 0, 1, 7, 1, 1, 0, 1, 1, 0, 1, 1, 16,
10), num = c(3, 2, 1, 3, 3, 0, 5, 0, 5, 4, 0, 2, 3, 3, 2,
6, 6, 6, 5, 3, 3, 2, 4, 8, 3, 3, 4, 4, 11, 6, 7, 3, 4, 9,
4, 2, 4, 6, 3, 4, 5, 5, 4, 5, 4, 6, 6, 4, 9, 2, 4, 4, 4,
5, 6, 5), adj = c(19, 9, 5, 10, 7, 12, 28, 20, 18, 19, 12,
1, 17, 16, 13, 10, 2, 29, 23, 19, 17, 1, 22, 16, 7, 2, 5,
3, 19, 17, 7, 35, 32, 31, 29, 25, 29, 22, 21, 17, 10, 7,
29, 19, 16, 13, 9, 7, 56, 55, 33, 28, 20, 4, 17, 13, 9, 5,
1, 56, 18, 4, 50, 29, 16, 16, 10, 39, 34, 29, 9, 56, 55,
48, 47, 44, 31, 30, 27, 29, 26, 15, 43, 29, 25, 56, 32, 31,
24, 45, 33, 18, 4, 50, 43, 34, 26, 25, 23, 21, 17, 16, 15,
9, 55, 45, 44, 42, 38, 24, 47, 46, 35, 32, 27, 24, 14, 31,
27, 14, 55, 45, 28, 18, 54, 52, 51, 43, 42, 40, 39, 29, 23,
46, 37, 31, 14, 41, 37, 46, 41, 36, 35, 54, 51, 49, 44, 42,
30, 40, 34, 23, 52, 49, 39, 34, 53, 49, 46, 37, 36, 51, 43,
38, 34, 30, 42, 34, 29, 26, 49, 48, 38, 30, 24, 55, 33, 30,
28, 53, 47, 41, 37, 35, 31, 53, 49, 48, 46, 31, 24, 49, 47,
44, 24, 54, 53, 52, 48, 47, 44, 41, 40, 38, 29, 21, 54, 42,
38, 34, 54, 49, 40, 34, 49, 47, 46, 41, 52, 51, 49, 38, 34,
56, 45, 33, 30, 24, 18, 55, 27, 24, 20, 18), sumNumNeigh = 234)
#INITIAL VALUES
list(tau = 1, alpha0 = 0, alpha1 = 0, b = c(0, 0,
0, 0, 0, NA, 0, NA, 0, 0, NA, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0))

WinBUGS在统计分析中的应用 第一部分完

## 《WinBUGS在统计分析中的应用（第一部分）》有74个想法

1. WinBugs，头一回听说啊。老齐在annpro，有空分享一下工作心得吧。

2. 刘思喆说道：

土了，还以为是 windows 下的 debug 工具。原来是这 ！！

3. 呵呵，能写点公式还是写点吧，看代码不知道做的问题是什么……

4. 齐韬说道：

我会在第二部分中将相关的理论部分加上，看来得要学习一下怎么嵌入LaTex了.

1. lemonc说道：

请问如何输入幂的形式的，比如a^b，但是中间那个符号不能识别

1. 卷卷说道：

”请问如何输入幂的形式的，比如a^b，但是中间那个符号不能识别“ 请问这个问题您解决了没，求教

5. 王化儒说道：

嗯，楼主辛苦了！跟着一点一点学，期待着空间数据分析。。。

6. 贴一首诗吧，刚BUGS team发过来的，非常有意思：

Each year you wait with bated breath,
The old WinBUGS key nearing death,
And will the brand new key appear
In time to join the festive cheer?
The waiting’s over – raise your glass,
And drink to rituals that pass.
Relax, sit back and have a chortle;
This time your WinBUGS key’s immortal.

7. xjx说道：

不知道在linux下能不能用?

1. tao liu说道：

另外，R里面有 R2bugs 的package，在数据不大的情况下，很方便。

8. hong说道：

大家好

麻烦问一下这个迭代次数如何选择.谢谢回答

9. DJ说道：

很好的教程，谢谢了。MCMC万岁！

10. 左伊秩訾说道：

前几天上课时听说了这个软件，真是及时雨！谢谢了！

11. icwei说道：

小弟现正在剑桥mrc-bsu做postdoc，具体的项目就是BUGS的开发以及在生物及医学方面应用。WinBUGS这个软件是我两个老板David Spiegehalter,Dave Lunn 和其它一些牛人共同开发的。我们现在正在从WinBUGS 转向openBUGS,目的是将它做成open source的软件以应用在更广的领域。 我现在正在开发BUGS中的WBDiff部分并将它应用在二型糖尿病的动态系统的数据分析中。有兴趣或有问题的同学可以和我联系：

chen.wei@mrc-bsu.cam.ac.uk

还有我们这里每年会举办3-4次BUGS的培训，2天的课程，在英国或能到英国出差的同学有兴趣的话可以参加，主讲人是David speigehalter 和 Dave Lunn。

1. DJ说道：

楼上的大牛人啊！

2. maple说道：

这位仁兄真牛

3. 海涛说道：

请问在
winbugs中我编译FRAILTY模型时提示
educational version cannot do this model 难道这个软件还有别的版本？？？我用的版本是1.4
我编译的模型是帮助文件Examples Volume 1中的最后一个文件Cox regression with random effects

12. zyn说道：

希望能和您联系,我这边看了gibbs sampling有不少问题，不知道您可否提供帮助？

13. phlissia说道：

弱弱的问一句，执行到第四步check model的时候，winbugs坐下角显示的不是“model is syntactically correct”，而是在alpha1 ~ dnorm(0, 1e-05)一句1e处显示e处应为”expected right parenthesis，在代码最后0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0))指示”invalid or unexpected token scanned”，
这要怎么弄呢？ 第一次接触这个软件，还不会使，请高手指点。

1. 数字格式不对。参见下面的回复。

To齐韬：你可能需要修改正文中的代码，我在OpenBUGS 3.0.7中运行你的代码会出错，错误信息和这里一样。我检查了一下，是数字格式的问题。

1. 王力说道：

请问就是把您给的代码复制进WINBUGS,到compile时,出现educational version cannot do this model.什么原因阿?对这个软件不甚掌握,论文却急用请指教,谢谢

2. 用OpenBUGS吧，我从不用WinBUGS，后者早已经停止开发了

14. peng说道：

yes I have the same problem with ”expected right parenthesis”

15. peng说道：

把 alpha1 ~ dnorm(0, 1e-05)
tau ~ dgamma(0.5, 5e-04)
替换成

alpha1 ~ dnorm(0, 0.00001)
tau ~ dgamma(0.5, 0.0005)

这样就可以了

1. 出错的原因是数字的格式不对，应该写作1.0E-05（1是整数，1.0才是浮点数；e要大写）

这里设置的正态分布实际上应该是想要一个flat prior，所以用dflat()就可以了。

1. 关于MH算法说道：

您好，请问进行MH算法需要另外编辑程序吗？

2. 不需要，要是MH都要编程，那WinBUGS就没有存在必要了……呃，确切地说，我说的是OpenBUGS……我从没用过WinBUGS，因为它已经确定没有继续的技术支持了，它的开发已经停止。OpenBUGS还在更新中。

3. yin说道：

谢老师，上面的问题急待回答，谢谢。

16. peng说道：

请问winbugs能计算递归吗？想用之解决一个随机变量自我n重卷积的问题。用递归来算，比如一个正态分布s~(mu,tau)，自我卷积n次求结果。该怎么弄呢？谢谢大虾们指教。

17. peng说道：

请问winbugs中有卷积函数吗？

18. yuanhang说道：

另一个可能性是你的模型中使用了数据没有提供的变量名，这种情况下，你要自己检查模型中的变量，如果有遗漏，就在数据的list()中加上。

19. smw说道：

model
{
for (i in 1 : I) {
for(j in 1 : J) {
for(k in 1:K) {
epsilon[i,j,k] ~ dnorm(0.0,tau[i])
y[i,j,k] <-alpha[i,k] + beta[i,k] * theta[j,k] + epsilon[i,j,k]
}
}
}

list(
I = 5,J = 5,K = 4,
y = structure(
.Data = c(
6.8,7.2,7.35,7.63,

.Dim = c(5,5,4)
)
)

编译出现multiple definitions of node y[1,1,1]
如果去掉数据部分的y定义就可以编译过去。
希望大家帮我解决下这个问题，万分感谢！

1. 去掉y即使模型语法正确，程序肯定也无法运行。

你再仔细看看你这里y的定义，c的回括号打到哪里去了？

还有，一共4个数字，请问该如何定义为一个5x5x4的数组？

2. tsdy2001说道：

请问，上面这位仁兄，你的问题解决了吗？我现在也遇到与你同样的问题。到底问题出在哪里呢？

20. 韩涵说道：
model
{
#distribution of Ys
###################
for (i in 1:N) {
Yisigma2[i,1,2] <- -rhoep[i]*exp(0.5*th[i,1]+0.5*th[i,2])/ysigmadet[i];             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[i]<-exp(0.5*th[i,1]);
sig2[i]<-exp(0.5*th[i,2]);
qmean[i]<-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)
1. copious说道：

这个是什么的程序呀 SO 长

21. 韩涵说道：

这是一个DC-MSV模型，大家能帮我看一下，到底是哪个初始值没有定义吗？
在编译完成之后，提示说有未定义的初始值，是怎么回事？
急急急，先谢过各位了！！！！！

1. 这是个体力活儿，自个儿把模型中所有的参数都列出来，然后对比一下是不是所有的参数都有初始值吧

1. whatever说道：

model;
{
mu ~ dnorm(1,0.4)
for( i in 2 : n ) {
v[i] ~ dnorm(vmean[i],ivd[i])I(0,50)
}
kt ~ dnorm( 0.0, 1.0)
k ~ dnorm( 0.0, 1.0)
# theta <- kt / k
for( i in 2 : n ) {
ksy[i] ~ dnorm(ksymean[i],tauy)
}
for( i in 2 : n ) {
ksv[i] ~ dexp(muv)
}
rhoj ~ dnorm( 0.0,4)
muv ~ dgamma(1,0.1)
for( i in 2 : n ) {
j[i] ~ dbern(lambda)
}
lambda ~ dbeta(1,20)
#ev ~ dnorm( 0.0, 1.0)
#ey ~ dnorm( 0.0, 1.0)
#e ~ dnorm( 0.0,1.0E-6)
sigy2 <- 1 / tauy
tauy ~ dgamma( 1,0.1)
muy ~ dnorm( 0.0,0.1)
for( i in 2 : n ) {
#rmean[i] <- mu + ksymean[i] * j[i]
rmean[i] <- mu + ksy[i] * j[i]
}
tauv ~ dgamma( 1, 0.1)
sigv2 <- 1 / tauv
for( i in 2 : n ) {
#vmean[i] <- v[i – 1] +( kt-k* v[i – 1]) + 1/muv * j[i] + rho * (y[i] – y[i – 1] – mu – ksymean[i] * j[i])
vmean[i] <- v[i – 1] +( kt-k* v[i – 1]) +ksv[i] * j[i] +sqrt(sigv2)* rho * (y[i] – y[i – 1] – mu – ksy[i] * j[i])
}
for( i in 2 : n ) {
#ivd[i] <- 1 / (rho*rho*sigy2*j[i]+(1-rho*rho)*sigv2*j[i]+v[i – 1]+1/pow(muv,2)*j[i])
ivd[i] <- tauv/((1-rho*rho)*v[i-1])
}
for( i in 2 : n ) {
#ird[i] <- 1 / (sigy2*j[i]+v[i – 1])
ird[i] <- 1/v[i-1]
}
for( i in 2 : n ) {
r[i] <- y[i]-y[i-1]
r[i] ~ dnorm(rmean[i],ird[i])
}
rho ~ dunif(-1,1)
for( i in 2 : n ) {
ksymean[i] <- muy + rhoj * ksv[i]
}
v[1] <- 1
ksv[1] <- 0.0
ksy[1] <- 0.0
j[1]<-1
}
list(y=c(3.274324143,3.287163251,3.284207325,3.273737212,3.279863062),
n=5)
list(mu=1,k=1,kt=1,tauv=1,rho=1,lambda=1,muy=1,tauy=1,muv=1,rhoj=1)

2. 刘宇佳说道：

您好，我也在做这个dcc-msv模型。代码部分也是这个。可是我到了compile这一步以后就一直提示“vairable N is not defined”.请问你遇到了这个问题吗？谢谢

22. Wangsen Lan说道：

每次启动winBUGS1.4.3都要弹出Licence Agreement，这是不是没安装好呢？

1. 用OpenBUGS！！WinBUGS早已经停止开发了。OpenBUGS开源、免费、还在开发更新中。

23. Wangsen Lan说道：

第七步时，它说edicational version cannot do this model，是不是必须注册了才行呀

1. hh说道：

请问现在注册网址是什么？上面楼主的网址好像有错误

24. copious说道：

大家好 OpenBUGS和WINbugs的程序一样吗

25. Zhangtao说道：

学《高级统计学》的时候，在介绍MCMC及EM的时候，老师介绍过这个软件，
非常好！就是我自己编不出程序来。我是统计的爱好老师与初学者。
另外，问一下谢老师，openbugs在那儿有下载的呢？
非常感谢！

1. 都这年头了，真的需要我把网址打出来么？……

2. openbugs说道：

26. yin说道：

请问谢老师，我的问题也是
它说edicational version cannot do this model，是什么原因呢？谢谢了

1. 前面回复了一遍又一遍：请用OpenBUGS。等我抽时间写一篇OpenBUGS的介绍，实际上它和WinBUGS并没有太大的不同。

1. 张为豪说道：

你好，我是WINBUGS的初学者，现在想用来做SV-T模型的MCMC估计，但总提示数据载入不了，还有初始值设置也不明白，请帮忙指导一下，真的很急，谢谢！！！

model{
for(i in 1:n)
{y[i]~dt(0,p[i],omega)
p[i]<-exp(-theta[i])
}
theta[1]~dnorm(mu,itau2)
for(j in 2:n)
{theta[j]~dnorm(theta2[j],itau2)
theta2[j]<-mu+phi*(theta[j-1]-mu)}
phi<-2*phi1-1
tau<-sqrt(1/itau2)
mu~dnorm(0,0.01)
itau2~dgamma(2.5,0.025)
phi1~dbeta(20,1.5)
omega~dchisqr(8)
}

list(N=236， y=c(-2.863137825，0.236862576，0.708743605，-0.410171084，-0.197900711，-0.017711569，-1.099961877，
-0.131942478，-0.628070303，0.329370517，-0.425979456，0.003865601，-1.562805432，-1.029900337，0.30855106，
-1.76405028，0.536937996，0.88903611，-0.478923481，-2.231280334，0.95179952，-0.240613751，-0.444776657，
0.238243297，2.343495649，-1.345814335，0.409521478，0.938909523，-0.541400821，-1.120087252，-0.465308369，
0.22970175，-0.442665133，0.208833427，-0.866786092，0.294963792，1.018515417，-0.296166322，-0.055762223，
-0.336635628，-0.431943542，1.38828862，-0.107762343，-0.502041508，0.094863359，-0.266547417，-0.188964816，
-2.292527804，-0.355015921，-0.569077156，0.017514077，-0.25962396，0.583479024，0.068167933，0.014009904，
1.352042569，0.177236882，-0.564620898，0.349787556，-0.707430539，-0.129763575，1.671152083，0.776163597，
-0.021446623，0.632287045，-0.008614125，0.209323915，-0.289051354，1.681336675，-0.155862607，-0.077731799，
0.82786481，-1.041519493，0.261535927，-0.412195369，1.112557672，-0.100371081，-1.144249282，0.040589822，
-0.189283934，0.434968595，1.094428192，0.290705338，-0.0138222，0.283080016，-0.782503525，-0.085722785，
0.129237661，-0.863186702，0.146550995，0.106793381，0.866076389，-0.290592628，-0.298612901，0.585978478，
0.060169353，0.871721164，0.067681259，-0.016429211，-0.691862319，-0.099393894，0.397879547，-0.0049166，
-0.686676288，-0.782562372，-0.017342403，0.106779034，-0.020327393，0.4941572，-0.184977933，0.045901852，
0.883238737，1.461385829，1.173875952，0.238537867，0.716623196，0.32768553，1.063529252，0.24979205，
1.218106841，0.438037488，-0.11434406，0.212401723，1.731100734，-0.575671103，-0.861522307，0.323921825，
-1.171227643，1.687658722，-0.670770608，-0.347428641，0.656840307，0.654451244，0.023281919，-0.413060303，
-0.236891092，-0.290287567，-2.747778165，0.143692637，-1.71000562，-0.935081269，0.692751369，
-0.411491889，0.928809898，-0.926643719，1.04862006，0.412139218，-0.193032523，-0.009653293，-0.79853627，
-0.100936572，0.309291433，0.102463117，-0.050883159，0.694099117，-0.842781601，
-0.269146402，0.441508007，1.405110928，0.033024447，-0.328643359，-0.323136664，-0.09309922，
-0.158225088，1.260303693，-0.789521669，-0.337560523，-0.291325249，-1.066279143，
-0.618890031，0.248983641，0.034759251，0.882100279，0.48269568，-0.066941817，
-0.184608038，-0.023517406，-0.672468729，0.253540536，0.191561488，-0.087120801，-0.686079288，
-1.793905963，0.036004183，1.022793487，-1.434730455，0.589690201，-0.217817872，
-0.145532668，0.401124556，0.627431273，0.224008103，0.485978832，0.082761656，
-0.662787178，0.939936177，0.179767043，1.528823088，-0.165954634，0.393731707，-0.033624367，
-0.370122603，0.936833337，-1.404577085，-0.053542798，0.249369071，0.293947579，0.548939249，-0.006808959，
-0.11640434，-0.288842871，0.870856137，0.706003843，0.278474099，-0.229181126，-0.72134986，
-0.524715596，0.113296712，-0.787005295，0.586313224，
-0.47130353，0.155244387，0.154718569，0.025373314，0.50625468，-0.025369976，0.419284333，-0.019841149，
-0.297545391，-0.131656389，-0.25774229，0.435472774，0.420723963，0.210803862，0.551788416，-0.27312702))

list(mu=1, tau=2)
2. 首先我想知道你了解这个模型吗？软件不是用来往里面乱塞东西的，尽管BUGS软件做得很简单（通常给定先验、给定数据，剩下的什么都不用管了），但你也得稍微明白一下这些代码的来龙去脉吧。比如你的模型中有个变量n，在数据中有吗？初始值中提供了一个变量tau，可是模型中它是参数吗？

27. yin说道：

谢谢，不好意思。不小心没看到前面有一个贴的回复。

28. 寇惠说道：

你好，我刚学WINBUGS软件，现在想用来做SV-N模型，但是在模型的编辑（compile）这一步时，软件总是显示array index is greater than array bound for y 不知为什么，麻烦请帮忙指导一下，真的很急，谢谢！！！

model
{
for(i in 1:N)
{y[i]~dnorm(0,p[i])
p[i]<-exp(-theta[i])
}
theta[1]~dnorm(mu,itau2)
for(j in 2:N)
{theta[j]~dnorm(theta2[j],itau2)
theta2[j]<-mu+phi*(theta[j-1]-mu)}
phi<-2*phi1-1
mu~dnorm(0,0.01)
itau2~dgamma(2.5,0.025)
phi1~dbeta(20,1.5)
}
list(N=1000,
y=c(-0.01332,0.215247,0.165196,0.268095,-0.08932,0.784046,0.194879,-0.55396,-0.06357,-0.75685,
-0.80325,-0.24164,-0.14004,-0.06396,-0.03859,0.03763,-0.41906,-0.2917,-0.01348,-0.01342,
-0.21579,0.189124,0.215324,0.778916,-0.60096,-0.97844,0.646815,-0.92657,0.341113,-0.08932,
-0.06394,-0.19079,-0.16536,0.393219,0.242184,0.064052,-0.08932,0.03867,-0.01245,-0.06368,
-0.06368,-0.01235,0.270668,0.764389,0.249,-0.42764,-0.01135,-0.45268,0.09219,0.300745,
-0.14142,-0.21945,0.275452,0.434102,-0.03683,-0.0368,0.015803,0.200342,0.016219,0.759039,
0.123898,-0.14267,-0.19551,-0.19539,-0.56526,-0.45794,0.252919,-0.74644,0.357049,-0.32588,
-0.0368,0.173698,0.306509,-0.27424,0.14849,0.361428,0.176777,0.204207,0.823436,0.072621,
-0.11633,-1.51045,-0.77916,-0.14219,0.414057,-0.61912,0.148742,0.255552,0.711109,-0.83659,
0.150263,-1.33421,-0.71908,0.093948,-0.97636,0.014628,-0.08932,-0.7889,-0.45015,-0.29492,
0.06484,0.348767,0.454509,0.248824,-1.28074,0.013718,0.323902,-0.03755,0.195911,0.248824,
-0.29755,0.197099,-0.08932,-0.11539,0.067213,0.38176,-0.06309,0.015692,0.252739,-0.03659,
-0.19475,0.042482,0.280657,0.149247,0.28292,-0.14258,0.123898,-0.06264,-0.03593,0.124525,
-0.14283,-0.16953,-0.08932,0.151481,0.286416,-0.08932,0.126021,-0.35843,-0.2236,0.018089,
-0.06245,-0.16991,-0.08932,-0.08932,0.233432,-0.11626,0.315482,0.154349,-0.17061,-0.00803,
-0.14352,-0.30583,-0.2244,0.045758,0.208499,-0.00794,-0.08932,-0.00788,0.346171,-0.55197,
-0.00783,-0.11649,0.209957,0.047011,0.40302,0.185253,0.351572,0.021206,-0.47564,-0.17191,
0.048362,-0.72111,-0.19879,0.266892,-0.28129,-0.11672,0.10259,0.130458,0.241257,0.408604,
0.188379,-0.58863,-0.03397,-0.14468,0.049125,0.104826,0.077391,-0.14492,-0.03372,0.35661,
-0.42396,0.049973,0.30174,-0.11731,0.247001,-0.03316,0.361129,0.136667,-0.1176,-0.14585,
-0.20228,0.08017,-0.54066,-0.20184,-0.03308,0.729776,-0.08932,0.280049,-0.57207,-0.96368,
-0.06124,-0.53775,0.359109,0.276515,0.193005,-0.65318,-0.65002,0.162603,0.388274,-0.03298,
-0.08932,0.136348,0.30683,0.393839,-0.37382,0.166687,0.195904,0.196719,-0.17522,0.082549,
-0.11799,-0.03198,-0.06064,-0.03193,-0.08932,0.313325,-0.49197,0.025555,0.198448,-0.29085,
-0.11808,0.227453,0.141692,-0.11823,-0.06042,-0.06041,-0.0604,-0.06039,-0.205,-0.69442,
-0.06059,0.025687,0.169932,0.112786,-0.0315,0.113312,-0.58072,-0.37726,0.314021,-0.37759,
0.343392,0.055333,-0.00243,0.317183,0.027126,0.290074,-0.08932,-0.00156,-0.14784,-0.08932,
-0.49796,-0.08932,-0.14756,0.582508,-0.84845,0.318723,0.320395,0.410459,0.028637,-0.20728,
0.590852,0.059156,-0.26747,0.237524,-0.2083,-0.14876,-0.23776,0.029408,-0.26736,0.118424,
-0.6522,-0.26642,-0.05983,-0.35446,-0.14815,0.14618,-0.05985,-0.53056,-0.11867,-0.26522,
0.115927,-0.47016,-0.23541,-0.38087,-0.14753,-0.69846,-0.08932,0.287326,-0.03125,-0.61076,
0.403076,0.844174,-0.26502,-0.26471,-0.78781,-0.08932,-0.35,-0.11824,-0.0604,-0.40701,-0.40601,
-0.06057,-0.4051,-0.2039,-0.3752,0.196556,-0.11795,-0.2609,-0.43159,-0.06084,-0.2885,
-0.62795,0.705459,0.110364,-0.71557,-0.06094,-0.1177,-0.06094,-0.08932,-0.74005,-0.62371,
-0.11737,-0.14539,-0.48093,-0.39594,0.049934,-0.42321,-0.20037,-0.36641,-0.06165,-0.22762,
0.27064,-0.25562,-0.06163,0.299244,0.300759,-0.14514,-0.5903,-0.56017,0.104289,-0.117,
-0.3105,0.353524,-0.00607,0.523323,-0.28466,0.357731,-0.00528,-0.42508,-0.14517,-0.86797,
0.021543,0.494743,0.217984,-0.31291,0.386403,-1.03852,0.691714,-0.64783,-0.39518,-0.61542,
-0.42017,-0.33674,-1.26308,-0.1707,-0.22481,-0.5486,-0.19708,-0.00851,-0.06237,-0.51966,
0.314067,-0.00845,-0.70769,-0.33025,-0.32967,0.204521,-1.10078,-0.1952,-0.32713,-0.32657,
0.438662,-0.27443,0.069324-0.14223,0.680627,-0.91215,-0.2478,0.572669,0.791271,0.071622,
-0.43771,-0.11607,1.311186,0.673001,-0.19858,-0.17119,-0.55197,-0.44167,-0.33253,0.316358,
-1.35499,0.151803,0.367729,-1.02804,0.28511,1.015364,0.45401,-0.84916,-1.45865,-1.01833,
0.26577,-1.961,-0.27197,0.171707,-1.85103,0.270853,1.783782,2.16811,-3.62531,-4.12899,
7.486337,-10.1731,-1.89421,0.125552,-2.59122,-1.0174,-0.02001,-0.57347,-0.73123,1.732474,
0.775472,0.145696,1.763023,-0.56756,-1.7923,-2.29282,0.324472,-1.12061,1.241889,1.823496,
-0.86352,0.944279,-0.84209,-0.27663,4.87394,-2.0368,-2.11779,-0.77208,0.216165,-3.44495,
-0.0438,1.401702,3.390926,-4.83265,-2.81318,0.378125,-1.24284,-1.4475,-0.5452,0.75897,
-0.3511,0.281736,1.120681,-1.08041,0.724853,0.420149,0.601489,1.148178,-0.08932,-0.08932,
-0.38321,0.772166,1.055641,0.072024,1.327823,0.992787,1.100318,1.308628,-0.62188,-0.90683,
1.06676,-2.10378,-0.82262,0.004995,1.264457,0.678805,0.00711,1.539893,0.64887,0.753937,
-0.98196,2.23317,-0.79426,-1.26155,-0.68265,0.826866,-2.78887,-0.83712,0.223592,-0.47431,
-0.85488,-0.42242,2.023095,-0.74221,-0.73797,-0.9002,-2.96907,0.095478,-0.6427,-1.95726,
-1.78997,-1.01689,1.193881,0.670628,0.789557,0.114598,-1.97659,0.804938,0.496265,-2.01332,
0.332573,-2.05039,-0.95833,0.344236,-0.58775,-1.24988,0.747503,-0.32606,-0.23968,1.033073,
-0.08932,2.394049,1.571026,-1.86087,-0.46647,0.132358,-0.5101,0.887812,-0.37902,-0.84304,
1.491123,0.834836,1.785019,1.164467,-2.89992,1.260041,-1.4614,0.754565,0.139983,-0.8896,
1.33272,1.282952,0.780965,-1.28694,0.825211,0.264661,1.411242,1.312344,1.20892,0.976496,
4.031302,-1.66091,-1.46025,0.492296,1.366828,-1.34238,-1.42722,1.121591,-1.75057,0.185644,
0.714905,3.219945,3.360211,-2.09393,1.591838,1.483805,0.513915,0.048286,1.409041,0.190633,
2.466183,2.828889,0.029197,0.148137,1.80064,-1.89009,0.05951,0.448314,0.210618,-1.19464,
2.133659,2.526614,-3.58232,-1.34607,0.626885,-1.30976,-0.08932,3.556209,-2.93283,2.020461,
1.074907,-2.43449,-2.32189,0.73835,-1.47475,0.057155,1.446897,0.388148,-0.80468,1.135854,
1.76174,1.765445,0.505266,2.196492,0.943612,2.474078,-4.93029,-1.84982,2.084334,-4.45138,
1.632738,1.505209,0.924307,-0.78729,0.863671,0.358535,0.489458,-0.34697,1.4993,1.69114,
-2.3913,0.628462,-0.15479,0.074425,-0.1876,2.162429,1.225378,-0.83307,0.451037,-0.19086,
-0.22455,2.096557,2.074811,-2.42596,0.014144,-0.26171,1.823447,0.121574,1.043033,0.696677,
0.378053,0.70663,0.092463,0.567855,0.756738,0.838501,1.375493,-0.08932,-0.20277,0.289323,
0.138555,-0.05129,0.981434,-0.24299,0.179754,0.373642,-0.43675,1.503254,-1.2958,1.588112,
-1.65063,0.376517,-0.16711,-0.67083,0.881748,2.560625,1.935039,-1.55098,1.372336,-2.31388,
1.808575,-8.60344,-0.98385,-0.86552,2.031688,2.039227,1.186551,3.55609,-0.89318,0.472703,
1.615957,-0.00739,-1.99723,-4.0315,0.453315,-6.23177,-0.78137,0.895548,2.134639,-1.90924,
1.655664,-0.94706,2.430082,2.300002,-0.24521,-2.20845,2.772436,3.83118,-1.42673,-0.77137,
1.807809,2.7627,0.79508,-1.76721,0.620345,-0.08932,2.327586,-0.04639,1.993183,0.394089,
2.137942,-0.71797,-1.24643,0.220893,-0.26671,2.924173,1.661068,1.408833,-0.04214,0.668616,
3.690198,-2.52847,-0.18566,-0.66542,0.390525,1.608512,0.745858,2.966263,1.345129,-0.14091,
-0.34685,0.013612,0.065277,0.947494,0.223831,1.703975,1.304055,-0.08932,2.479827,1.19239,
-1.37103,-0.80657,1.573688,-0.70231,-0.03375,0.188997,0.413613,0.078887,0.022974,0.756945,
-1.16006,-1.42564,0.632302,0.862266,1.612364,-1.67846,-0.033,0.079835,1.33156,-0.03206,
-1.90549,0.361637,-1.21293,-0.53525,-0.42247,1.756884,0.704333,0.882948,0.371508,-0.89439,
-0.31816,-0.37463,0.539431,0.082844,0.428963,0.60593,0.435308,0.203333,-0.90661,0.143507,
0.027296,0.085861,-0.5558,0.027092,0.085554,1.025386,-1.14567,0.0275,0.261966,0.912749,
2.06645,-0.33116,0.273655,-0.33145,0.943808,-0.02822,-0.15043,-0.57683,0.94949,0.89834,
-2.1158,0.459293,0.955234,-0.45924,1.398616,-0.5256,-1.81563,0.400875,-0.76273,-0.75822,
-2.60292,-0.73734,-0.90804,-0.55417,0.901074,0.556789,0.62032,0.148349,-2.26652,-0.26379,
0.20163,0.553767,-3.49121,-0.08932,-0.03262,-1.04897,-0.53775,1.432263,-0.65557,0.023672,
-0.25877,-0.14574,0.874131,-0.48716,0.024185,0.708858,0.139904,-0.43296,-0.48875,-0.26002,
1.399824,1.305239,0.262378,-0.90804,0.085554,1.202534,-0.56102,0.382377,-4.08638,-1.7227,
-0.31254,-1.58313,1.85141,-0.42471,2.453906,2.637856,1.272855,0.388862,0.210708,-2.99151,
-2.22628,2.457112,-0.38198,-1.42461,-2.36992,4.171541,3.564873,0.645525,0.527203,1.846424,
2.659333,1.412502,-0.08932,0.108241,0.43943,1.851187,1.820695,0.393937,3.002983,0.555149,
0.559329,1.73488,-2.70575,-0.51881,0.411932,0.414457,-0.01715,0.999618,1.455394,0.505475,
1.716596,1.133,-0.39631,1.14431,0.689498,0.61684,-0.32526,0.14662,0.860052,0.548639,
0.391834,-0.33019,0.071192,0.636195,-0.00838,0.723691,0.074077,-0.08932,-0.49732,-1.14239,
-0.16987,-0.00877,-0.57164,0.473607,0.395723,0.316687,0.809744,-0.00719,-0.25353,0.734406,
-0.00657,-0.08932,-1.64993,0.07381,0.565881,0.487523,1.074448,-1.17041,-1.89265,-0.97894,
1.861959,1.566345,0.412771,-0.00539,0.923345,0.676965,1.374329,0.084289,-1.29832,-0.60302,
0.338578,-0.08932,0.426586,0.515954,1.39607,1.329141,1.078181,-0.08932,-1.07809,1.715781,
0.092993,-0.08932,0.367926,0.831495,0.374286,-0.46038,0.281736,-0.18222,0.562753,-0.27606,
0.847892,-0.08932,-0.37141,1.043835,-0.18424,0.38619,-0.18461,-1.41384,-0.37088,-0.08932))

list(mu=0,itau1=0.02,phi1=0.975)

1. jindan514说道：

您好，请问你这个问题解决了吗？我最近在做用SV预测波动率的实证，遇到了同样的问题，难道theta[j]需要初始化吗？

29. 寇惠说道：

当我运行compile这一步是，软件界面上说educational version cannot do this model，是什么意思啊！有知道的麻烦告诉一声，有急用！谢谢

30. yinyin148说道：

请问用openbugs时，前面模型都没有问题，但是在load inits这步，总显示unable to generate initial values for node [011483COH] of type Grapht.Mixing,是什么意思？难道初值设定有问题吗？谢谢

1. yang说道：

请问你的问题解决了吗？我也碰到一样的问题

31. yanxueqing说道：

请教大家一下，如果炫耀载入的数据非常多，如何载入那么多数据啊，有没有简便方法？

1. yanxueqing说道：

更正一下，是需要载入的数据多

32. syljbuaa说道：

您好，winbugs软件有详细的说明吗

33. Ron说道：

等WinBUGS的update时，顺便把评论都看了一遍，觉得很多朋友对Bayes分析还是很有兴趣的，但椰丝儿们留错了地方，cos论坛的椰子板块会更靠谱。看到后来发现很多评论挂着一大长串代码的都在问一件事情：SV模型怎么做？这让我觉得很惊讶，不过想来也合理。

国内的金融方面最近几年一直很热，学生极多。很多学校会对金融的学生要求论文，而高质量的也是学生保研，各种评优的重要依据。最快的方法，就是用各种统计和计量方法来让自己的文章升级（记得有个神文写的是如何拿计量经济五大杀器力克CSSCI）。我接触的一些研究生，甚至是本科生，都在搞一些很高级的处理方法和模型（高级的意思是读一两本书都搞不清楚）。金融方面比较好用的，无非就是波动模型，因此在GARCH模型已经搞烂掉的现在，SV模型凭借它”拒人于千里之外“的难理解的特点越来越受到学生的青睐。这里的MCMC之所以被提如此之多，就是因为目前SV估计比较好的方法，就是它了。

有时候一想起这些就会有些郁闷，但觉得孩子们本身没错儿，靠自己的努力和智力，迅速学一些别人学不会的方法来争取自己的前程是很棒的。令我郁闷的是这种特别的规则而引发的浮躁心态以及对思考问题的态度。就拿SV来说，有个孩子学经济的，问想拿MCMC估计SV，但理解不了那个超长的函数，我就说这个叫做后验概率密度，他说什么叫后验密度？我说是通过贝叶斯原理得到的，他说贝叶斯不是统计学的吗？我要估计的是SV模型……尽最快速度达到目标，却忽略了积累与思索，这便得不偿失了。如果有负责的老师在这方面加以引导和帮助，我想这便是孩子们在大学，甚至是近十年的一大幸事。可惜的是这可遇而不可求，因为老师同样面临和孩子们一样的境况，只不过从保研变成了评职称而已，并且他们只能靠自己。

如果真的不想要培养研究人才，拜托请不要变相用这种方式折磨孩子们的思维了。至少在绝大多数工作中，老板不会让你没事就去研究MCMC……而在看过评论之后唯一欣慰的是，COS的编辑们很热心负责的，他们抱着对统计认真的态度，在帮助在这条路上迷乱的人，尽管这很有限，但也让很多人满载而归了。

深夜吐槽 尽请见谅~那我也来做点有用的吧。SV不多说，R上面也有包可以做，MCMC方法cos上也有文章介绍。如果你确实初学，没接触R也不想学非要用WinBUGS，这里有个paper
http://www.mysmu.edu/faculty/yujun/Research/YuEJ2000.pdf
我觉得这个东西，结合上面文章的内容，绝大多SV问题就应该解决了。另外有其他问题可以到cos论坛金融区去请教版主，版主在这方面可是很强的。

1. 同惊叹！
不知道Ron大侠是否愿意写一篇科普文投给COS主站呢？

34. 为什么一运算就死机了 ，winbgus能载入多少数据啊

35. 我载入了28080个数据 程序其他的都没问题 ，但是运行不出来结果 这有什么办法解决吗？

1. 说道：

请问winbugs如何像你一样载入大量数据？

36. 说道：

请问一个问题，我在看各种winBUGS的例题时，一般标准差都会有这么一句sigma<-sqrt(1/tau),其中tau是方差，但是不应该是sigma<-sqrt(tau)吗？

37. yang说道：

请问谢老师，用openbugs时，在load inits这步，总显示unable to generate initial values for node要怎么解决？

38. ranyouhua说道：

MODEL:
Y=a*T1*P1*W*E*P2/(1+P2/b)+c*T2+d

Input data：
T1= [0.9902 0.9696 0.9851 0.9341 0.9787 0.9502 0.8941 0.9987 0.9987 0.9897 0.9946 0.7887]
T2= [18.261 22.945 22.074 24.285 22.474 23.742 25.379 20.631 19.377 21.735 18.711 11.493]
P1= [0.5961 0.5777 0.5989 0.6213 0.6258 0.6478 0.6416 0.6554 0.6514 0.6553 1.0000 1.0000]
P2= [19265.5 18877 18052.5 15476 20851.5 21657 22012 6481 21738 20884 17597.5 18969]
W= [0.8755 0.8492 0.8791 0.9119 0.9176 0.9511 0.9417 0.9625 0.9552 0.9615 0.9525 0.9180]
E= [0.2927 0.2817 0.3708 0.4006 0.3684 0.4377 0.4311 0.4939 0.4847 0.4762 0.4602 0.4043]

Yo= [-6.23 -7.24 2.21 6.1 15.18 16.38 13.63 2.53 10.51 4.08 0.73 3.51]

Cost function:
J(p) = 1/2*[(Yo-Y(p))’Co ^(-1) (Yo-Y(p))+(p-pb) ‘ Pb^(-1) (p-pb)]
where Yo represents the data vector, Y(p) is the model output vector, Co and Pb the error covariance matrix on data and model parameters, respectively. p is the parameter vector. pb is a priori values of p.
pb = [0.06, 1600, 0.05, 1];

My question is how to estimate the parameter (a, b, c, d) using WinBUGS?