有限混合模型|Finite mixture models ----——
有限混合模型也叫做潜在分类类模型 （laten class models）或无监督学习模型 （unsupervised learning models）,在有限数量的潜在类别中提供了异质性的自然表示， 有限混合模型 适用于建立高灵活性的概率模型，能更好地解释 随机效应 、误差分布 和模型其他方面的参数选择的真实不确定性，可广泛用于单变量或多变量密度估计、分类和非参数回归 。
当在不同条件下测量随机变量会出现混合分布 。例如，一组成年人群中的身高分布实际是男性和女性的混合；精神分裂症患者在注意力任务上的反应时间可能是他们受到或未受到注意力延迟影响的一系列试验的混合。为了实施分层建模思想，使用样本构造拟合真实概率模型的一般原则，可以将这些分布构造为简单形式的混合。例如，将男性和女性身高作为单独的单变量分布建模——在此类问题中应用 混合模型（Mixture models） ，抽样单元的总体由若干子总体组成，每个子总体中应用一个相对简单的模型。
建立有限混合模型 的基本原则是引入未观测的指标——随机变量，通常将其标记为向量或矩阵$z$，指定每个特定观测所来自的混合成分，假设观测数据属于未观测总体的类，使用概率密度或回归模型的混合分布来建模拟合，拟合模型后，可以对每个观测值的类隶属度概率进行预测。
1.1 连续混合 Continuous mixture 有限混合 是连续混合 的一个特例，$p(y_i)=\int p(y_i|\theta)\lambda(\theta d\theta)$，每个可观测的 $y_i$ 是一个随机变量，其分布取决于 $\theta_i$ ， $\theta_i$ 的先验分布或总体分布由混合分布 $\lambda(\theta)$ 确定。连续混合 的计算方法与有限混合 的计算方法非常相似,下面给出自由度为$\nu$ 的 $t$ 分布的概率模型的建立来进一步说明:
给定位置参数 $\mu$ 、方差 $\sigma^2$ 和尺度参数 $z_i$ （类似于有限混合 中的指标变量）的可观测 $y_i$ 服从 $N(\mu,\sigma^2z_i)$ .
位置和方差是混合组成参数 θ, $z_i$ 是来自混合分布 的随机样本(符合$Inv-\chi^2$ 分布)。
对 $z_i$ 进行平均后， $y_i$ 的边际分布为$t_v(\mu,\sigma^2)$。
自由度$\nu$ （可能已知也可能未知）描述了连续混合 的混合分布 ，描述方式与多项式参数 $\lambda$ 对有限混合 的描述方式相同，$z_i$ 的后验分布可用于评估哪些观测值可能是异常值。
假设混合 的成分都来自同一参数族，具有不同的参数向量，将$y=（y_1，…，y_n）$的分布或 $y | x$ 的分布看为是$H$成分的混合（不知道每个特定观测来自哪个混合成分*(mixture component)*）。利用可以确定非混合模型的信息以简化模型，例如成人身高分布例子中的性别。对于$h=1，…，H$,$h$分量的分布$f_h（y_i |θ_h）$依赖于参数向量$θ_h$；$λ_h$表示组分 $h$ 中总体比例的参数，且 $\sum_{h=1}^{H}\lambda_h=1$ 。
$y$ 的抽样分布为:
$$
p(y_i|\theta,\lambda=\lambda_1f(y_i|\theta_1)+\lambda_2f(y_i|\theta_2)+...+\lambda_Hf(y_i|\theta_H)).\tag{1}
$$
概率为$\lambda=( \lambda_1,...,\lambda_H)$的混合分布可视为参数$\theta_H$上的离散先验分布，如果把这个先验分布（或混合分布）看作是对整个群体的 $\theta$ 的变化的表述，混合模型会类似于层次模型（hierarchical model） ——引入缺失的指标 $z_{ih}$ 后两者的相似性会进一步增强，
$$
z_{ih}=
\begin{cases}
1& \text{第$i$个单位来自第$h$个混合成分}\
0& \text{其它.}
\end{cases}
$$
给定 $\lambda$ ，向量$z_i=（z_{i1}，…，z_{iH}）$的分布为$Multin（1；\lambda_1，…，\lambda_H）$。其中，混合参数 $\lambda$ 是决定 $z$ 分布的超参数。可观测的 $y$ 和未观测指标的联合分布 $z$ 为：
$$
p(y,z|\theta,\lambda)=p(z|\lambda)p(y|z,\theta)=\prod_{i=1}^{n}\prod_{h=1}^{H}(\lambda_hf(y_i|\theta_h))^{z_ih}.\tag{2}
$$
其中，$\sum i=1$，并且混合成分的数量 $H$ 已知且固定。
混合似然 (mixture likelihood)参数的可辨识性 如果为多个模型参数选择获得相同的似然函数，那么模型中的参数将无法识别，所有的有限混合模型 在某种意义上都是不可识别的。所以要定义参数空间以消除歧义，例如，指定混合成分 的平均值为非递减顺序，或指定混合比例 $λ_H$ 为非递减。
有限混合模型 参数$(\theta,\lambda)$ 的先验分布是 $\theta$ 和 $\lambda$ 独立先验分布的乘积。如果向量$z_i=（z_{i1}，…，z_{iH}）$是参数 $\lambda$ 的多项式，则自然共轭先验分布为$Dirichlet$分布，$\lambda\sim Dirichlet(\alpha_1,...,\alpha_H)$。参数 $\alpha_h$ 的相对大小描述了 $\lambda$ 先验分布的均值，$\alpha_h$的总和是先验分布强度的度量。用 $\theta$ 表示由混合成分 中所有参数组成的向量，$\theta=（\theta_1，…，\theta_H）$。
先验分布 $p(\theta)$ 没有前提假设
混合分布的参数（如$t$ 分布中的自由度）需要超先验分布。
如果数据未表明数据中存在所有$ H $ 分量，则 $\lambda$ 的错误非信息性先验分布（对应于 $\alpha_i=0$ ）会使模型出错（如果对各成分的参数使用不适当的先验分布，更容易出现问题）。
混合成分（mixture components）的数量 有限混合模型 中包含的混合成分 $ H $ 的数量通常是不确定的。计算$ H $ 值较大的模型比较麻烦，所以考虑从一个数量较小的模型开始，并评估拟合的充分性，然后确定数据的某些特征是否未被在当前模型中。恰当的测试量应满足：其后验预测分布可用于确定当前的 成分数量 是否描述了观测数据的数据范围。可以通过 $y$ 对混合模型后验分布的平均来确定成分数量 。
假设样本中的 $i=1，…，n$ 项中的每一项属于 $H$ ，每个潜在类别在参数模型中具有一个或多个参数的不同值，设 $z_i\in{1,...,H}$ 表示 $i$ 项的子群指数，即潜在类状态。然后，$i$ 在 $z_i$ 条件下的分布为：
$$
y_i|z_i\sim f(\theta_{z_i},\phi)\tag{3}
$$
$f(\theta,\phi)$ 是抽样分布， $\phi$ 在不同类中不变化，$\theta_h$ 对应于子类 $h(h=1，…，H)$ 。
假设样本属于 $h$ 类的概率为 $Pr(z_i=h)=\pi_h$ ，在边缘化潜在类状态时获得以下似然：
$$
g(y|\pi,\theta,\phi)=\sum_{h=1}^{H}\pi_hf(y|\theta_h,\phi),
$$
上式相当于一个具有 $H$ 个成分的有限混合 。一般来说。
对于 $f(y|\theta_h,\phi)=N(y|\theta_h,\phi^2)$ ， 有$(y_i|z_i=h)\sim N(\theta_h,\phi^2)$ 。正态假设的局限性在于分布是单峰对称的，如果允许正态分布的均值在不同的子群中变化，在边缘化潜在子群时，就可以得到一个更加灵活的正态混合模型 ：
$$
g(y|\pi,\theta,\phi)=\sum_{h=1}^{H}\pi_hN(y|\theta_h,\phi^2),\tag{1.4}
$$
任何密度都可以用足够多的正态混合精确近似，$\phi_h$ 允许正态核的方差在各个成分（分量）间变化，于是我们可以得到位置尺度 的正态混合，可以提高精度并且允许采样变异性在不同群体中可以更真实地实现特征描述。
可以使用绘图/聚类获得 混合参数 和各组分相对比例的初始估计。图形或聚类方法可用于初步确定从各种 混合成分 中得出的观察结果（对观测值进行分类可以直接获得不同混合成分的参数估计值）。但是图形法——这种粗略估计的方法忽略了指标中的不确定性，可能高估混合成分 间的差异。下面简要介绍两种估计方法，详见[3.1案例：精神分裂症患者的反应时间问题](###3.1 案例：精神分裂症患者的反应时间问题)
使用 $EM$ 算法对指标变量进行平均，估计有限混合模型 的参数更为实用。
$E$ 步骤 $(y,z)$ 联合分布充分统计的期望。在 有限混合 中，这相当于用贝叶斯方法 计算指标变量的条件期望。
$M$ 步骤 $argmax\theta \in \phi$ ，求 $\theta$ 的估计值。
$Gibbs$ 采样 的起始值可通过对后验的近似（位于 $t_4$ 分布的混合）进行重采样获得。对于混合模型 ，$Gibbs$ 采样交替进行两个主要步骤是：
根据给定模型参数的指标分布获取抽样；
根据给定指标的模型参数中获取抽样。
给定指标后，使用共轭族作为先验分布，混合模型 简化为普通模型（可能是分层模型）。从指标的分布中获得抽样通常很简单——因为原本就是有限混合模型 中的多项式抽样。另外，通常可以在计算的迭代模拟阶段发现建模错误（例如是否不正确地应用不适当的先验密度）。
在零方差附近开始的 $Gibbs$ 序列 可能永远不会离开该区域。如果 $Gibbs$ 序列 由于成分的排列混叠而不能收敛，可识别性问题也可能变得明显。当 $Gibbs$ 采样 达到近似收敛，忽略抽取的指标来获得关于模型参数的后验推断。指标变量的后验分布包含关于每个观测值可能被抽取的成分信息。
有限混合模型 用于对观察结果进行模式识别、调整聚类、研究异质性、生存分析，普通的单一分布难以在混杂数据中挖掘出有效信息，混合模型 则有效地解决了这个问题，已被广泛应用于各个领域，如在医学领域应用泊松混合模型 ，在工程领域应用指数混合模型 ，正态混合模型的应用则更加广泛——等方差的正态密度混合可以用来近似任意连续分布，适用于处理多模态、偏态或不对称数据。
现在用开篇所说的精神分裂症患者的反应时间 案例来进一步说明简单 有限混合模型 的应用（“简单”是指有固定数量的混合成分具有指定的分布）。
受试者共17人，其中11名为非精神分裂症患者，6名为精神分裂症患者。设计一个有关反应时间的测试，对每名受试者测试30次。测试统计数据见图1、图2 。
图 1.11名非精神分裂症患者的反应时间对数（毫秒）
图 2.6名精神分裂症患者的反应时间对数（毫秒）
从图1、图2 中可以看出，精神分裂症患者的平均反应时间较高，并且有一些患者的反应时间比其他患者的反应时间变化要大得多。相关心理学理论指出精神分裂症患者在一些试验中存在注意力缺陷和运动反射迟钝，运动迟缓影响了所有试验，注意力缺陷仅影响部分试验，因此导致患者的反应相对较慢，。
非精神分裂症患者的反应时间由正常随机效应模型描述，其中第$j$人($j=1，…，11$ ) 的反应服从正态分布，方差 $\sigma_y^2$ 相同，均值 $\alpha_j$ 不同。
有限混合似然模型 为了反映注意缺陷，每名精神分裂症患者$j$($j=12，…，17$ ) 的反应时间是包含两组分的混合：
响应没有延迟的概率为 $1-\lambda$ ，呈正态分布，均值为 $\alpha_j$ ，方差为 $\sigma_y^2$ ；
响应延迟的概率为 $\lambda$ ，观测值平均值为 $\alpha_j+\tau$ ，方差为$\sigma_y^2$。
分层模型 对精神分裂症患者与非精神分裂症患者的 $\alpha=(\alpha_1，…，\alpha_{17})$ 成分进行比较，以确定精神分裂症患者运动反射迟滞的程度。
具体而言，个体间的变异是通过使均值 $\alpha_j$ 服从正态分布来模拟的，对于非精神分裂症患者，均值为 $\mu$ ，对于精神分裂症患者，均值为 $\mu+\beta$ ，每个分布的方差为$\sigma_y^2$。即，整体分布中 $\alpha_j$ 的均值为 $\mu+\beta S_j$ ，其中 $S_j$ 是一个可观测的的指标变量，如果 $j$ 是精神分裂症患者，则为1，反之为0。
主要关注的三个参数是：$\beta$，用于测量运动反射阻滞；$\lambda$ 表示精神分裂症患者反应延迟的比例；$\tau$表示注意力转移延迟的大小。
用指标变量表示的混合模型 假设 $y_{ij}$ 是个体 $j$ 的第 $i$ 个响应，则模型可以表示为：
$$
\begin{align}
y_{ij}|\alpha_j,z_{ij},\phi\sim& \ N(\alpha_j+\tau z_{ij},\sigma_y^2),\
\alpha_j|z,\phi\sim& \ N(\mu+\beta S_j,\sigma_{\alpha}^2),\
z_{ij}|\phi\sim& \ Bernoulli(\lambda S_j),
\end{align}
$$
其中 $\phi=(\sigma_{\alpha}^2,\beta,\lambda,\tau,\mu,\sigma_y^2)$ 和 $z_{ij}$ 是未观测到的指标变量，如果对受试者 $j$ 的测量 $i$ 来自延迟分量，则为1，如果来自未延迟分量，则为0。$\theta=(\alpha,\phi)$ 表示除z指标之外的所有参数。
$z_{ij}$ 不是建模所必需的，但简化了模型中的条件分布，使我们可以使用 $ECM$ 和 $Gibbs$ 来简化计算。本例中只有两种混合成分，只需要为每一个观测值 $y_{ij}$ 提供一个指标 $z_{ij}$ 。
模型中 $H=2$ ，混合概率 $\lambda_1=\lambda,\lambda_2=1-\lambda$ ，相应的混合指标为 $z_{ij1}=z_{ij}，z_{ij2}=1-z_{ij}$ 。
超先验分布 首先在 $\phi$ 上指定一个非信息性的均匀先验密度。
在模型中，限制 $\tau$ 、方差 $\sigma_{\alpha}^2$ 和 $\sigma_y^2$ 为正。此外，由于0或1的值与混合分布不一致，所有混合成分 $\lambda$ 在[0.001,0.999]范围内是均匀的。
在本例中，$\alpha_j$ 可以通过对受试者 $j$ 的观测均值进行粗略估计， $\sigma_y^2$ 可以通过非精神分裂症患者的平均样本方差进行估计。给定 $\alpha_j$ 的估计后，将 $\alpha_j$ 分为两组（即非精神分裂症患者和精神分裂症患者）估计超参数。
用非精神分裂症患者的 $\alpha_j$ 的估计值的均值来估计 $\mu$ ，用两组之间的平均差异来估计 $\beta$ ，通过组内估计 $\alpha_j$ 的方差来估计 $\sigma_{\alpha}^2$ 。根据图2 可得粗略估计 $\hat{\lambda}=\frac{1}{3}$ ，$\hat{\tau}=1.0$ ( $z_{ij}$ 在$ECM$和$Gibbs$采样中计算）。
从 $\phi$ 的简化分布 中随机抽取100个点，并将每个点作为$ECM$最大化算法的起点。简化分布 是通过在原始参数估计中添加随机性得到的，具体而言，为了从简化分布 中得到样本，我们首先设置了上述粗略估计的参数，然后将每个参数除以一个独立的 $\chi_1^2$ 随机变量，以确保100次抽样足够分散（为了覆盖参数空间）。
E步骤 $z$ 上进行平均，给定 $\theta^{old}$ ,即表达式 $E_{old}$ 是在分布 $p(z|\theta^{old},y)$ 上平均化 $z$ 。分层混合模型 中对数后验密度* 为：
$$
E_{old}(\log p(z,\theta|y)=constant+\sum_{j=1}^{17}\log \left( N(\alpha_j|\mu+\beta S_j,\sigma_{\alpha}^2) \right)\
+\sum_{j=1}^{17} \sum_{i=1}^{30}[\log \left(N(y_{ij|}\alpha_j,\sigma_y^2) \right) \left( 1-E_{old}(z{ij}) \right) +\log \left(N(y_{ij}|\alpha_j+\tau,\sigma_y^2) \right)E_{old}(z_{ij})]\
+\sum_{j=1}^{17} \sum_{i=1}^{30}S_j[\log(1-\lambda)(1-E_{old}(z_{ij}))+\log(\lambda)E_{old}(z_{ij})].
$$
接下来还要计算每个观测值 $(i,j)$ 的 $E_{old}(z_{ij})$ ,给定 $\theta_{old}$ 和 $y$ ，条件后验概率为:
$$
Pr(z_{ij}=0|\theta_{old},y)=&1-z_{ij}\
Pr(z_{ij}=1|\theta_{old},y)=&z_{ij}
$$
其中，
$$
z_{ij}=\frac{\lambda^{old}N(y_{ij}|\alpha_j^{old}+\tau^{old},(\sigma_y^{old})^2)}{(1-\lambda^{old})N(y_{ij}|\alpha_j^{old},(\sigma_y^{old})^2)+\lambda^{old}N(y_{ij}|\alpha_j^{old}+\tau^{old},(\sigma_y^{old})^2)}.\tag{5}
$$
M步骤 $\theta$ 以增加 $E_{old}(\log p(z,\theta|y))$ 。使用 $ECM$ ，一次改变一组成分，找到条件最大值。条件最大化步骤为：
计算表现出延迟反应的患者的试验比例来更新 $\lambda$ ，将每个精神分裂症患者30次试验的贡献相加：
$$
\lambda^{new}=\frac{1}{6·30}\sum_{j=12}^{17} \sum_{i=1}^{30}z_{ij}.
$$
将 $\alpha_j$ 的正态总体分布与受试者 $j$ 上30个数据的正态混合分布相结合，根据 $\theta$ 中其他参数的当前值更新 $\alpha_j$ ：
$$
\alpha_j^{new}=\frac{\frac{1}{\sigma_{\alpha}^2}(\mu+\beta S_j)+\sum_{i=1}^{30}\sigma_y^2(y_{ij}-z_{ij}\tau)}{\frac{1}{\sigma_{\alpha}^2}+\sum_{i=1}^{30}\sigma_y^2}.\tag{6}
$$
向量 $\alpha$ ，$\tau$ 和 $\sigma_y^2$ 的更新估计值来自患者反应时间的延迟分量：
$$
\tau^{new}&=&\frac{\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}(y_{ij}-\alpha_j)}{\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}}\
(\sigma_y^{new})^2&=&\frac{1}{17·30}\sum_{j=1}^{17}\sum_{i=1}^{30}(y_{ij}-\alpha_j-z_{ij}\tau^{new})^2.
$$
向量 $\alpha$ ，参数 $\mu$ 、$\beta$ 和 $\sigma_{\alpha}^2$ 的更新估计值遵循正态分布，$ECM$ 的条件模式满足:
$$
\mu^{new}=\frac{1}{17}\sum_{j=1}^{17}(\alpha_j-\beta^{new}S_j)\
\beta^{new}=\frac{1}{6}\sum_{j=12}^{17}(\alpha_j-\mu^{new})\
(\sigma_{\alpha}^{new})^2=\frac{1}{17}\sum_{j=1}^{17}(\alpha_j-\mu^{new}-\beta^{new}S_j)^2\tag{7}
$$
上式相当于：
$$
\mu^{new}=\frac{1}{11}\sum_{j=1}^{11}\alpha_j\
\beta^{new}=\frac{1}{6}\sum_{j=12}^{17}\alpha_j-\frac{1}{11}\sum_{j=1}^{11}\alpha_j,
$$
和$(7)$中的$(\sigma_{\alpha}^{new})^2$。
从100个起始点开始，进行100次$ECM$迭代后，可以发现$(\alpha,\phi)$的三个局部极大值：一个主模式(major mode)和两个 副模式 (minor mode)。副模式基本上是互不影响的，对应于混合参数 $\lambda$ 接近零的近似退化模型（near-degenerate models） ，与主模式的后验密度比小于$e^{-20}$。由上述几点来看，可以忽略副模式 ，并且就实际情况而言，目标分布可以可以看作是单峰的。
模型参数 $\theta$ 的边际后验分布（在指标 $z_{ij}$ 上求均值）是混合形式的：
$$
p(\theta|y)\varpropto \prod_{j=1}^{17}N(\alpha|\mu +\beta S_j,\sigma_{\alpha}^2)\times\
\sum_{j=1}^{17}\sum_{i=1}^{30}((1-\lambda S_j)N(y_{ij}|\alpha_j,\sigma_y^2)+\lambda S_jN(y_{ij}|\alpha_j+\tau,\sigma_y^2)).\tag{8}
$$
在运行$ECM$算法时计算式$(8)$，确保每一步的边际后验密度增加。一旦找到 模式 ，就为 $\theta$ 构造一个多变量 $t_4$ 近似，以主模式 为中心，尺度由模式 的数值计算的二阶导数矩阵确定。
从 $t_4$ 分布中抽取$S=2000$ 个 $\theta$ 的独立样本，作为参数向量 $\theta$ 重要性重采样的起始分布。
通过重要性重采样从以主模式 为中心的 $t_4$ 近似值中抽取一组10个起点，以创建$Gibbs$ 采样的起始分布，其中，平均值接近目标均值，方差大于目标方差。
$Gibbs$ 采样很容易应用，因为完整的条件后验分布 $-p(\theta|\alpha,z,y),p(\alpha|\phi,z,y)$ 和$p(z|\alpha,\phi,y)$ 具有标准形式，可以很容易地从中取样，所需步骤类似于用于查找后验分布模式的 $ECM$ 步骤：
对于 $z_{ij}(i=1,...,30,j=12,...,17)$ ，将$z_{ij}$ 转换为为独立的$Bernoulli(z_{ij})$ ，概率 $z_{ij}$ 见式$(5)$。对于非精神分裂症患者， $z_{ij}(j<12)$ 固定为0。
对每位受试者 $j$ ，从式 $(6)$ 中均值为 $\alpha^{new}$ 的正态分布中抽取 $\alpha_j$ （注意用$Bernoulli(z_{ij})$ 替换式中的 $z_{ij}$ ）。于是，现在是以 $z$ 为条件，而不是对其进行平均。$alpha_j$的正态条件分布的方差正好是式 $(6)$ 分母的倒数。
从 $Beta(h+1，180-h+1)$ 分布中抽取混合参数 $\lambda$ ，其中 $h=\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}$ 。在180项试验中，出现注意力缺失的试验数量精神分裂症患者。模拟受 $\lambda$ 限制在区间[0.001,0.999]的约束。
其余参数按照均值和方差未知的正态分布进行处理。从给定的$(\beta,\mu,\tau,\sigma_y^2,\sigma_{\alpha}^2)$ 的后验分布中提取数据$(\alpha,\lambda,z)$，首先从方差的边际后验分布中采样，然后从其他参数的条件后验分布中采样，有：
$$
\sigma_y^2|\alpha,\lambda,z\sim Inv-\chi^2 \left( 508,\frac{1}{508}\sum_{j=1}^{17}\sum_{i=1}^{30}(y_{ij}-\alpha_j-z_{ij\tau})^2 \right),
$$
和
$$
\sigma_{\alpha}^2|\alpha,\lambda,z\sim Inv-\chi^2 \left( 15,\frac{1}{15}\sum_{j=1}^{17}(\alpha_j-\mu-\beta S_j)^2 \right).
$$
然后，从正态分布中模拟 $\tau$ :
$$
\tau |\alpha,\lambda,z,\sigma_y^2,\sigma_{\alpha}^2\sim N \left( \frac{\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}(y_{ij}-\alpha_j)}{\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}},\frac{\sigma_y^2}{\sum_{j=12}^{17}\sum_{i=1}^{30}z_{ij}} \right).
$$
$\mu$ 的分布是正态的，取决于$\alpha、\beta$ 和$\sigma^2$:
$$
\mu |\alpha,\lambda,z,\beta,\sigma_y^2,\sigma_{\alpha}^2\sim N \left( \frac{1}{17}\sum_{j=1}^{17}(\alpha_j-\beta S_j),
\frac{1}{17} \sigma_{\alpha}^2\right).
$$
$\beta$ 也具有正态条件分布：
$$
\beta |\alpha,\lambda,z,\mu,\sigma_y^2,\sigma_{\alpha}^2\sim N \left( \frac{1}{6}\sum_{j=12}^{17}(\alpha_j-\mu),
\frac{1}{6} \sigma_{\alpha}^2\right).\tag{9}
$$
正如共轭模型中常见的那样，$Gibbs$采样的步骤只是$ECM$的随机版本（例如，方差来自相关的反标度 $\chi^2$ 分布，而不是被设置为后验）。
退化点 如果存在退化点——所有 $z_{ij}$ 均为零，则不能定义条件分布 $(9)$ 的均值和方差，因为 $\tau$ 具有错误的先验分布，且在 $\sum_{ij}z_{ij}=0$ 条件下，没有延迟反应，因此不存在关于 $\tau$ 的信息，后验分布也是错误的。在精神分裂症患者的反应时间试验中，这个退化点的后验概率极低，在模拟中从未出现，如果在其他实验中，数据表明存在 $\sum_{ij}z_{ij}=0$ ，那就有必要为 $\tau$ 指定一个有信息的先验分布。**
在上述混合模型示例中，我们计算了几个单变量的估计值：17个随机效应 $\alpha_j$ 及其标准偏差 $\sigma_{\alpha}$ 、移位参数 $\tau$ 和 $\beta$ 、观测值的标准偏差 $\sigma_y$ 、混合参数$\lambda$、标准偏差比率 $\sigma_{\alpha}/\sigma_y$ 以及对数后验密度。运行10个序列的200次模拟后计算所有标量估计的潜在缩减规模 $\hat{R}$ ，发现均小于1.1,因此，模拟结果较好。
表 1.反应时间实验——新旧混合模型下相关参数的后分位数。
表1 中的前三列显示了 $Gibbs$ 取样模拟参数的后验中位数和95%区间，其中，
$\lambda$ 表示反应延迟的精神分裂症患者观察延迟的概率； $\tau$ 表示对数化 的注意力延迟； $\beta$ 表示反应未延迟的患者观察的平均对数 响应时间减去非精神分裂症患者的平均对数 响应时间；$\omega$表示有注意延迟的精神分裂症患者的比例。
利用模型参数 $\theta$ 和混合分量 $z$ 的后验分布的结果，从后验预测分布中获得抽样：
为了对受试者$j$ 进行额外测量，后验模拟$\tilde{y_{ij}}$ 应以后验模拟中的单个参数$\alpha_j$ 为条件。
对于一个新的受试者$j$ 的测量，首先应在超参数 $\mu,\beta,\sigma_{\alpha}^2$ 条件下得到一个新参数 $\tilde{\alpha_j}$ ，然后以 $\tilde{\alpha_j}$ 为条件模拟测量结果$\tilde{y_j}$ 。
不同的预测方法服务于不同的目的，为了检验模型的拟合度，可以使用第一类预测，以便将观测数据与特定个体的测量值的后验预测分布进行比较。
确定测试数量以评估拟合不良的情况 选择该模型是为了准确地拟合研究中两组受试者的均值和方差，但在对个体的拟合方面仍然存在问题——图2 精神分裂症患者的反应时间直方图表明，个体间反应时间方差存在显著差异，现检验该模型是否能够反映释该特征，计算每个患者的30个对数反应时间 $y_{ij}$ 的标准偏差 $s_j(j=12,..,17)$ ，定义$T_{min}=min(s_j)$和 $T_{max}=(s_j)$ 。
从正态混合模型中模拟预测数据集，得到 $T_{min}$ 和 $T_{max}$ 的后验预测分布，对后验分布的参数进行1000次模拟抽样。对于1000个模拟数据集中的每一个 $y^{rep}$ ，定义两个测试量$T_{min}=min(y^{rep})$和 $T_{max}=(y^{rep})$ 。
图形法比较后验预测分布和观测值 1000个模拟值的散点图如图3 所示，X $y$ 不符合后验预测分布（$T_{min}$ 太小，$T_{max}$ 太大）。
图 3.模拟值的散点图
相比之下，$s_j$ 的均值等测试量对模型j检测没有用，因为它们本质上是模型参数 $\sigma_y^2$ 的充分统计量 ，模型拟合自然较好。
为了更准确地拟合数据，在模型中增加参数：
参数 $\omega$ 表示精神分裂症患者出现注意延迟的概率
参数 $\sigma_y^2$ 表示注意延迟的方差
引入指标 $W_j$ ，如果受试者$j$ 容易出现注意延迟，则 $W_j$ 为1，否则为0
$$
y_{ij}|z_{ij},\theta\sim N(\alpha_j+\tau z_{ij},(1-z_{ij})\sigma_y^2+z_{ij}\sigma_{y2}^2)\\
\alpha_j|z,S,\mu,\beta,\sigma_{\alpha}^2 \sim N(\mu +\beta S_j,\sigma_{\alpha}^2)\\
z_{ij}|S,W,\theta \sim $ Bernoulli(\lambda S_jW_j)\\
W_j|S_j,\theta \sim $ Bernoulli(\omega S_j).
$$
拟合新模型，在$Gibbs$采样中添加三个新步骤来更新$\omega、\sigma_y^2$ 和 $W$ ，使用从之前的后验模拟中随机选择的十个抽样作为$Gibbs$采样十次平行运行的起点。在新模型中，长度为500的10个模拟序列近似收敛，所有模型参数的估计潜在规模缩减小于1.1。
表1 的最后三列显示了新模型的测试结果，显示出与旧模型的显著差异：精神分裂症患者的观察延迟的比例更大，但平均延迟更短。
使用新的后验分布下相同测试量的后验预测模拟来检查模型改进后的拟合效果，结果如图4 所示,与图3 相比，后验预测分布更接近 X
图 4.模型改进后的模拟值散点图
[1]Andrew Gelman,John Caolin.Bayesian Data Analysis, Third Edition[M]
