I never looked at how to calculate the parameters for the fraction dispersing. The trick may be the betabinomial function in the VGAM library.
The first step is to amalgamate all the dispersing seeds from disperseLer:
nondispersers <- subset(disperseLer, Pot == 0, c("ID", "Seedlings"))
dispersers <- filter(disperseLer, Pot == 1) %>% group_by(ID) %>% summarise(dispersers = sum(Seedlings))
disperse_num <- merge(nondispersers, dispersers)
names(disperse_num)[2] <- "nondispersers"
disperse_num$dispersers <- round(2 * disperse_num$dispersers)So now we make an intercept-only model with vgam.
library(VGAM)
bbfit <- vglm(cbind(dispersers, nondispersers) ~ 1, betabinomial, data = disperse_num)
summary(bbfit)
Call:
vglm(formula = cbind(dispersers, nondispersers) ~ 1, family = betabinomial,
data = disperse_num)
Pearson residuals:
Min 1Q Median 3Q Max
logitlink(mu) -2.7263 -0.4860 0.05773 0.5152 2.162
logitlink(rho) -0.7111 -0.6557 -0.49997 0.1977 4.045
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept):1 0.05167 0.12312 0.420 0.675
(Intercept):2 -1.78155 0.21624 -8.239 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Names of linear predictors: logitlink(mu), logitlink(rho)
Log-likelihood: -238.6762 on 74 degrees of freedom
Number of Fisher scoring iterations: 3
No Hauck-Donner effect found in any of the estimatesCoef(bbfit) mu rho
0.5129148 0.1441117 The parameterization in vgam is somewhat unconventional. mu is simply the mean of the beta distribution, alpha/(alpha + beta), but rho is 1/(1 + alpha + beta) where alpha and beta are the conventional shape parameters.
To get to the variance parameter, we need to transform this. \[\begin{align} V &= \frac{\alpha \beta}{(\alpha + \beta)^2 (\alpha + \beta + 1)}\\ &= \mu (1-\mu) \rho. \end{align}\]
So we have
fd_mean <- Coef(bbfit)[1]
fd_stdv <- sqrt(fd_mean * (1 - fd_mean) * Coef(bbfit)[2])
fd_mean mu
0.5129148 fd_stdv mu
0.1897469