在Winbugs中拟合分布混合（高斯+均匀）

gbr 发表于 Dev

GBR

我正在尝试将混合分布模型拟合为值向量，混合需要包含2个高斯分布和1个均匀分布。我正在尝试在Winbugs中实现这一点。我发现了很多使用高斯混合的例子，但无法弄清楚如何添加制服。下面的代码粘贴程序当前是参数化的，以适应范围在零和一之间的值的向量，但是我得到“节点NSD [1]的多个定义”，因此看来我的结构还是错误的。有什么建议？

model{

   ## priors
    xmin~dunif(0,1)
    eps2~dunif(0,1)
    xmax<-min(xmin+eps2, 1)
    mu1~dunif(0,1)
    eps1~dunif(0,1)
    mu2<-min(mu1+eps1,1)

   sigma1 ~ dunif(0,.5)     
   sigma2 ~ dunif(0,.5)     
   tau1<-pow(sigma1,-2)
   tau2<-pow(sigma2,-2)
   alpha[1]<-1
   alpha[2]<-1
   alpha[3]<-1
   p.state[1:3]~ddirch(alpha[])

   for (t in 1:npts) {
     idx[t] ~ dcat(p.state[])   ##  idx is the latent variable and the parameter index
     x[t,1]~dnorm(mu1,tau1)
     x[t,2]~dnorm(mu2,tau2)
     x[t,3]~dunif(xmin,xmax) 

      NSD[t] <-x[t,idx[t]]    
      }
}

鲍姆斯

您可以尝试使用无信息的dnorm先验代替dunif先验，以便可以将先验建模NSD为~ dnorm(mu[idx[t]], tau[idx[t]])。但是，您需要截断，因此当您需要普通的先验时，可以为截断设置非常低/较高的边界。

也许是这样的：

model  {
  mu[1] ~ dunif(0, 1)
  mu[2] <- min(mu[1] + eps[1], 1)
  mu[3] <- 0.5
  eps[1] ~ dunif(0, 1)
  eps[2] ~ dunif(0, 1)
  sigma[1] ~ dunif(0,.5)     
  sigma[2] ~ dunif(0,.5)     
  tau[1] <- pow(sigma[1],-2)
  tau[2] <- pow(sigma[2],-2)
  tau[3] <- 0.000001
  left[1] <- -100 # something relatively very low
  left[2] <- -100 # something relatively very low
  left[3] ~ dunif(0, 1)
  right[1] <- 100 # something relatively very high
  right[2] <- 100 # something relatively very high
  right[3] <- min(left[3] + eps[2], 1)
  alpha[1] <- 1
  alpha[2] <- 1
  alpha[3] <- 1
  p.state[1:3] ~ ddirch(alpha[])

  for (t in 1:npts) {
    idx[t] ~ dcat(p.state[])
    NSD[t] ~ dnorm(mu[idx[t]], tau[idx[t]])T(left[idx[t]], right[idx[t]])  
  }
}

截断的模糊正态分布应大致等于均匀分布。我们可以比较adnorm(0, 0.000001)T(0, 1)和a样本的核密度dunif(0, 1)。在这里，我使用R中的JAGS，但WinBUGS的结果应该类似：

library(R2jags)
M <- '
model {
  y_tnorm ~ dnorm(0, 0.000001)T(0, 1)
  y_unif ~ dunif(0, 1)
}
'
out <- jags(list(), NULL, c('y_tnorm', 'y_unif'), textConnection(M), 1, 100000, 
            n.burnin=0, n.thin=1, DIC=FALSE)

plot(density(out$BUGSoutput$sims.matrix[, 'y_tnorm'], bw=0.001), main='')
lines(density(out$BUGSoutput$sims.matrix[, 'y_unif'], bw=0.001), col=2)
legend('bottomright', c('Truncated normal', 'Uniform'), bty='n', 
       col=1:2, lty=1, inset=0.05)

在此处输入图片说明

编辑

该模型似乎在JAGS中运行良好。

M <- 'model  {
  mu[1] ~ dunif(0, 1)
  mu[2] <- min(mu[1] + eps[1], 1)
  mu[3] <- 0.5
  eps[1] ~ dunif(0, 1)
  eps[2] ~ dunif(0, 1)
  sigma[1] ~ dunif(0,.5)     
  sigma[2] ~ dunif(0,.5)     
  tau[1] <- pow(sigma[1],-2)
  tau[2] <- pow(sigma[2],-2)
  tau[3] <- 0.000001
  left[1] <- -100 # something relatively very low
  left[2] <- -100 # something relatively very low
  left[3] ~ dunif(0, 1)
  right[1] <- 100 # something relatively very high
  right[2] <- 100 # something relatively very high
  right[3] <- min(left[3] + eps[2], 1)
  alpha[1] <- 1
  alpha[2] <- 1
  alpha[3] <- 1
  p.state[1:3] ~ ddirch(alpha[])

  for (t in 1:npts) {
    idx[t] ~ dcat(p.state[])
    NSD[t] ~ dnorm(mu[idx[t]], tau[idx[t]])T(left[idx[t]], right[idx[t]])  
  }
}'


d <- read.csv('NSD.csv')

library(R2jags)
jagsfit <- jags(list(NSD=d$NSD, npts=nrow(d)), NULL, 
                c('mu', 'sigma', 'eps', 'left', 'right', 'p.state'), 
                textConnection(M), 3, 50000)

我没有让它运行足够长的时间以使所有参数完全收敛，但是这里是对一些参数的初步了解。

##                  mean      sd       2.5%        25%        50%        75%      97.5%   Rhat n.eff
## deviance   -2650.2912 16.7002 -2667.7334 -2663.5577 -2656.7462 -2639.8387 -2610.2082 1.0054   450
## eps[1]         0.9514  0.0021     0.9472     0.9500     0.9514     0.9528     0.9556 1.0018  2500
## eps[2]         0.9100  0.0523     0.8438     0.8590     0.9018     0.9569     0.9975 1.0087   260
## left[1]     -100.0000  0.0000  -100.0000  -100.0000  -100.0000  -100.0000  -100.0000 1.0000     1
## left[2]     -100.0000  0.0000  -100.0000  -100.0000  -100.0000  -100.0000  -100.0000 1.0000     1
## left[3]        0.0021  0.0013     0.0001     0.0011     0.0021     0.0032     0.0043 1.0011 14000
## mu[1]          0.0008  0.0001     0.0007     0.0008     0.0008     0.0008     0.0009 1.0011 22000
## mu[2]          0.9522  0.0021     0.9480     0.9508     0.9522     0.9536     0.9564 1.0017  2600
## mu[3]          0.5000  0.0000     0.5000     0.5000     0.5000     0.5000     0.5000 1.0000     1
## p.state[1]     0.4721  0.0259     0.4217     0.4546     0.4721     0.4898     0.5227 1.0010 60000
## p.state[2]     0.3712  0.0265     0.3193     0.3532     0.3711     0.3890     0.4234 1.0017  2900
## p.state[3]     0.1567  0.0207     0.1189     0.1423     0.1558     0.1700     0.1999 1.0019  2300
## right[1]     100.0000  0.0000   100.0000   100.0000   100.0000   100.0000   100.0000 1.0000     1
## right[2]     100.0000  0.0000   100.0000   100.0000   100.0000   100.0000   100.0000 1.0000     1
## right[3]       0.9121  0.0522     0.8465     0.8610     0.9038     0.9589     0.9997 1.0087   260
## sigma[1]       0.0007  0.0000     0.0006     0.0007     0.0007     0.0007     0.0008 1.0010 60000
## sigma[2]       0.0238  0.0016     0.0210     0.0227     0.0238     0.0248     0.0272 1.0016  3200

本文收集自互联网，转载请注明来源。

如有侵权，请联系[email protected] 删除。

编辑于2021-02-19

我来说两句

0条评论

登录后参与评论

上一篇：在Powershell中无法在此对象上找到属性“值”

来自分类Dev

Related 相关文章

文章

在Winbugs中拟合分布混合（高斯+均匀）

在Winbugs中拟合分布混合（高斯+均匀）

Flexbox之间的空间均匀，项目分布均匀

均匀分布按钮

随机均匀分布

随机均匀分布

离散均匀圆形分布

衰变的块状均匀分布

Python适合均匀分布

箱线图均匀分布

在div中均匀分布按钮

将UIImageViews均匀分布在UIScrollView中

verilog中均匀分布的随机实数

数组中数字的均匀分布

从多集合中均匀采样

在MATLAB中强制均匀输出

从均匀分布的混合物中抽取随机数

文本不会包裹在均匀分布的、均匀宽度的 flexbox 中

Randint并不总是遵循均匀分布

在圆内均匀分布x点

NumPy：均匀分布的N维样本

给定比例的随机非均匀分布

证明随机产生的数是均匀分布的

Bootstrap 3、11列均匀分布？

如何以圆形均匀分布点

Python：舍入错误扭曲了均匀分布

如何生成均匀分布的整数序列？

均匀分布的变量之间的区别，Python

均匀分布html标签菜单的按钮

选择均匀分布的点算法

如何垂直均匀地分布动态div？