WhittleFactory

class WhittleFactory(*args)

Whittle estimator of a scalar ARMA Gaussian process.

Available constructors:

WhittleFactory()

WhittleFactory(p, q, invert)

WhittleFactory(indP, indQ, invertible)

Parameters:
pint

Order of the AR part of the ARMA(p,q) process of dimension d.

qint

Order of the MA part of the ARMA(p,q) process of dimension d.

invertiblebool, optional

Restrict the estimation to invertible ARMA processes.

By default: True.

indPIndices

All the p orders that will be investigated. Care: not yet implemented.

indQIndices

All the p orders that will be investigated. Care: not yet implemented.

Notes

We suppose here that the white noise is normal with zero mean and variance \sigma^2. It implies that the ARMA process estimated is normal.

For each order (p,q), the estimation of the coefficients (a_k)_{1\leq k\leq p}, (b_k)_{1\leq k\leq q} and the variance \sigma^2 is done using the Whittle estimator which is based on the maximization of the likelihood function in the frequency domain.

The principle is detailed hereafter for the case of a time series : in the case of a process sample, the estimator is similar except for the periodogram which is computed differently.

Let (t_i, \vect{X}_i)_{0\leq i \leq n-1} be a multivariate time series of dimension d from an ARMA(p,q) process.

The spectral density function of the ARMA(p,q) process writes :

f(\lambda, \vect{\theta}, \sigma^2) = \frac{\sigma^2}{2 \pi} \frac{|1 + b_1 \exp(-i \lambda) + \ldots + b_q \exp(-i q \lambda)|^2}{|1 + a_1 \exp(-i \lambda) + \ldots + a_p \exp(-i p \lambda)|^2} = \frac{\sigma^2}{2 \pi} g(\lambda, \vect{\theta})

where \vect{\theta} = (a_1, a_2,...,a_p,b_1,...,b_q) and \lambda is the frequency value.

The Whittle log-likelihood writes :

\log L_w(\vect{\theta}, \sigma^2) = - \sum_{j=1}^{m} \log f(\lambda_j,  \vect{\theta}, \sigma^2) - \frac{1}{2 \pi} \sum_{j=1}^{m} \frac{I(\lambda_j)}{f(\lambda_j,  \vect{\theta}, \sigma^2)}

where :

  • I is the non parametric estimate of the spectral density, expressed in the Fourier space (frequencies in [0,2\pi] instead of [-f_{max}, f_{max}]). OpenTURNS uses by default the Welch estimator.

  • \lambda_j is the j-th Fourier frequency, \lambda_j = \frac{2 \pi j}{n}, j=1 \ldots m with m the largest integer \leq \frac{n-1}{2}.

We estimate the (p+q+1) scalar coefficients by maximizing the log-likelihood function. The corresponding equations lead to the following relation :

\sigma^{2*} = \frac{1}{m} \sum_{j=1}^{m} \frac{I(\lambda_j)}{g(\lambda_j, \vect{\theta}^{*})}

where \vect{\theta}^{*} maximizes :

\log L_w(\vect{\theta}) = m (\log(2 \pi) - 1) - m \log\left( \frac{1}{m} \sum_{j=1}^{m} \frac{I(\lambda_j)}{g(\lambda_j, \vect{\theta})}\right) - \sum_{j=1}^{m} g(\lambda_j, \vect{\theta})

The Whitle estimation requires that :

  • the determinant of the eigenvalues of the companion matrix associated to the polynomial 1 + a_1X + \dots + a_pX^p are outside the unit disc, which guarantees the stationarity of the process;

  • the determinant of the eigenvalues of the companion matrix associated to the polynomial 1 + ba_1X + \dots + b_qX^q are outside the unit disc, which guarantees the invertibility of the process.

The criteria AIC, AIC_c (corrected AIC) and BIC are evaluated to help the model selection:

\begin{eqnarray*}
    AIC_c  & = &  -2\log L_w + 2(p + q + 1)\frac{m}{m - p - q - 2}\\
    AIC & = & -2\log L_w + 2(p + q + 1)\\
    BIC & = & -2\log L_w + 2(p + q + 1)\log(p + q + 1)
\end{eqnarray*}

where m is half the number of points of the time grid of the process sample (if the data are a process sample) or in a block of the time series (if the data are a time series).

The BIC criterion leads to a model that gives a better prediction. The AIC criterion selects the best model that fits the given data. The AIC_c criterion improves the previous one by penalizing a too high order that would artificially fit to the data.

Examples

Create a time series from a scalar ARMA(4,2) and a normal white noise:

>>> import openturns as ot
>>> myTimeGrid = ot.RegularGrid(0.0, 0.1, 100)
>>> myWhiteNoise = ot.WhiteNoise(ot.Triangular(-1.0, 0.0, 1.0), myTimeGrid)
>>> myARCoef = ot.ARMACoefficients([0.4, 0.3, 0.2, 0.1])
>>> myMACoef = ot.ARMACoefficients([0.4, 0.3])
>>> myARMAProcess = ot.ARMA(myARCoef, myMACoef, myWhiteNoise)
>>> myTimeSeries = myARMAProcess.getRealization()
>>> myProcessSample = myARMAProcess.getSample(10)

Estimate the ARMA process specifying the orders:

>>> myFactory_42 = ot.WhittleFactory(4, 2)

Check the default SpectralModelFactory:

>>> print(myFactory_42.getSpectralModelFactory())  

Set a particular spectral model: WelchFactory as SpectralModelFactory with the Hann filtering window:

>>> myFilteringWindow = ot.Hann()
>>> mySpectralFactory = ot.WelchFactory(myFilteringWindow, 4, 0)
>>> myFactory_42.setSpectralModelFactory(mySpectralFactory)
>>> print(myFactory_42.getSpectralModelFactory())  

Estimate the ARMA process specifying a range for the orders:

p = [1, 2, 4] and q = [4,5,6]:

>>> pIndices = [1, 2, 4]
>>> qIndices = [4, 5, 6]
>>> myFactory_Range = ot.WhittleFactory(pIndices, qIndices)

To get the quantified AICc, AIC and BIC criteria:

>>> myARMA_42, myCriterion = myFactory_42.buildWithCriteria(ot.TimeSeries(myTimeSeries))
>>> AICc, AIC, BIC = myCriterion[0:3]

Methods

build(*args)

Estimate the ARMA process.

buildWithCriteria(*args)

Estimate the ARMA process.

clearHistory()

Clear the history of the factory.

disableHistory()

Deactivate the history of all the estimated models.

enableHistory()

Activate the history of all the estimated models.

getClassName()

Accessor to the object's name.

getCurrentP()

Accessor to the current P order.

getCurrentQ()

Accessor to the current Q order.

getHistory()

Check whether the history mechanism is activated.

getInvertible()

Accessor to the invertible constraint.

getName()

Accessor to the object's name.

getP()

Accessor to the P orders.

getQ()

Accessor to the Q orders.

getSpectralModelFactory()

Accessor to the spectral factory.

getStartingPoints()

Accessor to the starting points for the optimization step.

hasName()

Test if the object is named.

isHistoryEnabled()

Check whether the history mechanism is activated.

setInvertible(invertible)

Accessor to the invertible constraint.

setName(name)

Accessor to the object's name.

setSpectralModelFactory(factory)

Accessor to the spectral factory.

setStartingPoints(startingPoints)

Accessor to the starting points for the optimization step.

__init__(*args)
build(*args)

Estimate the ARMA process.

Available usages:

build(myTimeSeries)

build(myProcessSample)

Parameters:
myTimeSeriesTimeSeries

One realization of the process.

myProcessSampleProcessSample

Several realizations of the process.

Returns:
myARMAARMA

The process estimated with the Whittle estimator.

Notes

The model selection is made using the spectral density built using the given data and theoretical spectral density of the ARMA process.

The best ARMA process is selected according to the corrected AIC criterion.

buildWithCriteria(*args)

Estimate the ARMA process.

Available usages:

buildWithCriteria(myTimeSeries)

buildWithCriteria(myProcessSample)

Parameters:
myTimeSeriesTimeSeries

One realization of the process.

myProcessSampleProcessSample

Several realizations of the process.

Returns:
myARMAARMA

The process estimated with the Whittle estimator.

criterionPoint

Result of the evaluation of the AICc, AIC and BIC criteria

Notes

The model selection is made using the spectral density built using the given data and theoretical spectral density of the ARMA process.

The best ARMA process is selected according to the corrected AIC criterion.

clearHistory()

Clear the history of the factory.

Notes

Clear the history of the factory.

disableHistory()

Deactivate the history of all the estimated models.

Notes

Deactivate the history mechanism which is the trace of all the tested models and their associated information criteria.

enableHistory()

Activate the history of all the estimated models.

Notes

Activate the history mechanism which is the trace of all the tested models and their associated information criteria.

By default, the history mechanism is activated.

getClassName()

Accessor to the object’s name.

Returns:
class_namestr

The object class name (object.__class__.__name__).

getCurrentP()

Accessor to the current P order.

Returns:
pint

Order of the AR part of the ARMA(p,q) process of dimension d.

getCurrentQ()

Accessor to the current Q order.

Returns:
qint

Order of the MA part of the ARMA(p,q) process of dimension d.

getHistory()

Check whether the history mechanism is activated.

Returns:
histMeca list of WhittleFactoryState

Returns the collection of all the states that have been built for the estimation.

getInvertible()

Accessor to the invertible constraint.

Returns:
invertiblebool

The initial AR coefficients used for the optimization algorithm.

getName()

Accessor to the object’s name.

Returns:
namestr

The name of the object.

getP()

Accessor to the P orders.

Returns:
indPIndices

All the p orders that will be investigated.

getQ()

Accessor to the Q orders.

Returns:
indQIndices

All the p orders that will be investigated.

getSpectralModelFactory()

Accessor to the spectral factory.

Returns:
initARCoeffSpectralModelFactory

The spectral factory used to estimate the spectral density based on the data.

getStartingPoints()

Accessor to the starting points for the optimization step.

Returns:
startPointsLista list of Point

Starting points for the optimization step, for each pair of orders that will be tested.

hasName()

Test if the object is named.

Returns:
hasNamebool

True if the name is not empty.

isHistoryEnabled()

Check whether the history mechanism is activated.

Returns:
histMecbool

Check whether the history mechanism is activated.

By default, the history mechanism is activated.

setInvertible(invertible)

Accessor to the invertible constraint.

Parameters:
invertiblebool

The initial AR coefficients used for the optimization algorithm.

setName(name)

Accessor to the object’s name.

Parameters:
namestr

The name of the object.

setSpectralModelFactory(factory)

Accessor to the spectral factory.

Parameters:
spectralModelFactSpectralModelFactory

The spectral factory used to estimate the spectral density based on the data.

setStartingPoints(startingPoints)

Accessor to the starting points for the optimization step.

Parameters:
startPointsLista list of Point

Starting points for the optimization step, for each pair of orders that will be tested.

Examples using the class

Estimate a scalar ARMA process

Estimate a scalar ARMA process