docs: def _ln_choose(n, k_agg):
      '''
      log binomial coefficient with extended gamma factorials. n and k_agg may be
      int or array - if both array, must be the same length.
  
      '''
      gammaln = special.gammaln
      return gammaln(n + 1) - (gammaln(k_agg + 1) + gammaln(n - k_agg + 1))