首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >PyMC中的负二项式混合物

PyMC中的负二项式混合物
EN

Stack Overflow用户
提问于 2014-06-29 10:40:07
回答 1查看 1.5K关注 0票数 1

我正试图在PyMC中加入一种负二项式混合物。我似乎做错了什么,因为预测看上去一点也不像输入数据。问题可能是负二项式参数的先验问题。有什么建议吗?

代码语言:javascript
复制
    from sklearn.cluster import KMeans
    import pymc as mc
    n = 3 #Number of components of the mixture
    ndata = len(data)

    dd = mc.Dirichlet('dd', theta=(1,)*n)
    category = mc.Categorical('category', p=dd, size=ndata)

    kme = KMeans(n) # This is not needed but it is to help convergence
    kme.fit(data[:,newaxis])
    alphas = mc.TruncatedNormal('alphas', kme.cluster_centers_[:,0], 0.1, a=0. ,b=100000 ,size=n)
    means = mc.TruncatedNormal('means', kme.cluster_centers_[:,0],0.1,a=0.0 ,b=100000, size=n)

    @mc.deterministic
    def mean(category=category, means=means):
        return means[category]

    @mc.deterministic
    def alpha(category=category, alphas=alphas):
        return alphas[category]

    obs = mc.NegativeBinomial('obs', mean, alpha, value=data, observed = True)

    predictive = mc.NegativeBinomial('predictive', mean, alpha)

    model = mc.Model({'dd': dd,
                  'category': category,
                  'alphas': alphas,
                  'means': means,
                  'predictive':predictive,
                  'obs': obs})

    mcmc = mc.MCMC( model )
    mcmc.sample( iter=n_samples, burn=int(n_samples*0.7))
EN

回答 1

Stack Overflow用户

回答已采纳

发布于 2014-07-17 21:44:05

您正确地实现了三种分布的混合贝叶斯估计,但是MCMC模型给出了错误的值。

问题是category没有足够快地收敛,而meansalphasdd中的参数在category决定哪个点属于哪个分布之前就偏离了好的值。

代码语言:javascript
复制
data = np.atleast_2d(list(mc.rnegative_binomial(100., 10., size=s)) +
    list(mc.rnegative_binomial(200., 1000., size=s)) +
    list(mc.rnegative_binomial(300., 1000., size=s))).T
nsamples = 10000

通过可视化category,您可以看到它的后缀是错误的:

代码语言:javascript
复制
G = [data[np.nonzero(np.round(mcmc.trace("category")[:].mean(axis=0)) == i)]
    for i in range(0,3) ]
plt.hist(G, bins=30, stacked = True)

期望最大化是稳定潜在变量的经典方法,但您也可以使用快速和脏的k-均值的结果来为MCMC提供初始值:

代码语言:javascript
复制
category = mc.Categorical('category', p=dd, size=ndata, value=kme.labels_)

然后,估计值收敛到看起来合理的值。

对于alpha的上级,您只需对所有这些都使用相同的发行版:

代码语言:javascript
复制
alphas = mc.Gamma('alphas', alpha=1, beta=.0001 ,size=n)

这个问题并不是关于负二项分布的;正态分布的Dirichlet-混合物也是这样失败的;它是由于具有高维范畴分布而导致的,而MCMC在优化方面是不有效的。

票数 3
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/24475293

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档