0

我想比较两组多项式 (k=3) 观测值的 Dirichlet 超参数。

我可以分开装,而且它们看起来很好。

第一组 (N=10) 观察是从 [.30 .60 .10] 的 alpha 向量模拟的,计数在 50 左右变化

import numpy as np
import pymc3 as pm
import pandas as pd
import theano.tensor as tt

#model 1
N1=10;
c1 = np.asarray([[16,29,4],
                 [16,29,6],
                 [14,30,4],
                 [16,29,3],
                 [16,31,5],
                 [13,29,5],
                 [15,32,5],
                 [15,29,6],
                 [17,31,6],
                 [13,29,4]]);

counts1 = np.asarray([49,51,48,48,52,47,52,50,54,46]);

with pm.Model() as model1:
    
    hyper_param1 = pm.HalfNormal('hyper_param1',10,shape=3);
    param1 = pm.Dirichlet('param1',a=hyper_param1,shape=(N1,3));
    y1 = pm.Multinomial('y1', n=counts1, p=param1, observed=c1);
    trace1 = pm.sample(10000, tune=10000, init='adapt_diag',target_accept=0.95);

后验和轨迹看起来不错,并且使用 hyper_param1 的每个维度的后验均值的狄利克雷分布正好位于它应该在的位置:

在此处输入图像描述

第二组 (N=11) 观察是从 [.30 .10 .60] 的 alpha 向量模拟的,计数再次在 50 左右变化:

#model 2
N2 = 11;
c2 = np.asarray([[15,4,31],
                 [15,3,30],
                 [15,5,31],
                 [16,6,31],
                 [14,5,29],
                 [16,5,30],
                 [15,6,30],
                 [13,6,30],
                 [14,6,31],
                 [16,6,29],
                 [15,7,29]]);

counts2 = np.asarray([50,48,51,53,48,51,51,49,51,51,51]);

with pm.Model() as model2:
    
    hyper_param2 = pm.HalfNormal('hyper_param2',10,shape=3);
    param2 = pm.Dirichlet('param2',a=hyper_param2,shape=(N2,3));
    y2 = pm.Multinomial('y2', n=counts2, p=param2, observed=c2);
    trace2 = pm.sample(10000, tune=10000, init='adapt_diag',target_accept=0.95);

同样,后验、轨迹和 Dirichlet 分布看起来很好: 在此处输入图像描述

我尝试使用通常在 PyMC3 中工作的索引格式将这些组合成一个层次模型:

## Combined model
c_comb = np.asarray([[16,29,4],
                [16,29,6],
                [14,30,4],
                [16,29,3],
                [16,31,5],
                [13,29,5],
                [15,32,5],
                [15,29,6],
                [17,31,6],
                [13,29,4],
                [15,4,31],
                [15,3,30],
                [15,5,31],
                [16,6,31],
                [14,5,29],
                [16,5,30],
                [15,6,30],
                [13,6,30],
                [14,6,31],
                [16,6,29],
                [15,7,29]]);

counts_comb = np.asarray([49,51,48,48,52,47,52,50,54,46,50,48,51,53,48,51,51,49,51,51,51]);
idx = np.asarray([0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1]).astype(int);

with pm.Model() as model_comb:
    
    hyper_param_comb1 = pm.HalfNormal('hyper_param_comb1',10,shape=3);
    hyper_param_comb2 = pm.HalfNormal('hyper_param_comb2',10,shape=3);

    param_comb1 = pm.Dirichlet('param_comb1',a=hyper_param_comb1,shape=(N1,3));
    param_comb2 = pm.Dirichlet('param_comb2',a=hyper_param_comb2,shape=(N2,3));
    
    param_comb = tt.concatenate((param_comb1,param_comb2),axis=0);
    
    y_comb = pm.Multinomial('y_comb', n=counts_comb, p=param_comb[idx], observed=c_comb);
    trace_comb = pm.sample(10000, tune=10000, init='adapt_diag',target_accept=0.95);

该模型几乎不适合。>100 的分歧,后验/轨迹是一团糟,狄利克雷分布相当缺乏信息:

在此处输入图像描述

我尝试增加对调整的接受要求(最高 0.99),增加调整样本,hyper_params 上的各种先验(gamma,StudentsT,normal),还尝试拟合参数向量的指数以防负值是不知何故出现并造成混乱......没有任何效果。要么我在模型设置上做一些愚蠢的事情,要么这是一些需要解决的偏离中心/中心问题?

PyMC3 版本:3.8 Python 版本:3.8 操作系统:Linux (PopOS v20.10) IDE:Spyder v4.2.1

4

0 回答 0