docs: def fit_mle(self, data, init_vals=(80, 80)):
          """%(super)s
          In addition to data, can take init_vals which allows the user to
          specify initial values for (alpha, theta) during the optimization.
  
          """
  
          if len(data) > 1:
              mu = np.mean(data)
              var = np.var(data)
              theta0 = var / mu
              alpha0 = mu / theta0
          else:
              alpha0 = init_vals[0]
              theta0 = init_vals[1]
  
          def mle(params):
              return -np.sum(np.log(self.pmf(data, params[0], params[1])))
  
          # Bounded fmin?
          alpha, theta = optim.fmin(mle, x0=[alpha0, theta0], disp=0)
  
          return alpha, theta