在 flexmix 中拟合二项式混合物 - 我做错了什么？

Question

我正在尝试理解 flexmix，特别是在尝试拟合最简单的可想象的二项式混合模型（两个仅截距模型的混合）时我做错了什么。

set.seed(42)
data=data.frame(class=rbinom(1000,size=1,prob=0.3)) %>% # 30% in class 1, 70% in class 0
  mutate(yes=map_dbl(class,~ifelse(.x,rbinom(1,1,prob=0.8),rbinom(1,1,prob=0.4)))) # class 1 = 80%, class 2 = 40%

# Algebraic
(mean(data$yes==1)-0.4)/(0.8-0.4)
# = 0.295 this is what the model should recover, right?


library(flexmix)
mo1=FLXMRglm(offset=rep(log(.4/.6), nrow(data)),family="binomial")
mo2=FLXMRglm(offset=rep(log(.8/.2), nrow(data)),family="binomial")
flexfit <- flexmix(cbind(yes,1-yes) ~ -1, data=data, k=2, model=list(mo1, mo2))
flexfit
summary(flexfit)

此代码返回属于第一类的所有数据点。我是否错误地设置了模型？还是我对混合模型的工作方式存在更根本的误解？

Answer 1

由于您没有回归量，我认为您需要在线性模型中指定一个截距，即cbind(yes, 1 - yes) ~ 1。

flexfit <- initFlexmix(cbind(yes, 1 - yes) ~ 1, data=data, k=2, model=list(mo1, mo2))
summary(flexfit)

将数据集拆分为 482 和 518 个单元的集群。 该模型的对数似然 = -692.499，AIC = 1390.998，BIC = 1419.537。

请注意，二项分布的混合通常无法识别。 第 2 组有 518 个单位并非巧合：它们正是 yes == 1 的观测数。 如果我们将标签传递给模型，我们显然会恢复确切的答案（使模型完全无用）：

data$label <- ifelse(data$class == 1, 2, 1)
flexfit <- initFlexmix(cbind(yes, 1 - yes) ~ 1 | label, data=data, k=2, model=list(mo1, mo2))
summary(flexfit)

我希望有人能证明我错了，但我认为没有办法通过仅观察 0/1 结果且没有强先验假设来恢复隐藏组。

这样想：你可以通过多种方式将你的二元观测分成两个簇（总数是第二类的斯特林数），每个分配都是对你的数据的完全合理的解释。最大似然估计是最合理的一种：通过观察到的相对频率估计隐藏类的概率，然后根据类将条件概率设置为0或1。

Pr(Y = 1 | π, θ, φ) = Pr(C = 1) Pr(Y = 1 | C = 1) + Pr(C = 2) Pr(Y = 1 | C = 2) = π θ + (1 - π) φ

将 π 设置为观察到的相对频率，φ = 1 和 θ = 0。 flexmix 模型比这更复杂，但考虑到我们只有一个截距，我确信它可以简化为我的示例。