from collections import OrderedDict
import contextlib
import warnings
import numpy as np
import pandas as pd
from scipy import stats
from statsmodels.base.data import PandasData
from statsmodels.tools._decorators import cache_readonly
from statsmodels.tools.eval_measures import aic, aicc, bic, hqic
from statsmodels.tools.numdiff import (
_get_epsilon,
approx_fprime,
approx_fprime_cs,
approx_hess_cs,
)
from statsmodels.tools.sm_exceptions import PrecisionWarning
from statsmodels.tools.tools import pinv_extended
from statsmodels.tools.validation import string_like
import statsmodels.tsa.base.tsa_model as tsbase
from statsmodels.tsa.statespace.tools import _safe_cond
from statsmodels.tsa.stattools._stattools import (
breakvar_heteroskedasticity_test,
)
class StateSpaceMLEModel(tsbase.TimeSeriesModel):
"""
Temporary base model for use with the exponential smoothing models; copies
the functionality needed from statespace.mlemodel.MLEModel
Parameters
----------
endog : array_like
The observed time-series process :math:`y`.
exog : array_like, optional
Array of exogenous regressors.
dates : array_like, optional
An array-like object of datetime objects. If a pandas object is given
for endog or exog, 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.
missing : str, optional
Available options are 'none', 'drop', and 'raise'. If 'none', no nan
checking is done. If 'drop', any observations with nans are dropped.
If 'raise', an error is raised. Default is 'none'.
**kwargs
Keyword arguments used to construct the state space representation
of the model and stored so that the model can be recreated.
"""
def __init__(
self, endog, exog=None, dates=None, freq=None, missing="none", **kwargs
):
# TODO: this was changed from the original, requires some work when
# using this as base class for state space and exponential smoothing
super().__init__(
endog=endog, exog=exog, dates=dates, freq=freq, missing=missing
)
# Store kwargs to recreate model
self._init_kwargs = kwargs
# Prepared the endog array: C-ordered, shape=(nobs x k_endog)
self.endog, self.exog = self.prepare_data(self.data)
self.use_pandas = isinstance(self.data, PandasData)
# Dimensions
self.nobs = self.endog.shape[0]
# Setup holder for fixed parameters
self._has_fixed_params = False
self._fixed_params = None
self._params_index = None
self._fixed_params_index = None
self._free_params_index = None
@staticmethod
def prepare_data(data):
raise NotImplementedError
def clone(self, endog, exog=None, **kwargs):
raise NotImplementedError
def _validate_can_fix_params(self, param_names):
for param_name in param_names:
if param_name not in self.param_names:
raise ValueError(
f'Invalid parameter name passed: "{param_name}".'
)
@property
def k_params(self):
"""The number of parameters in the model"""
return len(self.param_names)
@contextlib.contextmanager
def fix_params(self, params):
"""
Fix parameters to specific values (context manager)
Parameters
----------
params : dict
Dictionary describing the fixed parameter values, of the form
`param_name: fixed_value`. See the `param_names` property for valid
parameter names.
Examples
--------
>>> mod = sm.tsa.SARIMAX(endog, order=(1, 0, 1))
>>> with mod.fix_params({'ar.L1': 0.5}):
res = mod.fit()
"""
# Initialization (this is done here rather than in the constructor
# because param_names may not be available at that point)
if self._fixed_params is None:
self._fixed_params = {}
self._params_index = OrderedDict(
zip(self.param_names, np.arange(self.k_params), strict=True)
)
# Cache the current fixed parameters
cache_fixed_params = self._fixed_params.copy()
cache_has_fixed_params = self._has_fixed_params
cache_fixed_params_index = self._fixed_params_index
cache_free_params_index = self._free_params_index
# Validate parameter names and values
all_fixed_param_names = (
set(params.keys()) | set(self._fixed_params.keys())
)
self._validate_can_fix_params(all_fixed_param_names)
# Set the new fixed parameters, keeping the order as given by
# param_names
self._fixed_params.update(params)
self._fixed_params = OrderedDict(
[
(name, self._fixed_params[name])
for name in self.param_names
if name in self._fixed_params
]
)
# Update associated values
self._has_fixed_params = True
self._fixed_params_index = [
self._params_index[key] for key in self._fixed_params
]
self._free_params_index = list(
set(np.arange(self.k_params)).difference(self._fixed_params_index)
)
try:
yield
finally:
# Reset the fixed parameters
self._has_fixed_params = cache_has_fixed_params
self._fixed_params = cache_fixed_params
self._fixed_params_index = cache_fixed_params_index
self._free_params_index = cache_free_params_index
def fit_constrained(self, constraints, start_params=None, **fit_kwds):
"""
Fit the model with some parameters subject to equality constraints
Parameters
----------
constraints : dict
Dictionary of constraints, of the form `param_name: fixed_value`.
See the `param_names` property for valid parameter names.
start_params : array_like, optional
Initial guess of the solution for the loglikelihood maximization.
If None, the default is given by Model.start_params.
**fit_kwds : keyword arguments
fit_kwds are used in the optimization of the remaining parameters.
Returns
-------
results : Results instance
Examples
--------
>>> mod = sm.tsa.SARIMAX(endog, order=(1, 0, 1))
>>> res = mod.fit_constrained({'ar.L1': 0.5})
"""
with self.fix_params(constraints):
res = self.fit(start_params, **fit_kwds)
return res
@property
def start_params(self):
"""Starting parameters for maximum likelihood estimation"""
if hasattr(self, "_start_params"):
return self._start_params
else:
raise NotImplementedError
@property
def param_names(self):
"""List of human readable parameter names (for parameters actually included in the model)"""
if hasattr(self, "_param_names"):
return self._param_names
else:
try:
names = [f"param.{i:d}" for i in range(len(self.start_params))]
except NotImplementedError:
names = []
return names
@classmethod
def from_formula(
cls, formula, data, subset=None, drop_cols=None, *args, **kwargs
):
"""Not implemented for state space models"""
raise NotImplementedError
def _wrap_data(self, data, start_idx, end_idx, names=None):
# TODO: check if this is reasonable for statespace
# squeezing data: data may be:
# - m x n: m dates, n simulations -> squeeze does nothing
# - m x 1: m dates, 1 simulation -> squeeze removes last dimension
# - 1 x n: don't squeeze, already fine
# - 1 x 1: squeeze only second axis
if data.ndim > 1 and data.shape[1] == 1:
data = np.squeeze(data, axis=1)
if self.use_pandas:
if data.shape[0]:
_, _, _, index = self._get_prediction_index(start_idx, end_idx)
else:
index = None
if data.ndim < 2:
data = pd.Series(data, index=index, name=names)
else:
data = pd.DataFrame(data, index=index, columns=names)
return data
def _wrap_results(
self,
params,
result,
return_raw,
cov_type=None,
cov_kwds=None,
results_class=None,
wrapper_class=None,
):
if not return_raw:
# Wrap in a results object
result_kwargs = {}
if cov_type is not None:
result_kwargs["cov_type"] = cov_type
if cov_kwds is not None:
result_kwargs["cov_kwds"] = cov_kwds
if results_class is None:
results_class = self._res_classes["fit"][0]
if wrapper_class is None:
wrapper_class = self._res_classes["fit"][1]
res = results_class(self, params, result, **result_kwargs)
result = wrapper_class(res)
return result
def _score_complex_step(self, params, **kwargs):
# the default epsilon can be too small
# inversion_method = INVERT_UNIVARIATE | SOLVE_LU
epsilon = _get_epsilon(params, 2., None, len(params))
kwargs["transformed"] = True
kwargs["complex_step"] = True
return approx_fprime_cs(params, self.loglike, epsilon=epsilon,
kwargs=kwargs)
def _score_finite_difference(self, params, approx_centered=False,
**kwargs):
kwargs["transformed"] = True
return approx_fprime(params, self.loglike, kwargs=kwargs,
centered=approx_centered)
def _hessian_finite_difference(self, params, approx_centered=False,
**kwargs):
params = np.array(params, ndmin=1)
warnings.warn("Calculation of the Hessian using finite differences"
" is usually subject to substantial approximation"
" errors.",
PrecisionWarning,
stacklevel=3,
)
if not approx_centered:
epsilon = _get_epsilon(params, 3, None, len(params))
else:
epsilon = _get_epsilon(params, 4, None, len(params)) / 2
hessian = approx_fprime(params, self._score_finite_difference,
epsilon=epsilon, kwargs=kwargs,
centered=approx_centered)
# TODO: changed this to nobs_effective, has to be changed when merging
# with statespace mlemodel
return hessian / (self.nobs_effective)
def _hessian_complex_step(self, params, **kwargs):
"""Hessian matrix computed by second-order complex-step differentiation on the `loglike` function"""
# the default epsilon can be too small
epsilon = _get_epsilon(params, 3., None, len(params))
kwargs["transformed"] = True
kwargs["complex_step"] = True
hessian = approx_hess_cs(
params, self.loglike, epsilon=epsilon, kwargs=kwargs)
# TODO: changed this to nobs_effective, has to be changed when merging
# with statespace mlemodel
return hessian / (self.nobs_effective)
class StateSpaceMLEResults(tsbase.TimeSeriesModelResults):
r"""
Class to hold results from fitting a state space model
Parameters
----------
model : Model instance
The fitted model instance
params : ndarray
Fitted parameters
scale : float, optional
An optional scale parameter for the model. Default is 1.0.
Attributes
----------
model : Model instance
A reference to the model that was fit.
nobs : float
The number of observations used to fit the model.
params : ndarray
The parameters of the model.
"""
def __init__(self, model, params, scale=1.0):
self.data = model.data
self.endog = model.data.orig_endog
super().__init__(model, params, None, scale=scale)
# Save the fixed parameters
self._has_fixed_params = self.model._has_fixed_params
self._fixed_params_index = self.model._fixed_params_index
self._free_params_index = self.model._free_params_index
# TODO: seems like maybe self.fixed_params should be the dictionary
# itself, not just the keys?
if self._has_fixed_params:
self._fixed_params = self.model._fixed_params.copy()
self.fixed_params = list(self._fixed_params.keys())
else:
self._fixed_params = None
self.fixed_params = []
self.param_names = [
f"{name} (fixed)" if name in self.fixed_params else name
for name in (self.data.param_names or [])
]
# Dimensions
self.nobs = self.model.nobs
self.k_params = self.model.k_params
self._rank = None
@cache_readonly
def nobs_effective(self):
"""The number of observations used to fit the model, after adjustments"""
raise NotImplementedError
@cache_readonly
def df_resid(self):
"""The residual degrees of freedom"""
return self.nobs_effective - self.df_model
@cache_readonly
def aic(self):
"""The Akaike Information Criterion"""
return aic(self.llf, self.nobs_effective, self.df_model)
@cache_readonly
def aicc(self):
"""The Akaike Information Criterion with small sample correction"""
return aicc(self.llf, self.nobs_effective, self.df_model)
@cache_readonly
def bic(self):
"""The Bayes Information Criterion"""
return bic(self.llf, self.nobs_effective, self.df_model)
@cache_readonly
def fittedvalues(self):
"""The predicted values of the model"""
# TODO
raise NotImplementedError
@cache_readonly
def hqic(self):
"""The Hannan-Quinn Information Criterion"""
# return (-2 * self.llf +
# 2 * np.log(np.log(self.nobs_effective)) * self.df_model)
return hqic(self.llf, self.nobs_effective, self.df_model)
@cache_readonly
def llf(self):
"""The value of the log-likelihood function evaluated at `params`"""
raise NotImplementedError
@cache_readonly
def mae(self):
"""The mean absolute error"""
return np.mean(np.abs(self.resid))
@cache_readonly
def mse(self):
"""The mean squared error"""
return self.sse / self.nobs
@cache_readonly
def pvalues(self):
"""
The p-values associated with the z-statistics of the coefficients.
Note that the coefficients are assumed to have a Normal distribution.
"""
pvalues = np.zeros_like(self.zvalues) * np.nan
mask = np.ones_like(pvalues, dtype=bool)
mask[self._free_params_index] = True
mask &= ~np.isnan(self.zvalues)
pvalues[mask] = stats.norm.sf(np.abs(self.zvalues[mask])) * 2
return pvalues
@cache_readonly
def resid(self):
"""The model residuals"""
raise NotImplementedError
@cache_readonly
def sse(self):
"""The sum of squared errors"""
return np.sum(self.resid ** 2)
@cache_readonly
def zvalues(self):
"""The z-statistics for the coefficients"""
return self.params / self.bse
def _get_prediction_start_index(self, anchor):
"""Returns a valid numeric start index for predictions/simulations"""
if anchor is None or anchor == "start":
iloc = 0
elif anchor == "end":
iloc = self.nobs
else:
iloc, _, _ = self.model._get_index_loc(anchor)
if isinstance(iloc, slice):
iloc = iloc.start
iloc += 1 # anchor is one before start of prediction/simulation
if iloc < 0:
iloc = self.nobs + iloc
if iloc > self.nobs:
raise ValueError("Cannot anchor simulation outside of the sample.")
return iloc
def _cov_params_approx(
self, approx_complex_step=True, approx_centered=False
):
evaluated_hessian = self.nobs_effective * self.model.hessian(
params=self.params,
transformed=True,
includes_fixed=True,
method="approx",
approx_complex_step=approx_complex_step,
approx_centered=approx_centered,
)
# TODO: Case with "not approx_complex_step" is not hit in
# tests as of 2017-05-19
if len(self.fixed_params) > 0:
mask = np.ix_(self._free_params_index, self._free_params_index)
if len(self.fixed_params) < self.k_params:
(tmp, singular_values) = pinv_extended(evaluated_hessian[mask])
else:
tmp, singular_values = np.nan, [np.nan]
neg_cov = np.zeros_like(evaluated_hessian) * np.nan
neg_cov[mask] = tmp
else:
(neg_cov, singular_values) = pinv_extended(evaluated_hessian)
self.model.update(self.params, transformed=True, includes_fixed=True)
if self._rank is None:
self._rank = np.linalg.matrix_rank(np.diag(singular_values))
return -neg_cov
@cache_readonly
def cov_params_approx(self):
"""
The variance / covariance matrix, computed using the numerical
Hessian approximated by complex step or finite differences methods.
"""
return self._cov_params_approx(
self._cov_approx_complex_step, self._cov_approx_centered
)
def test_serial_correlation(self, method, lags=None):
"""
Ljung-Box test for no serial correlation of standardized residuals
Null hypothesis is no serial correlation.
Parameters
----------
method : {'ljungbox', 'boxpierce', None}
The statistical test for serial correlation. If None, an attempt is
made to select an appropriate test.
lags : int, array_like of int, or None, optional
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 min(10, nobs//5) for
non-seasonal time series and min (2*m, nobs//5) for seasonal time
series.
Returns
-------
output : ndarray
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
Ljung-Box test for serial correlation.
Notes
-----
For statespace models: let `d` = max(loglikelihood_burn, nobs_diffuse);
this test is calculated ignoring the first `d` residuals.
Output is nan for any endogenous variable which has missing values.
"""
if method is None:
method = "ljungbox"
method = string_like(
method, "method", options=("ljungbox", "boxpierce")
)
if self.standardized_forecasts_error is None:
raise ValueError("Cannot compute test statistic when standardized"
" forecast errors have not been computed.")
from statsmodels.stats.diagnostic import acorr_ljungbox
if hasattr(self, "loglikelihood_burn"):
d = np.maximum(self.loglikelihood_burn, self.nobs_diffuse)
# This differs from self.nobs_effective because here we want to
# exclude exact diffuse periods, whereas self.nobs_effective
# only excludes explicitly burned (usually approximate diffuse)
# periods.
nobs_effective = self.nobs - d
else:
nobs_effective = self.nobs_effective
output = []
# Default lags for acorr_ljungbox is 40, but may not always have
# that many observations
if lags is None:
seasonal_periods = getattr(self.model, "seasonal_periods", 0)
if seasonal_periods:
lags = min(2 * seasonal_periods, nobs_effective // 5)
else:
lags = min(10, nobs_effective // 5)
cols = [2, 3] if method == "boxpierce" else [0, 1]
for i in range(self.model.k_endog):
if hasattr(self, "filter_results"):
x = self.filter_results.standardized_forecasts_error[i][d:]
else:
x = self.standardized_forecasts_error
results = acorr_ljungbox(
x, lags=lags, boxpierce=(method == "boxpierce")
)
output.append(np.asarray(results)[:, cols].T)
output = np.c_[output]
return output
def test_heteroskedasticity(self, method, alternative="two-sided",
use_f=True):
r"""
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. The null hypothesis
is of no heteroskedasticity.
Parameters
----------
method : {'breakvar', 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 : {'increasing', 'decreasing', 'two-sided'}, optional
This specifies the alternative for the p-value calculation. Default
is two-sided.
use_f : bool, 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 : ndarray
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.
See Also
--------
statsmodels.tsa.stattools.breakvar_heteroskedasticity_test
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 :math:`h = [T/3]`, the test statistic is:
.. math::
H(h) = \sum_{t=T-h+1}^T \tilde v_t^2
\Bigg / \sum_{t=d+1}^{d+1+h} \tilde v_t^2
where :math:`d` = max(loglikelihood_burn, nobs_diffuse) (usually
corresponding to diffuse initialization under either the approximate
or exact approach).
This statistic can be tested against an :math:`F(h,h)` distribution.
Alternatively, :math:`h H(h)` is asymptotically distributed according
to :math:`\chi_h^2`; this second test can be applied by passing
`asymptotic=True` as an argument.
See section 5.4 of [1]_ for the above formula and discussion, as well
as additional details.
TODO
- Allow specification of :math:`h`
References
----------
.. [1] Harvey, Andrew C. 1990. *Forecasting, Structural Time Series*
*Models and the Kalman Filter.* Cambridge University Press.
"""
if method is None:
method = "breakvar"
_ = string_like(method, "method", options=("breakvar",))
if self.standardized_forecasts_error is None:
raise ValueError("Cannot compute test statistic when standardized"
" forecast errors have not been computed.")
# Store some values
if hasattr(self, "filter_results"):
resid = self.filter_results.standardized_forecasts_error
d = np.maximum(self.loglikelihood_burn, self.nobs_diffuse)
# This differs from self.nobs_effective because here we want to
# exclude exact diffuse periods, whereas self.nobs_effective
# only excludes explicitly burned (usually approximate diffuse)
# periods.
nobs_effective = self.nobs - d
else:
resid = np.asarray(self.standardized_forecasts_error)
if resid.ndim == 1:
resid = resid[np.newaxis, :]
nobs_effective = self.nobs_effective
d = 0
h = int(np.round(nobs_effective / 3))
test_statistics = []
p_values = []
for i in range(self.model.k_endog):
_het_result = breakvar_heteroskedasticity_test(
resid[i, d:], subset_length=h, alternative=alternative,
use_f=use_f
)
test_statistics.append(_het_result.statistic)
p_values.append(_het_result.pvalue)
output = np.c_[test_statistics, p_values]
return output
def test_normality(self, method):
"""
Test for normality of standardized residuals
Null hypothesis is normality.
Parameters
----------
method : {'jarquebera', 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.
Returns
-------
output : ndarray
An array with `(Jarque-Bera, p-value, skew, kurtosis)` for each
endogenous variable, sized `(k_endog, 4)`.
See Also
--------
statsmodels.stats.stattools.jarque_bera
The Jarque-Bera test of normality.
Notes
-----
For statespace models: let `d` = max(loglikelihood_burn, nobs_diffuse);
this test is calculated ignoring the first `d` residuals.
In the case of missing data, the maintained hypothesis is that the
data are missing completely at random. This test is then run on the
standardized residuals excluding those corresponding to missing
observations.
"""
if method is None:
method = "jarquebera"
_ = string_like(method, "method", options=("jarquebera",))
if self.standardized_forecasts_error is None:
raise ValueError("Cannot compute test statistic when standardized"
" forecast errors have not been computed.")
from statsmodels.stats.stattools import jarque_bera
if hasattr(self, "loglikelihood_burn"):
d = np.maximum(self.loglikelihood_burn, self.nobs_diffuse)
else:
d = 0
output = []
for i in range(self.model.k_endog):
if hasattr(self, "filter_results"):
resid = self.filter_results.standardized_forecasts_error[
i, d:
]
else:
resid = self.standardized_forecasts_error
mask = ~np.isnan(resid)
output.append(jarque_bera(resid[mask]))
return np.array(output)
def summary(
self,
alpha=0.05,
start=None,
title=None,
model_name=None,
display_params=True,
):
"""
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.
title : str, optional
The title used for the summary table.
model_name : str, optional
The name of the model used. Default is to use model class name.
display_params : bool, optional
Whether or not to display the parameters table. Default is True.
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
"""
from statsmodels.iolib.summary import Summary
# Model specification results
model = self.model
if title is None:
title = "Statespace Model Results"
if start is None:
start = 0
if self.model._index_dates:
ix = self.model._index
d = ix[start]
sample = [f"{d.month:02d}-{d.day:02d}-{d.year:02d}"]
d = ix[-1]
sample += ["- " + f"{d.month:02d}-{d.day:02d}-{d.year:02d}"]
else:
sample = [str(start), " - " + str(self.nobs)]
# Standardize the model name as a list of str
if model_name is None:
model_name = model.__class__.__name__
# Diagnostic tests results
try:
het = self.test_heteroskedasticity(method="breakvar")
except Exception: # FIXME: catch something specific
het = np.array([[np.nan] * 2])
try:
lb = self.test_serial_correlation(method="ljungbox")
except Exception: # FIXME: catch something specific
lb = np.array([[np.nan] * 2]).reshape(1, 2, 1)
try:
jb = self.test_normality(method="jarquebera")
except Exception: # FIXME: catch something specific
jb = np.array([[np.nan] * 4])
# Create the tables
if not isinstance(model_name, list):
model_name = [model_name]
top_left = [("Dep. Variable:", None), ("Model:", [model_name[0]])]
top_left.extend(
("", ["+ " + model_name[i]]) for i in range(1, len(model_name))
)
top_left += [
("Date:", None),
("Time:", None),
("Sample:", [sample[0]]),
("", [sample[1]]),
]
top_right = [
("No. Observations:", [self.nobs]),
("Log Likelihood", [f"{self.llf:#5.3f}"]),
]
if hasattr(self, "rsquared"):
top_right.append(("R-squared:", [f"{self.rsquared:#8.3f}"]))
top_right += [
("AIC", [f"{self.aic:#5.3f}"]),
("BIC", [f"{self.bic:#5.3f}"]),
("HQIC", [f"{self.hqic:#5.3f}"]),
]
if hasattr(self, "filter_results"):
if (
self.filter_results is not None
and self.filter_results.filter_concentrated
):
top_right.append(("Scale", [f"{self.scale:#5.3f}"]))
else:
top_right.append(("Scale", [f"{self.scale:#5.3f}"]))
if hasattr(self, "cov_type"):
top_left.append(("Covariance Type:", [self.cov_type]))
def format_str(array):
return [", ".join([f"{i:.2f}" for i in array])]
diagn_left = [
("Ljung-Box (Q):", format_str(lb[:, 0, -1])),
("Prob(Q):", format_str(lb[:, 1, -1])),
("Heteroskedasticity (H):", format_str(het[:, 0])),
("Prob(H) (two-sided):", format_str(het[:, 1])),
]
diagn_right = [
("Jarque-Bera (JB):", format_str(jb[:, 0])),
("Prob(JB):", format_str(jb[:, 1])),
("Skew:", format_str(jb[:, 2])),
("Kurtosis:", format_str(jb[:, 3])),
]
summary = Summary()
summary.add_table_2cols(
self, gleft=top_left, gright=top_right, title=title
)
if len(self.params) > 0 and display_params:
summary.add_table_params(
self, alpha=alpha, xname=self.param_names, use_t=False
)
summary.add_table_2cols(
self, gleft=diagn_left, gright=diagn_right, title=""
)
# Add warnings/notes, added to text format only
etext = []
if hasattr(self, "cov_type") and "description" in self.cov_kwds:
etext.append(self.cov_kwds["description"])
if self._rank < (len(self.params) - len(self.fixed_params)):
cov_params = self.cov_params()
if len(self.fixed_params) > 0:
mask = np.ix_(self._free_params_index, self._free_params_index)
cov_params = cov_params[mask]
etext.append(
"Covariance matrix is singular or near-singular,"
f" with condition number {_safe_cond(cov_params):6.3g}. Standard errors may be"
" unstable."
)
if etext:
etext = [
f"[{i + 1}] {text}" for i, text in enumerate(etext)
]
etext.insert(0, "Warnings:")
summary.add_extra_txt(etext)
return summary
