我给楼主说下逆高斯怎么弄，winbugs里，有些分布没有预设，需要我们自己另外定义。
# --------------- INVERSE GAUSSIAN MODEL ------------------------------------
# --------------- ZEROS TRICK ------------------------------------
model{
C <- 10000
for (i in 1:n) {
zeros <- 0
zeros ~ dpois( zeros.mean)
zeros.mean <- -l + C
l <- 0.5*( log( lambda ) - log(2*3.14) - 3*log(y) ) - 0.5* lambda * pow( (y-mu)/mu, 2 )/y
log(mu) <- beta[1] + beta[2]*x1 + beta[3]*x2 + beta[4]*x3 + beta[5]*x4
}
# priors
for (j in 1:5){ beta[j] ~ dnorm( 0.0 , 0.001 ) }
lambda ~ dgamma( 0.01, 0.01)
s2 <- 1/lambda
}
INITS
list( beta=c(0,0,0,0,0), lambda=1 )
# --------------- ONES TRICK ------------------------------------
model{
C <- 10
for (i in 1:n) {
ones <- 1
ones ~ dbern( ones.p)
ones.p <- exp( l - C )
# write here the expression of the log-likelihood for i observation
l <- 0.5*( log( lambda ) - log(2*3.14) - 3*log(y) ) - 0.5* lambda * pow( (y-mu)/mu, 2 )/y
log(mu) <- beta[1] + beta[2]*x1 + beta[3]*x2 + beta[4]*x3 + beta[5]*x4
}
# priors
for (j in 1:5){ beta[j] ~ dnorm( 0.0 , 0.001 ) }
lambda ~ dgamma( 0.01, 0.01)
s2 <- 1/lambda
}
在自己的程序里这样运用如下：
#model
{
C <- 10
for (i in 1:n) {
ones <- 1
ones ~ dbern( ones.p)
ones.p <- exp( l - C )
# write here the expression of the log-likelihood for i observation
l <- 0.5*( log( lambda ) - log(2*3.14) - 3*log(y) ) - 0.5* lambda * pow( (y-mu)/mu, 2 )/y
lambda<-2
mu<-4
delta2<-y
y,就是我们要用的逆高斯分布，我这里期望值：4，delta的平方：2，lambda与delta平方实际室友关系的，但是，转化后lambda还是为2，转化没有什么意义，相当于没变，在用逆高斯分布时，你只需要改这两个值，然后引用y
扫码加好友，拉您进群
全部版块
我的主页
收藏
