Maximum Likelihood Estimation

Stan provides optimization algorithms which find modes of the density specified by a Stan program. Three different algorithms are available: a Newton optimizer, and two related quasi-Newton algorithms, BFGS and L-BFGS. The L-BFGS algorithm is the default optimizer. Newton’s method is the least efficient of the three, but has the advantage of setting its own stepsize.

Optimize configuration

  • algorithm: Algorithm to use. One of: “BFGS”, “LBFGS”, “Newton”.
  • init_alpha: Line search step size for first iteration.
  • iter: Total number of iterations.
  • data: Values for all data variables in the model, specified either as a dictionary with entries matching the data variables, or as the path of a data file in JSON or Rdump format.
  • seed: The seed for random number generator.
  • inits: Specifies how the sampler initializes parameter values.
  • output_dir: A path or file name which will be used as the basename for the CmdStan output files.
  • save_diagnostics: Whether or not to save diagnostics.

All of these arguments are optional; when unspecified, the CmdStan defaults will be used.

Example: estamate MLE for model bernoulli.stan by optimization

In this example we use the CmdStan example model bernoulli.stan and data file

The CmdStanModel class method optimize returns a CmdStanMLE object which provides properties to retrieve the estimate of the penalized maximum likelihood estaimate of all model parameters:

  • column_names
  • optimized_params_dict
  • optimized_params_np
  • optimized_params_pd

In the following example, we instantiate a model and do optimization using the default CmdStan settings:

import os
from cmdstanpy.model import CmdStanModel
from cmdstanpy.utils import cmdstan_path

# instantiate, compile bernoulli model
bernoulli_path = os.path.join(cmdstan_path(), 'examples', 'bernoulli', 'bernoulli.stan')
bernoulli_model = CmdStanModel(stan_file=bernoulli_path)

# run CmdStan's optimize method, returns object `CmdStanMLE`
bern_data = os.path.join(cmdstan_path(), 'examples', 'bernoulli', '')
bern_mle = bernoulli_model.optimize(data=bernoulli_data)