docs: def pdf(self):
          r"""
          Generate the vector of probabilities for the Beta-binomial
          (n, a, b) distribution.
  
          The Beta-binomial distribution takes the form
  
          .. math::
              p(k \,|\, n, a, b) =
              {n \choose k} \frac{B(k + a, n - k + b)}{B(a, b)},
              \qquad k = 0, \ldots, n,
  
          where :math:`B` is the beta function.
  
          Parameters
          ----------
          n : scalar(int)
              First parameter to the Beta-binomial distribution
          a : scalar(float)
              Second parameter to the Beta-binomial distribution
          b : scalar(float)
              Third parameter to the Beta-binomial distribution
  
          Returns
          -------
          probs: array_like(float)
              Vector of probabilities over k
  
          """
          n, a, b = self.n, self.a, self.b
          k = np.arange(n + 1)
          probs = binom(n, k) * beta(k + a, n - k + b) / beta(a, b)
          return probs