MLEModel

class dismalpy.ssm.mlemodel.MLEModel(endog, k_states, exog=None, dates=None, freq=None, **kwargs)[source]

State space model for maximum likelihood estimation

Parameters:

endog : array_like

The observed time-series process \(y\)

k_states : int

The dimension of the unobserved state process.

exog : array_like, optional

Array of exogenous regressors, shaped nobs x k. Default is no exogenous regressors.

dates : array-like of datetime, optional

An array-like object of datetime objects. If a Pandas object is given for endog, it is assumed to have a DateIndex.

freq : str, optional

The frequency of the time-series. A Pandas offset or ‘B’, ‘D’, ‘W’, ‘M’, ‘A’, or ‘Q’. This is optional if dates are given.

**kwargs

Keyword arguments may be used to provide default values for state space matrices or for Kalman filtering options. See Representation, and KalmanFilter for more details.

Notes

This class wraps the state space model with Kalman filtering to add in functionality for maximum likelihood estimation. In particular, it adds the concept of updating the state space representation based on a defined set of parameters, through the update method or updater attribute (see below for more details on which to use when), and it adds a fit method which uses a numerical optimizer to select the parameters that maximize the likelihood of the model.

The start_params update method must be overridden in the child class (and the transform and untransform methods, if needed).

Attributes

ssm (KalmanFilter) Underlying state space representation.
filter(params, transformed=True, cov_type=None, cov_kwds=None, return_ssm=False, **kwargs)

Kalman filtering

Parameters:

params : array_like

Array of parameters at which to evaluate the loglikelihood function.

transformed : boolean, optional

Whether or not params is already transformed. Default is True.

return_ssm : boolean,optional

Whether or not to return only the state space output or a full results object. Default is to return a full results object.

cov_type : str, optional

See MLEResults.fit for a description of covariance matrix types for results object.

cov_kwds : dict or None, optional

See MLEResults.get_robustcov_results for a description required keywords for alternative covariance estimators

**kwargs

Additional keyword arguments to pass to the Kalman filter. See KalmanFilter.filter for more details.

fit(*args, **kwargs)

Fits the model by maximum likelihood via Kalman filter.

Parameters:

start_params : array_like, optional

Initial guess of the solution for the loglikelihood maximization. If None, the default is given by Model.start_params.

transformed : boolean, optional

Whether or not start_params is already transformed. Default is True.

method : str, optional

The method determines which solver from scipy.optimize is used, and it can be chosen from among the following strings:

  • ‘newton’ for Newton-Raphson, ‘nm’ for Nelder-Mead
  • ‘bfgs’ for Broyden-Fletcher-Goldfarb-Shanno (BFGS)
  • ‘lbfgs’ for limited-memory BFGS with optional box constraints
  • ‘powell’ for modified Powell’s method
  • ‘cg’ for conjugate gradient
  • ‘ncg’ for Newton-conjugate gradient
  • ‘basinhopping’ for global basin-hopping solver

The explicit arguments in fit are passed to the solver, with the exception of the basin-hopping solver. Each solver has several optional arguments that are not the same across solvers. See the notes section below (or scipy.optimize) for the available arguments and for the list of explicit arguments that the basin-hopping solver supports.

cov_type : str, optional

The cov_type keyword governs the method for calculating the covariance matrix of parameter estimates. Can be one of:

  • ‘opg’ for the outer product of gradient estimator
  • ‘oim’ for the observed information matrix estimator, calculated using the method of Harvey (1989)
  • ‘cs’ for the observed information matrix estimator, calculated using a numerical (complex step) approximation of the Hessian matrix.
  • ‘delta’ for the observed information matrix estimator, calculated using a numerical (complex step) approximation of the Hessian along with the delta method (method of propagation of errors) applied to the parameter transformation function transform_params.
  • ‘robust’ for an approximate (quasi-maximum likelihood) covariance matrix that may be valid even in the presense of some misspecifications. Intermediate calculations use the ‘oim’ method.
  • ‘robust_cs’ is the same as ‘robust’ except that the intermediate calculations use the ‘cs’ method.

cov_kwds : dict or None, optional

See MLEResults.get_robustcov_results for a description required keywords for alternative covariance estimators

maxiter : int, optional

The maximum number of iterations to perform.

full_output : boolean, optional

Set to True to have all available output in the Results object’s mle_retvals attribute. The output is dependent on the solver. See LikelihoodModelResults notes section for more information.

disp : boolean, optional

Set to True to print convergence messages.

callback : callable callback(xk), optional

Called after each iteration, as callback(xk), where xk is the current parameter vector.

return_params : boolean, optional

Whether or not to return only the array of maximizing parameters. Default is False.

optim_hessian : {‘opg’,’oim’,’cs’}, optional

The method by which the Hessian is numerically approximated. ‘opg’ uses outer product of gradients, ‘oim’ uses the information matrix formula from Harvey (1989), and ‘cs’ uses second-order complex step differentiation. This keyword is only relevant if the optimization method uses the Hessian matrix.

**kwargs

Additional keyword arguments to pass to the optimizer.

Returns:

MLEResults

See also

statsmodels.base.model.LikelihoodModel.fit, MLEResults

from_formula(formula, data, subset=None)[source]

Not implemented for state space models

hessian(params, *args, **kwargs)[source]

Hessian matrix of the likelihood function, evaluated at the given parameters

Parameters:

params : array_like

Array of parameters at which to evaluate the hessian.

*args, **kwargs

Additional arguments to the loglike method.

Returns:

hessian : array

Hessian matrix evaluated at params

Notes

This is a numerical approximation.

Both *args and **kwargs are necessary because the optimizer from fit must call this function and only supports passing arguments via *args (for example scipy.optimize.fmin_l_bfgs).

impulse_responses(params, steps=1, impulse=0, orthogonalized=False, cumulative=False, **kwargs)[source]

Impulse response function

Parameters:

params : array_like

Array of model parameters.

steps : int, optional

The number of steps for which impulse responses are calculated. Default is 1. Note that the initial impulse is not counted as a step, so if steps=1, the output will have 2 entries.

impulse : int or array_like

If an integer, the state innovation to pulse; must be between 0 and k_posdef-1. Alternatively, a custom impulse vector may be provided; must be shaped k_posdef x 1.

orthogonalized : boolean, optional

Whether or not to perform impulse using orthogonalized innovations. Note that this will also affect custum impulse vectors. Default is False.

cumulative : boolean, optional

Whether or not to return cumulative impulse responses. Default is False.

**kwargs

If the model is time-varying and steps is greater than the number of observations, any of the state space representation matrices that are time-varying must have updated values provided for the out-of-sample steps. For example, if design is a time-varying component, nobs is 10, and steps is 15, a (k_endog x k_states x 5) matrix must be provided with the new design matrix values.

Returns:

impulse_responses : array

Responses for each endogenous variable due to the impulse given by the impulse argument. A (steps + 1 x k_endog) array.

Notes

Intercepts in the measurement and state equation are ignored when calculating impulse responses.

information(params)

Fisher information matrix of model

Returns -Hessian of loglike evaluated at params.

initialize()

Initialize (possibly re-initialize) a Model instance. For instance, the design matrix of a linear model may change and some things must be recomputed.

initialize_statespace(**kwargs)

Initialize the state space representation

Parameters:

**kwargs

Additional keyword arguments to pass to the state space class constructor.

loglike(params, transformed=True, **kwargs)[source]

Loglikelihood evaluation

Parameters:

params : array_like

Array of parameters at which to evaluate the loglikelihood function.

transformed : boolean, optional

Whether or not params is already transformed. Default is True.

**kwargs

Additional keyword arguments to pass to the Kalman filter. See KalmanFilter.filter for more details.

See also

update
modifies the internal state of the model to reflect new params

Notes

[R2] recommend maximizing the average likelihood to avoid scale issues; this is done automatically by the base Model fit method.

References

[R2](1, 2) Koopman, Siem Jan, Neil Shephard, and Jurgen A. Doornik. 1999. Statistical Algorithms for Models in State Space Using SsfPack 2.2. Econometrics Journal 2 (1): 107-60. doi:10.1111/1368-423X.00023.
loglikeobs(params, transformed=True, **kwargs)[source]

Loglikelihood evaluation

Parameters:

params : array_like

Array of parameters at which to evaluate the loglikelihood function.

transformed : boolean, optional

Whether or not params is already transformed. Default is True.

**kwargs

Additional keyword arguments to pass to the Kalman filter. See KalmanFilter.filter for more details.

See also

update
modifies the internal state of the Model to reflect new params

Notes

[R3] recommend maximizing the average likelihood to avoid scale issues; this is done automatically by the base Model fit method.

References

[R3](1, 2) Koopman, Siem Jan, Neil Shephard, and Jurgen A. Doornik. 1999. Statistical Algorithms for Models in State Space Using SsfPack 2.2. Econometrics Journal 2 (1): 107-60. doi:10.1111/1368-423X.00023.
observed_information_matrix(params, **kwargs)[source]

Observed information matrix

Parameters:

params : array_like, optional

Array of parameters at which to evaluate the loglikelihood function.

**kwargs

Additional keyword arguments to pass to the Kalman filter. See KalmanFilter.filter for more details.

Notes

This method is from Harvey (1989), which shows that the information matrix only depends on terms from the gradient. This implementation is partially analytic and partially numeric approximation, therefore, because it uses the analytic formula for the information matrix, with numerically computed elements of the gradient.

References

Harvey, Andrew C. 1990. Forecasting, Structural Time Series Models and the Kalman Filter. Cambridge University Press.

opg_information_matrix(params, **kwargs)[source]

Outer product of gradients information matrix

Parameters:

params : array_like, optional

Array of parameters at which to evaluate the loglikelihood function.

**kwargs

Additional arguments to the loglikeobs method.

References

Berndt, Ernst R., Bronwyn Hall, Robert Hall, and Jerry Hausman. 1974. Estimation and Inference in Nonlinear Structural Models. NBER Chapters. National Bureau of Economic Research, Inc.

param_names

(list of str) List of human readable parameter names (for parameters actually included in the model).

predict(params, exog=None, *args, **kwargs)

After a model has been fit predict returns the fitted values.

This is a placeholder intended to be overwritten by individual models.

prepare_data()[source]

Prepare data for use in the state space representation

score(params, *args, **kwargs)[source]

Compute the score function at params.

Parameters:

params : array_like

Array of parameters at which to evaluate the score.

*args, **kwargs

Additional arguments to the loglike method.

Returns:

score : array

Score, evaluated at params.

Notes

This is a numerical approximation, calculated using first-order complex step differentiation on the loglike method.

Both *args and **kwargs are necessary because the optimizer from fit must call this function and only supports passing arguments via *args (for example scipy.optimize.fmin_l_bfgs).

score_obs(params, **kwargs)[source]

Compute the score per observation, evaluated at params

Parameters:

params : array_like

Array of parameters at which to evaluate the score.

*args, **kwargs

Additional arguments to the loglike method.

Returns:

score : array (nobs, k_vars)

Score per observation, evaluated at params.

Notes

This is a numerical approximation, calculated using first-order complex step differentiation on the loglikeobs method.

set_conserve_memory(conserve_memory=None, **kwargs)[source]

Set the memory conservation method

By default, the Kalman filter computes a number of intermediate matrices at each iteration. The memory conservation options control which of those matrices are stored.

Parameters:

conserve_memory : integer, optional

Bitmask value to set the memory conservation method to. See notes for details.

**kwargs

Keyword arguments may be used to influence the memory conservation method by setting individual boolean flags.

Notes

This method is rarely used. See the corresponding function in the KalmanFilter class for details.

set_filter_method(filter_method=None, **kwargs)[source]

Set the filtering method

The filtering method controls aspects of which Kalman filtering approach will be used.

Parameters:

filter_method : integer, optional

Bitmask value to set the filter method to. See notes for details.

**kwargs

Keyword arguments may be used to influence the filter method by setting individual boolean flags. See notes for details.

Notes

This method is rarely used. See the corresponding function in the KalmanFilter class for details.

set_inversion_method(inversion_method=None, **kwargs)[source]

Set the inversion method

The Kalman filter may contain one matrix inversion: that of the forecast error covariance matrix. The inversion method controls how and if that inverse is performed.

Parameters:

inversion_method : integer, optional

Bitmask value to set the inversion method to. See notes for details.

**kwargs

Keyword arguments may be used to influence the inversion method by setting individual boolean flags. See notes for details.

Notes

This method is rarely used. See the corresponding function in the KalmanFilter class for details.

set_smoother_output(smoother_output=None, **kwargs)[source]

Set the smoother output

The smoother can produce several types of results. The smoother output variable controls which are calculated and returned.

Parameters:

smoother_output : integer, optional

Bitmask value to set the smoother output to. See notes for details.

**kwargs

Keyword arguments may be used to influence the smoother output by setting individual boolean flags.

Notes

This method is rarely used. See the corresponding function in the KalmanSmoother class for details.

set_stability_method(stability_method=None, **kwargs)[source]

Set the numerical stability method

The Kalman filter is a recursive algorithm that may in some cases suffer issues with numerical stability. The stability method controls what, if any, measures are taken to promote stability.

Parameters:

stability_method : integer, optional

Bitmask value to set the stability method to. See notes for details.

**kwargs

Keyword arguments may be used to influence the stability method by setting individual boolean flags. See notes for details.

Notes

This method is rarely used. See the corresponding function in the KalmanFilter class for details.

simulate(params, nsimulations, measurement_shocks=None, state_shocks=None, initial_state=None)[source]
Simulate a new time series following the state space model
Parameters:

params : array_like

Array of model parameters.

nsimulations : int

The number of observations to simulate. If the model is time-invariant this can be any number. If the model is time-varying, then this number must be less than or equal to the number

measurement_shocks : array_like, optional

If specified, these are the shocks to the measurement equation, :math:`

arepsilon_t`. If unspecified, these are automatically

generated using a pseudo-random number generator. If specified, must be shaped nsimulations x k_endog, where k_endog is the same as in the state space model.

state_shocks : array_like, optional

If specified, these are the shocks to the state equation, \(\eta_t\). If unspecified, these are automatically generated using a pseudo-random number generator. If specified, must be shaped nsimulations x k_posdef where k_posdef is the same as in the state space model.

initial_state : array_like, optional

If specified, this is the state vector at time zero, which should be shaped (k_states x 1), where k_states is the same as in the state space model. If unspecified, but the model has been initialized, then that initialization is used. If unspecified and the model has not been initialized, then a vector of zeros is used. Note that this is not included in the returned simulated_states array.

Returns:

simulated_obs : array

An (nsimulations x k_endog) array of simulated observations.

simulation_smoother(**kwargs)

Retrieve a simulation smoother for the statespace model.

Parameters:

simulation_output : int, optional

Determines which simulation smoother output is calculated. Default is all (including state and disturbances).

simulation_smooth_results_class : class, optional

Default results class to use to save output of simulation smoothing. Default is SimulationSmoothResults. If specified, class must extend from SimulationSmoothResults.

prefix : string

The prefix of the datatype. Usually only used internally.

**kwargs

Additional keyword arguments, used to set the simulation output. See set_simulation_output for more details.

Returns:

SimulationSmoothResults

smooth(params, transformed=True, cov_type=None, cov_kwds=None, return_ssm=False, **kwargs)

Kalman smoothing

Parameters:

params : array_like

Array of parameters at which to evaluate the loglikelihood function.

transformed : boolean, optional

Whether or not params is already transformed. Default is True.

return_ssm : boolean,optional

Whether or not to return only the state space output or a full results object. Default is to return a full results object.

cov_type : str, optional

See MLEResults.fit for a description of covariance matrix types for results object.

cov_kwds : dict or None, optional

See MLEResults.get_robustcov_results for a description required keywords for alternative covariance estimators

**kwargs

Additional keyword arguments to pass to the Kalman filter. See KalmanFilter.filter for more details.

start_params

(array) Starting parameters for maximum likelihood estimation.

transform_jacobian(unconstrained)[source]

Jacobian matrix for the parameter transformation function

Parameters:

unconstrained : array_like

Array of unconstrained parameters used by the optimizer.

Returns:

jacobian : array

Jacobian matrix of the transformation, evaluated at unconstrained

See also

transform_params

Notes

This is a numerical approximation.

transform_params(unconstrained)[source]

Transform unconstrained parameters used by the optimizer to constrained parameters used in likelihood evaluation

Parameters:

unconstrained : array_like

Array of unconstrained parameters used by the optimizer, to be transformed.

Returns:

constrained : array_like

Array of constrained parameters which may be used in likelihood evalation.

Notes

This is a noop in the base class, subclasses should override where appropriate.

untransform_params(constrained)[source]

Transform constrained parameters used in likelihood evaluation to unconstrained parameters used by the optimizer

Parameters:

constrained : array_like

Array of constrained parameters used in likelihood evalution, to be transformed.

Returns:

unconstrained : array_like

Array of unconstrained parameters used by the optimizer.

Notes

This is a noop in the base class, subclasses should override where appropriate.

update(params, transformed=True)[source]

Update the parameters of the model

Parameters:

params : array_like

Array of new parameters.

transformed : boolean, optional

Whether or not params is already transformed. If set to False, transform_params is called. Default is True.

Returns:

params : array_like

Array of parameters.

Notes

Since Model is a base class, this method should be overridden by subclasses to perform actual updating steps.

class dismalpy.ssm.mlemodel.MLEResults(model, params, smoother_results, cov_type='opg', cov_kwds=None, **kwargs)[source]

Class to hold results from fitting a state space model.

Parameters:

model : MLEModel instance

The fitted model instance

params : array

Fitted parameters

filter_results : KalmanFilter instance

The underlying state space model and Kalman filter output

Attributes

model (Model instance) A reference to the model that was fit.
filter_results (KalmanFilter instance) The underlying state space model and Kalman filter output
nobs (float) The number of observations used to fit the model.
params (array) The parameters of the model.
scale (float) This is currently set to 1.0 and not used by the model or its results.
aic()[source]

(float) Akaike Information Criterion

bic()[source]

(float) Bayes Information Criterion

conf_int(alpha=0.05, cols=None, method='default')

Returns the confidence interval of the fitted parameters.

Parameters:

alpha : float, optional

The significance level for the confidence interval. ie., The default alpha = .05 returns a 95% confidence interval.

cols : array-like, optional

cols specifies which confidence intervals to return

method : string

Not Implemented Yet Method to estimate the confidence_interval. “Default” : uses self.bse which is based on inverse Hessian for MLE “hjjh” : “jac” : “boot-bse” “boot_quant” “profile”

Returns:

conf_int : array

Each row contains [lower, upper] limits of the confidence interval for the corresponding parameter. The first column contains all lower, the second column contains all upper limits.

Notes

The confidence interval is based on the standard normal distribution. Models wish to use a different distribution should overwrite this method.

Examples

>>> import statsmodels.api as sm
>>> data = sm.datasets.longley.load()
>>> data.exog = sm.add_constant(data.exog)
>>> results = sm.OLS(data.endog, data.exog).fit()
>>> results.conf_int()
array([[-5496529.48322745, -1467987.78596704],
       [    -177.02903529,      207.15277984],
       [      -0.1115811 ,        0.03994274],
       [      -3.12506664,       -0.91539297],
       [      -1.5179487 ,       -0.54850503],
       [      -0.56251721,        0.460309  ],
       [     798.7875153 ,     2859.51541392]])
>>> results.conf_int(cols=(2,3))
array([[-0.1115811 ,  0.03994274],
       [-3.12506664, -0.91539297]])
cov_params(r_matrix=None, column=None, scale=None, cov_p=None, other=None)

Returns the variance/covariance matrix.

The variance/covariance matrix can be of a linear contrast of the estimates of params or all params multiplied by scale which will usually be an estimate of sigma^2. Scale is assumed to be a scalar.

Parameters:

r_matrix : array-like

Can be 1d, or 2d. Can be used alone or with other.

column : array-like, optional

Must be used on its own. Can be 0d or 1d see below.

scale : float, optional

Can be specified or not. Default is None, which means that the scale argument is taken from the model.

other : array-like, optional

Can be used when r_matrix is specified.

Returns:

cov : ndarray

covariance matrix of the parameter estimates or of linear combination of parameter estimates. See Notes.

Notes

(The below are assumed to be in matrix notation.)

If no argument is specified returns the covariance matrix of a model (scale)*(X.T X)^(-1)

If contrast is specified it pre and post-multiplies as follows (scale) * r_matrix (X.T X)^(-1) r_matrix.T

If contrast and other are specified returns (scale) * r_matrix (X.T X)^(-1) other.T

If column is specified returns (scale) * (X.T X)^(-1)[column,column] if column is 0d

OR

(scale) * (X.T X)^(-1)[column][:,column] if column is 1d

cov_params_cs()[source]

(array) The variance / covariance matrix. Computed using the numerical Hessian computed without using parameter transformations.

cov_params_delta()[source]

(array) The variance / covariance matrix. Computed using the numerical Hessian computed using parameter transformations and the Delta method (method of propagation of errors).

cov_params_oim()[source]

(array) The variance / covariance matrix. Computed using the method from Harvey (1989).

cov_params_opg()[source]

(array) The variance / covariance matrix. Computed using the outer product of gradients method.

cov_params_robust()[source]

(array) The QMLE variance / covariance matrix. Alias for cov_params_robust_oim

cov_params_robust_cs()[source]

(array) The QMLE variance / covariance matrix. Computed using the numerical Hessian computed without using parameter transformations as the evaluated hessian.

cov_params_robust_oim()[source]

(array) The QMLE variance / covariance matrix. Computed using the method from Harvey (1989) as the evaluated hessian.

f_test(r_matrix, cov_p=None, scale=1.0, invcov=None)

Compute the F-test for a joint linear hypothesis.

This is a special case of wald_test that always uses the F distribution.

Parameters:

r_matrix : array-like, str, or tuple

  • array : An r x k array where r is the number of restrictions to test and k is the number of regressors. It is assumed that the linear combination is equal to zero.
  • str : The full hypotheses to test can be given as a string. See the examples.
  • tuple : A tuple of arrays in the form (R, q), q can be either a scalar or a length k row vector.

cov_p : array-like, optional

An alternative estimate for the parameter covariance matrix. If None is given, self.normalized_cov_params is used.

scale : float, optional

Default is 1.0 for no scaling.

invcov : array-like, optional

A q x q array to specify an inverse covariance matrix based on a restrictions matrix.

Returns:

res : ContrastResults instance

The results for the test are attributes of this results instance.

See also

statsmodels.stats.contrast.ContrastResults, wald_test, t_test, patsy.DesignInfo.linear_constraint

Notes

The matrix r_matrix is assumed to be non-singular. More precisely,

r_matrix (pX pX.T) r_matrix.T

is assumed invertible. Here, pX is the generalized inverse of the design matrix of the model. There can be problems in non-OLS models where the rank of the covariance of the noise is not full.

Examples

>>> import numpy as np
>>> import statsmodels.api as sm
>>> data = sm.datasets.longley.load()
>>> data.exog = sm.add_constant(data.exog)
>>> results = sm.OLS(data.endog, data.exog).fit()
>>> A = np.identity(len(results.params))
>>> A = A[1:,:]

This tests that each coefficient is jointly statistically significantly different from zero.

>>> print(results.f_test(A))
<F contrast: F=330.28533923463488, p=4.98403052872e-10,
 df_denom=9, df_num=6>

Compare this to

>>> results.fvalue
330.2853392346658
>>> results.f_pvalue
4.98403096572e-10
>>> B = np.array(([0,0,1,-1,0,0,0],[0,0,0,0,0,1,-1]))

This tests that the coefficient on the 2nd and 3rd regressors are equal and jointly that the coefficient on the 5th and 6th regressors are equal.

>>> print(results.f_test(B))
<F contrast: F=9.740461873303655, p=0.00560528853174, df_denom=9,
 df_num=2>

Alternatively, you can specify the hypothesis tests using a string

>>> from statsmodels.datasets import longley
>>> from statsmodels.formula.api import ols
>>> dta = longley.load_pandas().data
>>> formula = 'TOTEMP ~ GNPDEFL + GNP + UNEMP + ARMED + POP + YEAR'
>>> results = ols(formula, dta).fit()
>>> hypotheses = '(GNPDEFL = GNP), (UNEMP = 2), (YEAR/1829 = 1)'
>>> f_test = results.f_test(hypotheses)
>>> print(f_test)
fittedvalues()[source]

(array) The predicted values of the model. An (nobs x k_endog) array.

forecast(steps=1, **kwargs)[source]

Out-of-sample forecasts

Parameters:

steps : int, str, or datetime, optional

If an integer, the number of steps to forecast from the end of the sample. Can also be a date string to parse or a datetime type. However, if the dates index does not have a fixed frequency, steps must be an integer. Default

**kwargs

Additional arguments may required for forecasting beyond the end of the sample. See FilterResults.predict for more details.

Returns:

forecast : array

Array of out of sample forecasts. A (steps x k_endog) array.

get_forecast(steps=1, **kwargs)[source]

Out-of-sample forecasts

Parameters:

steps : int, str, or datetime, optional

If an integer, the number of steps to forecast from the end of the sample. Can also be a date string to parse or a datetime type. However, if the dates index does not have a fixed frequency, steps must be an integer. Default

**kwargs

Additional arguments may required for forecasting beyond the end of the sample. See FilterResults.predict for more details.

Returns:

forecast : array

Array of out of sample forecasts. A (steps x k_endog) array.

get_prediction(start=None, end=None, dynamic=False, **kwargs)[source]

In-sample prediction and out-of-sample forecasting

Parameters:

start : int, str, or datetime, optional

Zero-indexed observation number at which to start forecasting, i.e., the first forecast is start. Can also be a date string to parse or a datetime type. Default is the the zeroth observation.

end : int, str, or datetime, optional

Zero-indexed observation number at which to end forecasting, i.e., the last forecast is end. Can also be a date string to parse or a datetime type. However, if the dates index does not have a fixed frequency, end must be an integer index if you want out of sample prediction. Default is the last observation in the sample.

dynamic : boolean, int, str, or datetime, optional

Integer offset relative to start at which to begin dynamic prediction. Can also be an absolute date string to parse or a datetime type (these are not interpreted as offsets). Prior to this observation, true endogenous values will be used for prediction; starting with this observation and continuing through the end of prediction, forecasted endogenous values will be used instead.

**kwargs

Additional arguments may required for forecasting beyond the end of the sample. See FilterResults.predict for more details.

Returns:

forecast : array

Array of out of in-sample predictions and / or out-of-sample forecasts. An (npredict x k_endog) array.

hqic()[source]

(float) Hannan-Quinn Information Criterion

impulse_responses(steps=1, impulse=0, orthogonalized=False, cumulative=False, **kwargs)[source]

Impulse response function

Parameters:

steps : int, optional

The number of steps for which impulse responses are calculated. Default is 1. Note that the initial impulse is not counted as a step, so if steps=1, the output will have 2 entries.

impulse : int or array_like

If an integer, the state innovation to pulse; must be between 0 and k_posdef-1. Alternatively, a custom impulse vector may be provided; must be shaped k_posdef x 1.

orthogonalized : boolean, optional

Whether or not to perform impulse using orthogonalized innovations. Note that this will also affect custum impulse vectors. Default is False.

cumulative : boolean, optional

Whether or not to return cumulative impulse responses. Default is False.

**kwargs

If the model is time-varying and steps is greater than the number of observations, any of the state space representation matrices that are time-varying must have updated values provided for the out-of-sample steps. For example, if design is a time-varying component, nobs is 10, and steps is 15, a (k_endog x k_states x 5) matrix must be provided with the new design matrix values.

Returns:

impulse_responses : array

Responses for each endogenous variable due to the impulse given by the impulse argument. A (steps + 1 x k_endog) array.

Notes

Intercepts in the measurement and state equation are ignored when calculating impulse responses.

kalman_gain

Kalman gain matrices

llf()[source]

(float) The value of the log-likelihood function evaluated at params.

llf_obs()[source]

(float) The value of the log-likelihood function evaluated at params.

load(fname)

load a pickle, (class method)

Parameters:

fname : string or filehandle

fname can be a string to a file path or filename, or a filehandle.

Returns:

unpickled instance

loglikelihood_burn()[source]

(float) The number of observations during which the likelihood is not evaluated.

plot_diagnostics(variable=0, lags=10, fig=None, figsize=None)[source]

Diagnostic plots for standardized residuals of one endogenous variable

Parameters:

variable : integer, optional

Index of the endogenous variable for which the diagnostic plots should be created. Default is 0.

lags : integer, optional

Number of lags to include in the correlogram. Default is 10.

fig : Matplotlib Figure instance, optional

If given, subplots are created in this figure instead of in a new figure. Note that the 2x2 grid will be created in the provided figure using fig.add_subplot().

figsize : tuple, optional

If a figure is created, this argument allows specifying a size. The tuple is (width, height).

See also

statsmodels.graphics.gofplots.qqplot, statsmodels.graphics.tsaplots.plot_acf

Notes

Produces a 2x2 plot grid with the following plots (ordered clockwise from top left):

  1. Standardized residuals over time
  2. Histogram plus estimated density of standardized residulas, along with a Normal(0,1) density plotted for reference.
  3. Normal Q-Q plot, with Normal reference line.
  4. Correlogram
predict(start=None, end=None, dynamic=False, **kwargs)[source]

In-sample prediction and out-of-sample forecasting

Parameters:

start : int, str, or datetime, optional

Zero-indexed observation number at which to start forecasting, i.e., the first forecast is start. Can also be a date string to parse or a datetime type. Default is the the zeroth observation.

end : int, str, or datetime, optional

Zero-indexed observation number at which to end forecasting, i.e., the last forecast is end. Can also be a date string to parse or a datetime type. However, if the dates index does not have a fixed frequency, end must be an integer index if you want out of sample prediction. Default is the last observation in the sample.

dynamic : boolean, int, str, or datetime, optional

Integer offset relative to start at which to begin dynamic prediction. Can also be an absolute date string to parse or a datetime type (these are not interpreted as offsets). Prior to this observation, true endogenous values will be used for prediction; starting with this observation and continuing through the end of prediction, forecasted endogenous values will be used instead.

**kwargs

Additional arguments may required for forecasting beyond the end of the sample. See FilterResults.predict for more details.

Returns:

forecast : array

Array of out of in-sample predictions and / or out-of-sample forecasts. An (npredict x k_endog) array.

pvalues()[source]

(array) The p-values associated with the z-statistics of the coefficients. Note that the coefficients are assumed to have a Normal distribution.

remove_data()

remove data arrays, all nobs arrays from result and model

This reduces the size of the instance, so it can be pickled with less memory. Currently tested for use with predict from an unpickled results and model instance.

Warning

Since data and some intermediate results have been removed calculating new statistics that require them will raise exceptions. The exception will occur the first time an attribute is accessed that has been set to None.

Not fully tested for time series models, tsa, and might delete too much for prediction or not all that would be possible.

The list of arrays to delete is maintained as an attribute of the result and model instance, except for cached values. These lists could be changed before calling remove_data.

resid()[source]

(array) The model residuals. An (nobs x k_endog) array.

save(fname, remove_data=False)

save a pickle of this instance

Parameters:

fname : string or filehandle

fname can be a string to a file path or filename, or a filehandle.

remove_data : bool

If False (default), then the instance is pickled without changes. If True, then all arrays with length nobs are set to None before pickling. See the remove_data method. In some cases not all arrays will be set to None.

Notes

If remove_data is true and the model result does not implement a remove_data method then this will raise an exception.

simulate(nsimulations, measurement_shocks=None, state_shocks=None, initial_state=None)[source]
Simulate a new time series following the state space model
Parameters:

nsimulations : int

The number of observations to simulate. If the model is time-invariant this can be any number. If the model is time-varying, then this number must be less than or equal to the number

measurement_shocks : array_like, optional

If specified, these are the shocks to the measurement equation, :math:`

arepsilon_t`. If unspecified, these are automatically

generated using a pseudo-random number generator. If specified, must be shaped nsimulations x k_endog, where k_endog is the same as in the state space model.

state_shocks : array_like, optional

If specified, these are the shocks to the state equation, \(\eta_t\). If unspecified, these are automatically generated using a pseudo-random number generator. If specified, must be shaped nsimulations x k_posdef where k_posdef is the same as in the state space model.

initial_state : array_like, optional

If specified, this is the state vector at time zero, which should be shaped (k_states x 1), where k_states is the same as in the state space model. If unspecified, but the model has been initialized, then that initialization is used. If unspecified and the model has not been initialized, then a vector of zeros is used. Note that this is not included in the returned simulated_states array.

Returns:

simulated_obs : array

An (nsimulations x k_endog) array of simulated observations.

summary(alpha=0.05, start=None, title=None, model_name=None, display_params=True)[source]

Summarize the Model

Parameters:

alpha : float, optional

Significance level for the confidence intervals. Default is 0.05.

start : int, optional

Integer of the start observation. Default is 0.

model_name : string

The name of the model used. Default is to use model class name.

Returns:

summary : Summary instance

This holds the summary table and text, which can be printed or converted to various output formats.

See also

statsmodels.iolib.summary.Summary

t_test(r_matrix, cov_p=None, scale=None, use_t=None)

Compute a t-test for a each linear hypothesis of the form Rb = q

Parameters:

r_matrix : array-like, str, tuple

  • array : If an array is given, a p x k 2d array or length k 1d array specifying the linear restrictions. It is assumed that the linear combination is equal to zero.
  • str : The full hypotheses to test can be given as a string. See the examples.
  • tuple : A tuple of arrays in the form (R, q). If q is given, can be either a scalar or a length p row vector.

cov_p : array-like, optional

An alternative estimate for the parameter covariance matrix. If None is given, self.normalized_cov_params is used.

scale : float, optional

An optional scale to use. Default is the scale specified by the model fit.

use_t : bool, optional

If use_t is None, then the default of the model is used. If use_t is True, then the p-values are based on the t distribution. If use_t is False, then the p-values are based on the normal distribution.

Returns:

res : ContrastResults instance

The results for the test are attributes of this results instance. The available results have the same elements as the parameter table in summary().

See also

tvalues
individual t statistics
f_test
for F tests

patsy.DesignInfo.linear_constraint

Examples

>>> import numpy as np
>>> import statsmodels.api as sm
>>> data = sm.datasets.longley.load()
>>> data.exog = sm.add_constant(data.exog)
>>> results = sm.OLS(data.endog, data.exog).fit()
>>> r = np.zeros_like(results.params)
>>> r[5:] = [1,-1]
>>> print(r)
[ 0.  0.  0.  0.  0.  1. -1.]

r tests that the coefficients on the 5th and 6th independent variable are the same.

>>> T_test = results.t_test(r)
>>> print(T_test)
<T contrast: effect=-1829.2025687192481, sd=455.39079425193762,
t=-4.0167754636411717, p=0.0015163772380899498, df_denom=9>
>>> T_test.effect
-1829.2025687192481
>>> T_test.sd
455.39079425193762
>>> T_test.tvalue
-4.0167754636411717
>>> T_test.pvalue
0.0015163772380899498

Alternatively, you can specify the hypothesis tests using a string

>>> from statsmodels.formula.api import ols
>>> dta = sm.datasets.longley.load_pandas().data
>>> formula = 'TOTEMP ~ GNPDEFL + GNP + UNEMP + ARMED + POP + YEAR'
>>> results = ols(formula, dta).fit()
>>> hypotheses = 'GNPDEFL = GNP, UNEMP = 2, YEAR/1829 = 1'
>>> t_test = results.t_test(hypotheses)
>>> print(t_test)
test_heteroskedasticity(method, alternative='two-sided', use_f=True)[source]

Test for heteroskedasticity of standardized residuals

Tests whether the sum-of-squares in the first third of the sample is significantly different than the sum-of-squares in the last third of the sample. Analogous to a Goldfeld-Quandt test.

Parameters:

method : string {‘breakvar’} or None

The statistical test for heteroskedasticity. Must be ‘breakvar’ for test of a break in the variance. If None, an attempt is made to select an appropriate test.

alternative : string, ‘increasing’, ‘decreasing’ or ‘two-sided’

This specifies the alternative for the p-value calculation. Default is two-sided.

use_f : boolean, optional

Whether or not to compare against the asymptotic distribution (chi-squared) or the approximate small-sample distribution (F). Default is True (i.e. default is to compare against an F distribution).

Returns:

output : array

An array with (test_statistic, pvalue) for each endogenous variable. The array is then sized (k_endog, 2). If the method is called as het = res.test_heteroskedasticity(), then het[0] is an array of size 2 corresponding to the first endogenous variable, where het[0][0] is the test statistic, and het[0][1] is the p-value.

Notes

The null hypothesis is of no heteroskedasticity. That means different things depending on which alternative is selected:

  • Increasing: Null hypothesis is that the variance is not increasing throughout the sample; that the sum-of-squares in the later subsample is not greater than the sum-of-squares in the earlier subsample.
  • Decreasing: Null hypothesis is that the variance is not decreasing throughout the sample; that the sum-of-squares in the earlier subsample is not greater than the sum-of-squares in the later subsample.
  • Two-sided: Null hypothesis is that the variance is not changing throughout the sample. Both that the sum-of-squares in the earlier subsample is not greater than the sum-of-squares in the later subsample and that the sum-of-squares in the later subsample is not greater than the sum-of-squares in the earlier subsample.

For \(h = [T/3]\), the test statistic is:

\[H(h) = \sum_{t=T-h+1}^T \tilde v_t^2 \left / \sum_{t=d+1}^{d+1+h} \tilde v_t^2 \right .\]

where \(d\) is the number of periods in which the loglikelihood was burned in the parent model (usually corresponding to diffuse initialization).

This statistic can be tested against an \(F(h,h)\) distribution. Alternatively, \(h H(h)\) is asymptotically distributed according to \(\chi_h^2\); this second test can be applied by passing asymptotic=True as an argument.

See section 5.4 of [R4] for the above formula and discussion, as well as additional details.

TODO

  • Allow specification of \(h\)

References

[R4](1, 2) Harvey, Andrew C. 1990. Forecasting, Structural Time Series Models and the Kalman Filter. Cambridge University Press.
test_normality(method)[source]

Test for normality of standardized residuals.

Null hypothesis is normality.

Parameters:

method : string {‘jarquebera’} or None

The statistical test for normality. Must be ‘jarquebera’ for Jarque-Bera normality test. If None, an attempt is made to select an appropriate test.

See also

statsmodels.stats.stattools.jarque_bera

Notes

If the first d loglikelihood values were burned (i.e. in the specified model, loglikelihood_burn=d), then this test is calculated ignoring the first d residuals.

test_serial_correlation(method, lags=None)[source]

Ljung-box test for no serial correlation of standardized residuals

Null hypothesis is no serial correlation.

Parameters:

method : string {‘ljungbox’,’boxpierece’} or None

The statistical test for serial correlation. If None, an attempt is made to select an appropriate test.

lags : None, int or array_like

If lags is an integer then this is taken to be the largest lag that is included, the test result is reported for all smaller lag length. If lags is a list or array, then all lags are included up to the largest lag in the list, however only the tests for the lags in the list are reported. If lags is None, then the default maxlag is 12*(nobs/100)^{1/4}

Returns:

output : array

An array with (test_statistic, pvalue) for each endogenous variable and each lag. The array is then sized (k_endog, 2, lags). If the method is called as ljungbox = res.test_serial_correlation(), then ljungbox[i] holds the results of the Ljung-Box test (as would be returned by statsmodels.stats.diagnostic.acorr_ljungbox) for the `i`th endogenous variable.

See also

statsmodels.stats.diagnostic.acorr_ljungbox

Notes

If the first d loglikelihood values were burned (i.e. in the specified model, loglikelihood_burn=d), then this test is calculated ignoring the first d residuals.

tvalues()

Return the t-statistic for a given parameter estimate.

wald_test(r_matrix, cov_p=None, scale=1.0, invcov=None, use_f=None)

Compute a Wald-test for a joint linear hypothesis.

Parameters:

r_matrix : array-like, str, or tuple

  • array : An r x k array where r is the number of restrictions to test and k is the number of regressors. It is assumed that the linear combination is equal to zero.
  • str : The full hypotheses to test can be given as a string. See the examples.
  • tuple : A tuple of arrays in the form (R, q), q can be either a scalar or a length p row vector.

cov_p : array-like, optional

An alternative estimate for the parameter covariance matrix. If None is given, self.normalized_cov_params is used.

scale : float, optional

Default is 1.0 for no scaling.

invcov : array-like, optional

A q x q array to specify an inverse covariance matrix based on a restrictions matrix.

use_f : bool

If True, then the F-distribution is used. If False, then the asymptotic distribution, chisquare is used. If use_f is None, then the F distribution is used if the model specifies that use_t is True. The test statistic is proportionally adjusted for the distribution by the number of constraints in the hypothesis.

Returns:

res : ContrastResults instance

The results for the test are attributes of this results instance.

See also

statsmodels.stats.contrast.ContrastResults, f_test, t_test, patsy.DesignInfo.linear_constraint

Notes

The matrix r_matrix is assumed to be non-singular. More precisely,

r_matrix (pX pX.T) r_matrix.T

is assumed invertible. Here, pX is the generalized inverse of the design matrix of the model. There can be problems in non-OLS models where the rank of the covariance of the noise is not full.

wald_test_terms(skip_single=False, extra_constraints=None, combine_terms=None)

Compute a sequence of Wald tests for terms over multiple columns

This computes joined Wald tests for the hypothesis that all coefficients corresponding to a term are zero.

Terms are defined by the underlying formula or by string matching.

Parameters:

skip_single : boolean

If true, then terms that consist only of a single column and, therefore, refers only to a single parameter is skipped. If false, then all terms are included.

extra_constraints : ndarray

not tested yet

combine_terms : None or list of strings

Each string in this list is matched to the name of the terms or the name of the exogenous variables. All columns whose name includes that string are combined in one joint test.

Returns:

test_result : result instance

The result instance contains table which is a pandas DataFrame with the test results: test statistic, degrees of freedom and pvalues.

Examples

>>> res_ols = ols("np.log(Days+1) ~ C(Duration, Sum)*C(Weight, Sum)",
                  data).fit()
>>> res_ols.wald_test_terms()
<class 'statsmodels.stats.contrast.WaldTestResults'>
                                          F                P>F  df constraint  df denom
Intercept                        279.754525  2.37985521351e-22              1        51
C(Duration, Sum)                   5.367071    0.0245738436636              1        51
C(Weight, Sum)                    12.432445  3.99943118767e-05              2        51
C(Duration, Sum):C(Weight, Sum)    0.176002      0.83912310946              2        51
>>> res_poi = Poisson.from_formula("Days ~ C(Weight) * C(Duration)",
                                   data).fit(cov_type='HC0')
>>> wt = res_poi.wald_test_terms(skip_single=False,
                                 combine_terms=['Duration', 'Weight'])
>>> print(wt)
                            chi2             P>chi2  df constraint
Intercept              15.695625  7.43960374424e-05              1
C(Weight)              16.132616  0.000313940174705              2
C(Duration)             1.009147     0.315107378931              1
C(Weight):C(Duration)   0.216694     0.897315972824              2
Duration               11.187849     0.010752286833              3
Weight                 30.263368  4.32586407145e-06              4
zvalues()[source]

(array) The z-statistics for the coefficients.