Maximum likelihood estimation of k of the negative binomial

# Input observed freq dist f=0:8 # parasites per host nobs=c(164,147,92,43,25,11,3,0,1) # hosts sum(nobs) # total no of hosts # reconstruct data into variable 'Parasites' Parasites=rep(f, nobs) # obtain their sum, mean and variance (tot=sum(Parasites)) (m=mean(Parasites)) (v=var(Parasites)) (kapr=m^2/(v-m)) # get approximate k # find expected freqs dnbinom(f, s=kapr, mu=m)*tot # Now for a less biased estimate library(MASS) # Load MASS library # Obtain max. lik. est. of k (kmle=theta.ml(Parasites,m)[1]) # find expected freqs dnbinom(f, s=kmle, mu=m)*tot