SVCJ程序 求WINBUGS程序

shatian 发表于 2012-11-5 21:40
老师，我看不太懂原来抓的程序，还有为什么我对mu的参数设置是不对的，我在相关论文里看到确实是这样的啊 ...

WinBUGS uses “precision” as a parameter in specifying
a Normal distribution instead of variance!!!
  precision = 1/variance
  dnorm (0, 0.01) is the same as a Normal distribution with
                            mean 0 and variance 1/0.01 = 100.

y你已经由Data输入,
#Data
list(k=5038;y=c(0.200393,0,0.901378,-0.06532,....))

v贴文时被截断了,代码如下,是vector

1. for(i in 1:k){
   
2. vmean[i+1]<-alpha+(1+beta)*v[i];
   
3. vtheta2[i+1]<-sigma2*v[i];
   
4. }

复制代码

epoh 发表于 2012-11-6 10:33
WinBUGS uses “precision” as a parameter in specifying
a Normal distribution instead of varia ...

老师，你好，什么叫做v贴文时被截断了？

epoh 发表于 2012-11-5 10:37
svj程序有些错误需要更正
1.在winbugs是 mu ~ dnorm(1, 0.04)而不是 mu ~ dnorm(0,25)
              ...

老师，我还想问问，我在一些程序里看到v ~ dnorm(vmean,ivd)I(0,)的语句，后面的I(0,)是什么意思？

shatian 发表于 2012-11-6 16:32
老师，我还想问问，我在一些程序里看到v ~ dnorm(vmean,ivd)I(0,)的语句，后面的I(0,)是什么意思？

Range constraints: manual14.pdf,page13/60
  interval restrictions coded with the I(. , .)

shatian 发表于 2012-11-6 16:32
老师，我还想问问，我在一些程序里看到v ~ dnorm(vmean,ivd)I(0,)的语句，后面的I(0,)是什么意思？

1.在c碟,新建文件夹"Bugs",
  作为working.directory
  放进两文件 svcjl.bug,(your code)
                    shatian.txt(your data)   shatian.txt (47.97 KB)

2012-11-16 10:43:12 上传

2.注意你的winbugs安装路径, 若不是"D:/WinBUGS14/"
  请依你所安装修改.
#####
3.in R 运行底下程序:
library(R2WinBUGS)
#Data
n=5038
yy=read.table(file="c:/Bugs/shatian.txt")
r=yy[,1]
summary(r)
data=list("n","r")
#initial values 请自行补上,符号应跟你的不同,请自行更正  
inits = function() {list(mu=,k=,kt=,tauv=,rho=,muy=,lambda=,tauy=,rhoj=,muv=)}
  
parameters <- c("mu","theta","k","sigv","muy","rhoj","sigy","muv","rho","lambda")
#小样先看结果
SVCJ.sim <- bugs(data, inits, parameters.to.save=parameters,"svcj.bug",n.chains=1,
               n.thin=1,n.iter=5000,n.burnin=1500,debug=TRUE,DIC=TRUE,
               bugs.directory="D:/WinBUGS14/",working.directory = "c:/Bugs/")
attach.bugs(SVCJ.sim)
print(SVCJ.sim,digits=4)


所有运行的结果都存在working.directory = "c:/Bugs"
coda1.txt
codaIndex.txt
(以上两文件工package"coda"继续分析)
data.txt(所有数据)
log.odc (所有运行结果)
log.txt
...
...

你也可以先参考我上传的范例
Three parameter logistic model
  https://bbs.pinggu.org/thread-1221672-1-1.html

epoh 发表于 2012-11-16 10:42
1.在c碟,新建文件夹"Bugs",
  作为working.directory
  放进两文件 svcjl.bug,(your code)

老师，我svcj模型出了点问题，你能帮我看看吗

shatian 发表于 2012-11-16 21:40
老师，我svcj模型出了点问题，你能帮我看看吗

你应该要参考The impact of jumps in volatility and returns
AppendixA Posterior distributions for Jumps and Volatility
定义出Xiv,Xiy,J
Xiv ~ dexp(muv)           #jump sizes in volatility
Xmean <- muy+rhoj*Xiv
Xiy ~ dnorm(Xmean,tauy)   #jump sizes in returns
J~ dbern(lambda)          #Bernoulli random variables with intensity lambda

epoh老师您好:很冒昧打扰您..我这几天参考网上的参数来试用一些资料..结果一直无法成功..一是MODEL一直出现--expect a comma....-二是load data 会出现--defined variable name...不知那里出问题..我参考网上和您推荐的文章亦想不出来.请您可否指点迷津..谢谢老师.
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 ) {
   rmean[i] <-  mu + ksymean[i] * j[i]
        rmean[i] <-  mu + ksy[i] * j[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])
  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
}


n = length(y)
ksy0=rep(0,n)
ksv0=rep(0,n)
J0=rep(0,n)
v0=rep(0,n)
#data

list (n=180,v=c(
-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, ev=1,ey=1,e=1,tauv=1,theta=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-11-25 12:39
epoh老师您好:很冒昧打扰您..我这几天参考网上的参数来试用一些资料..结果一直无法成功..一是MODEL一直出现 ...

依你的bug code,你要给的数据是y,由y算出 r
  for( i in 2 : n ) {
      r <-  y-y[i-1]
      r ~ dnorm(rmean,ird)
   }
但你没给y,反而给v=c(-0.560.08,-0.95,-0.60,.261,..),v是要求的数据

重复定义
  vmean <- v[i - 1] +( kt-k* v[i - 1]) + 1/muv * j.....
  vmean <- v[i - 1] +( kt-k* v[i - 1]) +ksv * j ....
也重复定义
  rmean <-  mu + ksymean * j
  rmean <-  mu + ksy * j

你的bog code,的确是乱凑出来的.
建议你应由SV,SVJ,SVCJ逐渐了解.

如果你只是要跑出结果,楼主网上下载提供的bug code就ok了
那个程序我看过,是可以运行的

EPOH老师您好..感谢您的指导.我将您对我的缺点做一一修正结果仍是无法运作.调整如下:
model;
{

kt ~ dnorm( 0.0, 1.0)
   k ~ dnorm( 0.0, 1.0)
theta<-alpha/k;
mu ~ dnorm(1,0.04)
   for( i in 2 : n ) {
      v[i] ~ dt(vmean[i],ivd[i])I(0,)
   
isigmav2~dgamma(2.5,0.1);
sigmav2<-1/isigmav2;
muy~dnorm(0,0.01);
isigmay2~dgamma(10,40);
sigmay2<-1/isigmay2;
isigay2 ~ dgamma( 5.0,20)
lambda ~ dbeta(2,40)
   ev ~ dnorm( 0.0, 1.0)
   ey ~ dnorm( 0.0, 1.0)
   e ~ dnorm( 0.0,1.0E-6)
rhoj~dnorm(0,2);
imuv~dgamma(10,20);
muv<-1/imuv;
lamda~dbeta(2,40);
for( i in 2 : n ) {
      ksv[i] ~ dexp(muv)
v0mean<-(k*theta+lamda/muv);
v0theta2<-(2*lamda-pow(lamda,2)/(pow(muv,2)));
v0~dt(v0mean,v0theta2);
vmean[i-1]<-(k*theta+(1-k)*v0+lamda/muv);
vtheta2[i-1]<-(sigmav2*v0+(2*lamda-pow(lamda,2))/(pow(muv,2)));
v[i-1]~dt(vmean[i-1],vtheta2[i-1]);
muystar<-muy+rhoj*muy;
ymean[i-1]<-mu+lamda*muystar;
ytheta2[i-1]<-v0+pow(muystar,2)*lamda*(1-lamda)+sigmay2*lamda;
y[i-1]~dt(ymean[i-1],ytheta2[i-1]);
for (i in 2:N){
vmean[i]<-(k*theta+(1-k)*v[i-i]+lamda/muv);
vtheta2[i]<-(sigmav2*v[i-1]+(2*lamda-pow(lamda,2))/(pow(muv,2)));
v[i]~dt(vmean[i],vtheta2[i]);}
for (i in 2:N){
ymean[i]<-mu+lamda*muystar;
ytheta2[i]<-v0+pow(muystar,2)*lamda*(1-lamda)+sigmay2*lamda;
y[i]~dnorm(ymean[i],ytheta2[i]);
}
}
#data
list (n=180)
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 )
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



还是无法成功运作..请您再帮学生指导..感谢您