请选择 进入手机版 | 继续访问电脑版
楼主: shatian
16916 71

[问答] SVCJ程序 求WINBUGS程序 [推广有奖]

ccw9796 发表于 2012-11-27 08:12:40 |显示全部楼层 |坛友微信交流群
對不住..該把第四行改成theta <- kt / k

使用道具

epoh 发表于 2012-11-27 09:12:16 |显示全部楼层 |坛友微信交流群
ccw9796 发表于 2012-11-27 08:12
對不住..該把第四行改成theta
你下载错了,不是16楼的程序,
16楼的程序是楼主依自己模型自己写的
还在修改中
你应该下载:新建 Microsoft Office Word 文档.doc

你既然要做这个模型,
应该尝试自己把公式写出来,
然后依公式写程序.

使用道具

sun999bd 发表于 2012-11-28 13:59:03 |显示全部楼层 |坛友微信交流群
Sharing is caring!!!

使用道具

ccw9796 发表于 2012-12-2 17:37:48 |显示全部楼层 |坛友微信交流群
老師您好啊..前兩天有事無法上線.今日就把您所建議的部分做一修正.並加諸風險補償觀念進來..但仍無法達陣成功--expect a comma..所以附上請您多指正
model;
{      
   mu ~ dnorm(1,0.04)
   for( i in 2 : n ) {
      v[i] ~ dt(vmean[i],ivd[i])I(0,);
   }
   kt ~ dnorm( 0.0, 1.0);
   k ~ dnorm( 0.0, 1.0);
  theta <- kt / k
   for( i in 2 : n ) {
      ksy[i] ~ dt(ksymean[i],tauy);
   }
   for( i in 2 : n ) {
      ksv[i] ~ dexp(muv);
   }
   rhoj ~ dnorm( 0.0,4);
   muv ~ dgamma(20,10);
   for( i in 2 : n ) {
      j[i] ~ dbern(lambda);
   }
   lambda ~ dbeta(2,40);
   ev ~ dnorm( 0.0, 1.0);
   ey ~ dnorm( 0.0, 1.0);
   e ~ dnorm( 0.0,1.0E-6);
   sigy2 <- 1 / tauy;
   tauy ~ dgamma( 5.0,20);
   muy ~ dnorm( 0.0,0.01);
   for( i in 2 : n ) {
    ymean[i] <-  y[i - 1] +( kt-k* y[i - 1]) +ksv[i] * j[i] +sqrt(sigv2)* rho * (y[i] - y[i - 1] -d*exp(theta[i]) - ksy[i] * j[i])
        yisig2[i] <- 1 / (rho*rho*sigy2*j[i]+(1-rho*rho)*sigv2*j[i]+v[i - 1]+1/pow(muv,2)*j[i]);
        y[i]~dt(ymean[i],yisig2[i],);
         r[i] <-  y[i]-y[i-1] ;
         rmean[i] <-  mu + ksy[i] * j[i];
      r[i] ~ dnorm(rmean[i],ird[i]);
   }
   tauv ~ dgamma( 2.5, 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])
   }
   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]);
       
   }
   for( i in 2 : n ) {
      ird[i] <- 1 / (sigy2*j[i]+v[i - 1]);
   }
   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;
}
#data
list (n=180)
n = length(y)
ksy0=rep(0,n)
ksv0=rep(0,n)
J0=rep(0,n)
V0=rep(0,n)


Y[]
-0.560.08-0.95-0.60.261.071.48-3.881.19-0.94-0.20.781.74-1.72.16-0.35-0.46-0.21-0.6-0.4-0.032.171.32-1.78-1.721.69-2.560.370.6-1.340.021.96-0.55-0.18-1.31.040.360.410.20.120.35-0.230.09-1.171.34-0.61.24-0.610.71-3.381.45-0.01-0.640.070.250.35-0.971.49-2.11-0.5403.1-1.410.2-0.470.350.32-0.69-0.49-0.931.141.320.33-2.611.340.080.42-0.87-0.430.40.360.17-0.420.420.060.630.260.05-0.030.821.6-0.790.24-0.4-0.860.09-1.131.450.36-0.660.660.080.021.17-0.41-0.42-0.24-0.290.46-1.69-0.570.820.360.14-0.94-1.35-0.670.041.430.43-0.321.4-0.9-0.17-0.70.15-0.461.941.52.22-0.53-1.47-0.390.920.08-0.06-0.171.260.84-1.13-1.58-0.71-0.45-0.670.21-2.811.18-0.310.37-1.89-0.351.51-0.06-1.010.271.370.512.08-1.080.08-0.580.430.12-1.780.760.270.96-0.28-0.140.230.07-0.090.350.29-0.620.010.03-0.11-0.230.99


END
#init



list(mu=1,k=1,kt=1,tauv=1,rho=1,lambda=1,muy=1,tauy=1,muv=1,rhoj=1, theta=1,omega=1 )
ksy=c(0,0,0,0,0,0,0,0,0,0,0,0
),
ksv=c(0,0,0,0,0,0,0,0,0,0,0,0
),
J=c(0,0,0,0,0,0,0,0,0,0,0,0
),
v = c

使用道具

ccw9796 发表于 2012-12-3 14:39:11 |显示全部楼层 |坛友微信交流群
epoh 发表于 2012-11-27 09:12
你下载错了,不是16楼的程序,
16楼的程序是楼主依自己模型自己写的
还在修改中
老師您好.請您有空幫學生修正svcj的模型.已照您的建議模型bugs去做修改.但仍無法成功達陣.萬請你不吝指正..謝謝

使用道具

shatian 发表于 2012-12-5 13:09:03 |显示全部楼层 |坛友微信交流群
epoh 发表于 2012-11-27 09:12
你下载错了,不是16楼的程序,
16楼的程序是楼主依自己模型自己写的
还在修改中
老师您好,又有问题要打搅您了,我现在用的是OPENBUGS软件,里面有个DIC功能,但我运行不起来,能点出来,但是node中输入参数名称等,都无法set,想请教下这个该怎么解决

使用道具

epoh 发表于 2012-12-5 13:38:35 |显示全部楼层 |坛友微信交流群
shatian 发表于 2012-12-5 13:09
老师您好,又有问题要打搅您了,我现在用的是OPENBUGS软件,里面有个DIC功能,但我运行不起来,能点出来, ...
你一定要照我15楼的方式作,
不仅可以取得DIC,
且模型后的分析:
normal probability plots,Jump Sizes,Jump Probability
更必须要储存上需要的数据.
已有 1 人评分学术水平 热心指数 信用等级 收起 理由
ywh19860616 + 5 + 5 + 5 谢谢,学习了

总评分: 学术水平 + 5  热心指数 + 5  信用等级 + 5   查看全部评分

使用道具

shatian 发表于 2012-12-5 16:33:47 |显示全部楼层 |坛友微信交流群
epoh 发表于 2012-11-16 10:42
1.在c碟,新建文件夹"Bugs",
  作为working.directory
  放进两文件 svcjl.bug,(your code)
老师,我还是弄不对,我想问问*.bug文件是怎么产生的,就是你说的在Bugs文件夹中要放入一个svcj1.bug文件,我生成不了

使用道具

epoh 发表于 2012-12-5 16:48:38 |显示全部楼层 |坛友微信交流群
shatian 发表于 2012-12-5 16:33
老师,我还是弄不对,我想问问*.bug文件是怎么产生的,就是你说的在Bugs文件夹中要放入一个svcj1.bug文件 ...
就是将你的svcj winbugs code 存成svcj.bug
譬如底下的code.
如果你不会存成svcj.bug, 存成svcj.txt也可以
{
   mu ~ dnorm(1,0.04)
   kt ~ dnorm( 0.0, 1.0)
   k ~ dnorm( 0.0, 1.0)
   theta <- kt / k
   rho ~ dunif(-1,1)
   rhoj ~ dnorm( 0,0.25)
   muv ~ dgamma(20,10)
   tauy ~ dgamma( 5.0,20)
   sigy2<- 1/tauy
   sigy<-sqrt(sigy2)
   tauv ~ dgamma( 2.5, 0.1)
   sigv2 <- 1 / tauv
   sigv<-sqrt(sigv2)
   muy ~ dnorm( 0.0,0.01)
   lambda ~ dbeta(2,40)
   ...
   ...
  
}

使用道具

shatian 发表于 2012-12-5 18:48:06 |显示全部楼层 |坛友微信交流群
epoh 发表于 2012-12-5 16:48
就是将你的svcj winbugs code 存成svcj.bug
譬如底下的code.
如果你不会存成svcj.bug, 存成svcj.txt也可 ...
老师,按你说的方法做了,但是我的这个程序winbugs好像运行不成功,我当时是用OPENBUGS产生的,请问R怎么调用Openbugs呢?

使用道具

您需要登录后才可以回帖 登录 | 我要注册

本版微信群
加好友,备注cda
拉您进交流群

京ICP备16021002-2号 京B2-20170662号 京公网安备 11010802022788号 论坛法律顾问:王进律师 知识产权保护声明   免责及隐私声明

GMT+8, 2024-3-29 20:46