Speed up GLM mode

[jr-leary7]

• GLM mode is the fastest already but i imagine there are some speedups to be had
• for example, switching to using the speedglm package instead of MASS for fitting the NB GLM might be an easy improvement
  • speedglm appears to be orders of magnitude faster than MASS::glm.nb() based on some preliminary testing
  • however, this appears to partially stem from the fact that speedglm does not iteratively estimate the overdispersion parameter $\theta$
  • whether this matters is unclear, so do some testing using marge2() and see if results meaningfully differ

[jr-leary7]

commit 4dc10cd added a basic version of this to marge2() for testing