python 矩量法_矩量法：β二项分布

我想我已经能够用python准确地编写这个代码：import numpy as np

males = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12], dtype=float)

fams = np.array([3, 24, 104, 286, 670, 1033, 1343, 1112, 829, 478, 181, 45, 7], dtype=float)

n = 12

k = 1

m_1 = (sum(fams*males**k))/(sum(fams))

k = 2

m_2 = (sum(fams*males**k))/(sum(fams))

alpha = (n*m_1-m_2)/(n*(m_2/m_1-m_1-1)+m_1)

beta = (n-m_1)*(n-m_2/m_1)/(n*(m_2/m_1-m_1-1)+m_1)

print "n =", n

print "m_1 =", m_1

print "m_2 =", m_2

print "alpha =", alpha

print "beta =", beta

哪些输出：

^{pr2}$

这与示例的结果相同。但是，如果我使用R包VGAM，它使用最大似然，α和β的估计值完全不同x = c(0,1,2,3,4,5,6,7,8,9,10,11,12)

y = c(3,24,104,286,670,1033,1343,1112,829,478,181,45,7)

library("VGAM")

fit=vglm(cbind(x, y) ~ 1, betabinomialff, trace = TRUE)

Coef(fit)

shape1 shape2

0.4241806 4.9069560

我做错什么了吗？在