如何测试/证明数据是零膨胀的？

机器算法验证 r 分布 引导程序 零通胀
2022-03-18 19:57:58

我有一个问题，我认为应该很简单，但无法弄清楚。我正在研究种子授粉，我有成簇开花的植物（n = 36），我从每株植物中取样 3 个花簇，从每个簇中抽取 6 个种子荚（每株植物总共有 18 个种子荚）。一个豆荚可以有 0 到最多 4 颗授粉的种子。因此，数据是计数的，具有上限。我发现平均约 10% 的种子被授粉，但在给定植物上的任何地方都在 1 - 30% 之间，因此在分散的数据中，当然，在 3 个植物上缺少 4 个簇重复，所以不是完全对称的.

我要问的问题是，这些数据是否支持这种植物需要传粉者进行种子集结的想法。

我发现一个豆荚中种子数量的分布看起来比 0 授粉种子豆荚（16 个中的 6-9 个豆荚）和 3 和 4 个授粉种子豆荚（每个 2-4 个）更多如果种群中的种子只是随机授粉，则可以预期。基本上，我认为这是零膨胀数据的经典示例，首先昆虫要么访问要么根本不访问花朵（一个零生成器），如果它访问，然后在另一个分布中为 0-4 个种子授粉。另一种假设是植物部分自交，然后预计每颗种子都有相同的授粉概率（该数据表明大约 0.1 的机会，这意味着同一豆荚中的两颗种子有 0.01 的机会等） .

但我只是想证明数据最适合一种或另一种分布，而不是实际上对数据进行 ZIP 或 ZINB。我认为我使用的任何方法都应该考虑到授粉种子的实际数量和每株植物上取样的豆荚数量。我想出的最好的办法是做一些引导式的事情，我只是将给定植物的授粉种子数量随机分配到我采样的种子荚数量中，做 10,000 次，看看它的可能性有多大给定植物的实验数据来自该随机分布。

我只是觉得这应该比蛮力引导要容易得多，但是经过几天的思考和搜索，我放弃了。我不能只与泊松分布进行比较，因为它是上限，它不是二项式的，因为我需要以某种方式首先生成预期分布。有什么想法吗？而且我正在使用 R，所以那里的建议（特别是如何最优雅地将 10,000 个随机分布的 n 个球生成到 16 个盒子中，每个盒子最多可以包含 4 个球）将是最受欢迎的。

添加于 2012 年 9 月 7 日 首先，感谢大家的关注和帮助。阅读答案让我想重新措辞一下我的问题。我要说的是我有一个假设（现在我认为它是零）种子在豆荚中随机授粉，而我的替代假设是具有至少 1 个授粉种子的种子豆荚更有可能有多个授粉种子比随机过程预期的要多。我提供了三个工厂的真实数据作为例子来说明我在说什么。第一列是豆荚中授粉种子的数量，第二列是具有该种子数量的豆荚的频率。

植物 1（共 3 粒种子：4% 授粉）

种子数 ::pod.freq

0::16

1::1

2::1

3::0

4::0

植物 2（共 19 粒种子：26% 授粉）

num.seeds::pod.freq

0::12

1::1

2::1

3::0

4::4

植物 3（共 16 粒种子：22% 授粉）

num.seeds::pod.freq

0::9

1::4

2::3

3::2

4::0

在植物#1 中，18 个豆荚中只有 3 颗种子授粉，一个豆荚有一颗种子，一个豆荚有两颗种子。考虑在豆荚中随机添加一颗种子的过程，前两颗种子各自进入自己的豆荚，但对于第三颗种子，已经有一颗种子的豆荚中有 6 个可用点，而 16 个豆荚中有 64 个点没有种子，所以这里有 2 个种子的豆荚的最高概率是 6/64= 0.094。这有点低，但不是很极端，所以我想说这种植物符合所有种子随机授粉的假设，授粉发生的几率约为 4%。但对我来说，植物 2 看起来要极端得多，有 4 个豆荚完全授粉，但 12 个豆荚什么都没有。我不太确定如何直接计算这种分布的几率（因此我的引导想法），但我猜如果每颗种子有约 25% 的授粉机会非常低，这种分布的几率是随机发生的。植物＃3我真的不知道，我认为随机分布的0和3比预期的要多，但我的直觉是这种种子数量的分布比植物＃2的分布更有可能，并且可能不是那么不可能。但显然我想确定，并且跨越所有工厂。我认为随机分布的 0 和 3 比预期的要多，但我的直觉是，这种种子数量的分布比植物#2 的分布更有可能，而且可能不太可能。但显然我想确定，并且跨越所有工厂。我认为随机分布的 0 和 3 比预期的要多，但我的直觉是，这种种子数量的分布比植物#2 的分布更有可能，而且可能不太可能。但显然我想确定，并且跨越所有工厂。

最后，我希望写一个声明，如“授粉种子在种子荚中的分布符合（或不符合）植物不仅仅是部分自相容的假设，而是需要访问传粉者来影响种子集。（统计检验结果）。” 这实际上只是我的前瞻性部分的一部分，我正在谈论下一步要进行哪些实验，所以我并不急于这是一回事，但如果可能的话，我想自己知道。如果我不能用这些数据做我想做的事情，我也想知道！

一开始我确实问了一个相当广泛的问题，因为我很好奇是否有任何好的测试来显示数据是否应该首先进入零膨胀模型。我看到的所有例子似乎都在说——“看，这里有很多零，对此有一个合理的解释，所以让我们使用零膨胀模型”。这就是我现在在这个论坛上正在做的事情，但我在上一章有一个经验，我使用泊松 glm 来处理计数数据，我的一位主管说“不，glms 太复杂和不必要，这个数据应该进入列联表”，然后向我发送了由他们昂贵的统计数据包生成的海量列联表的数据转储，该数据包为我的所有因素 + 交互提供了相同的 p 值到三位有效数字！所以，我试图保持统计数据清晰简单，并确保我对它们的理解足以有力地捍卫我的选择，我觉得我现在无法为零膨胀模型做到这一点。我已经使用了准二项式（用于整株植物以摆脱假复制）和上述数据的混合模型来比较处理并回答我的主要实验问题，两者似乎都做同样的工作，但我也会今晚和 ZINB 一起玩，看看它的表现如何。我在想，如果我能明确地证明这些数据首先是强烈聚集的（或零膨胀），然后为这种情况的发生提供一个很好的生物学原因，我会更好地设置随后拔出 ZINB，而不是只需将一个与准二项式/混合模型进行比较并争论，因为它给出了更好的结果，这就是我应该使用的。对于零膨胀模型，我觉得我现在做不到。我已经使用了准二项式（用于整株植物以摆脱假复制）和上述数据的混合模型来比较处理并回答我的主要实验问题，两者似乎都做同样的工作，但我也会今晚和 ZINB 一起玩，看看它的表现如何。我在想，如果我能明确地证明这些数据首先是强烈聚集的（或零膨胀），然后为这种情况的发生提供一个很好的生物学原因，我会更好地设置随后拔出 ZINB，而不是只需将一个与准二项式/混合模型进行比较并争论，因为它给出了更好的结果，这就是我应该使用的。对于零膨胀模型，我觉得我现在做不到。我已经使用了准二项式（用于整株植物以摆脱假复制）和上述数据的混合模型来比较处理并回答我的主要实验问题，两者似乎都做同样的工作，但我也会今晚和 ZINB 一起玩，看看它的表现如何。我在想，如果我能明确地证明这些数据首先是强烈聚集的（或零膨胀），然后为这种情况的发生提供一个很好的生物学原因，我会更好地设置随后拔出 ZINB，而不是只需将一个与准二项式/混合模型进行比较并争论，因为它给出了更好的结果，这就是我应该使用的。我已经使用了准二项式（用于整株植物以摆脱假复制）和上述数据的混合模型来比较处理并回答我的主要实验问题，两者似乎都做同样的工作，但我也会今晚和 ZINB 一起玩，看看它的表现如何。我在想，如果我能明确地证明这些数据首先是强烈聚集的（或零膨胀），然后为这种情况的发生提供一个很好的生物学原因，我会更好地设置随后拔出 ZINB，而不是只需将一个与准二项式/混合模型进行比较并争论，因为它给出了更好的结果，这就是我应该使用的。我已经使用了准二项式（用于整株植物以摆脱假复制）和上述数据的混合模型来比较处理并回答我的主要实验问题，两者似乎都做同样的工作，但我也会今晚和 ZINB 一起玩，看看它的表现如何。我在想，如果我能明确地证明这些数据首先是强烈聚集的（或零膨胀），然后为这种情况的发生提供一个很好的生物学原因，我会更好地设置随后拔出 ZINB，而不是只需将一个与准二项式/混合模型进行比较并争论，因为它给出了更好的结果，这就是我应该使用的。

但我不想过多地分散我的主要问题，我如何确定我的数据是否真的比随机分布的预期更零膨胀？就我而言，这个问题的答案是我真正感兴趣的，模型证明可能带来的好处是一个奖励。

再次感谢您的所有时间和帮助！

干杯， BWGIA

3个回答

对我来说，这似乎是一个相对简单（非线性）的混合模型。您将种子荚嵌套在嵌套在植物中的簇中，并且您可以在每个阶段拟合具有随机效应的二项式模型：

    library(lme4)
    binre <- lmer( pollinated ~ 1 + (1|plant) + (1|cluster), data = my.data, family = binomial)

或者如果你有协变量。如果花朵自花授粉，那么您可能会看到一些轻微的影响，因为植物自身的生存能力存在自然变异性。但是，如果响应中的大多数变异性是由集群变异性驱动的，那么您将有更强的证据表明昆虫授粉可能只访问植物上的选定集群。理想情况下，您需要随机效应的非参数分布而不是高斯分布：点质量为零，没有昆虫访问，点质量为正值——这本质上是 Michael Chernick 所考虑的混合模型。您可以将其与GLLAMM Stata 包配合使用，如果这在 R 中是不可能的，我会感到惊讶。

可能对于一个干净的实验，您可能希望将植物放在里面，或者至少在没有昆虫进入的地方，看看有多少种子会授粉。这可能会以更严格的方法回答你所有的问题。

在我看来，这是每种昆虫的混合分布。在概率 p 下，昆虫确实以 1-p 的概率着陆，它着陆并分发 0 到 4 粒种子。但是如果你没有关于昆虫是否落在植物上的信息，你就无法区分得到 0 的两种方式。所以你可以让 p 是 0 的概率，然后你就有多项分布（p1，p2， p3, p4) 其中 pi 是给定昆虫授粉受约束 p1+p2+p3+p4=1 的 i 种子的概率。该模型有五个未知数 p、p1、p2、p3、p4，每个 i 的约束条件为 0=0。有了足够的数据，您应该能够使用受限的最大似然方法来估计这些参数。

这是您问题最后一部分的答案，即如何快速生成传粉者假设所需的数据：

n = 16
max = 4
p1 = 0.1
p2 = 0.9
Y1 = rbinom(10000*n,1,p1)
Y2 = matrix(Y1*rbinom(10000*n,4,p2),ncol=16)

You can also use rzibinom() in package VGAM. Although I'm not sure what you want to do with it. You have 2 free parameters, p1 and p2, that need to be estimated. Why not use a zero inflated binomial model to estimate them from the data?

You should look at package VGAM, which fits ZIB models among others. In fact, you can get the expected distribution for a ZIB from the VGAM function dzibinom(), which you could use to compare your observed distribution with if you know the parameters of visitation and pollination. Again, you really should fit the ZIB model.

如果您的部分自交假设仅适用于昆虫授粉，那么预期分布只是二项式，您可以使用二项式家庭 glm 或可能具有植物 id 作为随机效应的 glmm 来估计参数。但是，如果它们可以部分自我并接受昆虫授粉，那么您将需要两个二项分布的混合。在那种情况下，我会使用 OpenBUGS 或 JAGS 来研究使用 MCMC 来拟合模型。

将这两个模型拟合到数据后，您可以使用 AIC 或 BIC 或您选择的其他一些指标比较模型以查看哪个模型更适合。

其它你可能感兴趣的问题
上一篇Cox 比例风险模型和涉及较高案例交互时的系数解释 下一篇SVM 类权重的先验选择