Parametric Estimators for Long Memory

ParametricEstimators

This module contains functions for estimating the parameters of the fractional differenced process à la ARFIMA and the CSA process.

Author

J. Eduardo Vera-Valdés

ar_var_vals(T::Int, AR::AbstractArray)

Computes the autocovariance function of the stationary AR(p) process with unit innovation variance at lags 0, 1, ..., T-1, solving the Yule-Walker equations. Internal function used by arfima_var_vals.

arfima_cor_vals(T::Int, d::Real, AR::AbstractArray, MA::AbstractArray)

Computes the autocorrelation function of the stationary ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter, d∈(-1/2,1/2).
• AR::AbstractArray: The AR coefficients, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ.
• MA::AbstractArray: The MA coefficients, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.

Output

• cors::Array: The autocorrelation function of the ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Notes

This function uses arfimavarvals() to compute the autocovariance function and then normalizes by the variance (first computed value).

Examples

julia> arfima_cor_vals(10, 0.3, [0.5], [0.2])

arfima_exact_gen(T::Int, d::Real, AR::AbstractArray, MA::AbstractArray; μ=0, σ=1)

Generate a stationary ARFIMA(p,d,q) process of length T with exactly the theoretical autocovariances, using the Davies-Harte method.

Arguments

• T::Int: The length of the time series.
• d::Real: The fractional differencing parameter, d∈(-1/2,1/2).
• AR::AbstractArray: The AR coefficients, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ, which must be stationary.
• MA::AbstractArray: The MA coefficients, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.

Optional arguments

• μ::Real: The mean of the time series.
• σ::Real: The standard deviation of the innovations.

Output

• x::Array: The time series.

Notes

In contrast to arfima_gen, which applies the fractional filter truncated at the start of the sample (a nonstationary, so-called type II, process), this function draws from the stationary process, with the autocovariances of arfima_var_vals. Note the order of the arguments, (d, AR, MA) as in arfima_var_vals and arfima_forecast, while arfima_gen takes (μ, AR, d, MA). An error is thrown if the circulant embedding is not nonnegative definite. See davies_harte for details.

Examples

julia> arfima_exact_gen(100, 0.3, [0.5], [0.2])

arfima_llk(θ::AbstractArray, x::AbstractArray, p::Int, q::Int)

Computes the log-likelihood of the ARFIMA(p,d,q) process given the data x.

Arguments

• θ::AbstractArray: The 1 + p + q parameters on the unconstrained scale: d on the logit scale, mapped to -1/2 + exp(θ₁)/(1+exp(θ₁)) ∈ (-1/2, 1/2), followed by the partial autocorrelations of the AR and MA polynomials on the atanh scale (Monahan, 1984).
• x::AbstractArray: The data (zero mean).
• p::Int: The order of the autoregressive polynomial.
• q::Int: The order of the moving average polynomial.

Output

• llk::Real: The negative concentrated log-likelihood (up to constants, divided by T), suitable for minimization.

Notes

This function computes the exact concentrated log-likelihood function of the ARFIMA(p,d,q) process, using the autocovariances of arfima_var_vals and the Durbin-Levinson algorithm, see durbin_levinson. The parametrization imposes stationarity and invertibility.

Examples

julia> arfima_llk([0.4, 0.2, -0.1], randn(100), 1, 1)

arfima_mle_est(x::AbstractArray, p::Int, q::Int)

Computes the exact maximum likelihood estimates of the parameters of the ARFIMA(p,d,q) model given the data x.

Arguments

• x::AbstractArray: The data.
• p::Int: The order of the autoregressive polynomial.
• q::Int: The order of the moving average polynomial.

Output

• d::Real: The maximum likelihood estimate of the fractional differencing parameter.
• AR::Array: The maximum likelihood estimates of the AR coefficients φ, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ.
• MA::Array: The maximum likelihood estimates of the MA coefficients θ, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.
• σ::Real: The maximum likelihood estimate of the standard deviation of the innovations.

Notes

The model is φ(L)(1-L)^d x_t = θ(L)ε_t, the same parametrization used by arfima_gen and arfima_whittle_est. The series is demeaned with the sample mean, and the Gaussian likelihood is evaluated exactly with the autocovariances of arfima_var_vals and the Durbin-Levinson algorithm (Sowell, 1992). Stationarity and invertibility are imposed by optimizing over d∈(-1/2,1/2) and over the partial autocorrelations of the AR and MA polynomials (Monahan, 1984). Hence, nonstationary data should be differenced before estimation. This function uses the Optim package to minimize the log-likelihood function, started at the Whittle estimates of arfima_whittle_est. With p = q = 0 it estimates the same model as fi_mle_est. As usual for the exact maximum likelihood estimator with an estimated mean, d is biased downward in small samples, more so with AR terms that compete with d for the low-frequency dynamics (in simulations of an ARFIMA(1,0.3,0) with φ₁ = 0.5, the mean estimate of d is about 0.19 for T = 500 and 0.27 for T = 2000). In such cases the likelihood can also have several local maxima. The Whittle estimator arfima_whittle_est does not depend on the mean and has a smaller bias.

References

Sowell, F. (1992). Maximum likelihood estimation of stationary univariate fractionally integrated time series models. Journal of Econometrics, 53(1-3), 165-188.

Examples

julia> arfima_mle_est(arfima_gen(500, 0, [0.5], 0.3, []), 1, 0)

arfima_params(θ::AbstractArray, p::Int, q::Int)

Maps the unconstrained parameters θ to d∈(-1/2,1/2), stationary AR coefficients and invertible MA coefficients. Internal function used by arfima_whittle_est.

arfima_spec(w::AbstractArray, d::Real, AR::AbstractArray, MA::AbstractArray)

Computes the spectral density of the ARFIMA(p,d,q) model at the frequencies w, up to the factor σ²/(2π). Internal function used by arfima_whittle_est.

arfima_var_matrix(T::Int, d::Real, AR::AbstractArray, MA::AbstractArray)

Constructs the autocovariance matrix of the stationary ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter, d∈(-1/2,1/2).
• AR::AbstractArray: The AR coefficients, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ.
• MA::AbstractArray: The MA coefficients, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.

Output

• V::Array: The autocovariance matrix of the ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Examples

julia> arfima_var_matrix(10, 0.3, [0.5], [0.2])

arfima_var_vals(T::Int, d::Real, AR::AbstractArray, MA::AbstractArray; σ=1)

Computes the autocovariance function of the stationary ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter, d∈(-1/2,1/2).
• AR::AbstractArray: The AR coefficients φ, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ, which must be stationary.
• MA::AbstractArray: The MA coefficients θ, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.

Optional arguments

• σ::Real: The standard deviation of the innovations.

Output

• acv::Array: The autocovariance function of the ARFIMA(p,d,q) process at lags 0, 1, ..., T-1.

Notes

The model is φ(L)(1-L)^d x_t = θ(L)ε_t, the same parametrization used by arfima_gen, arfima_whittle_est and arfima_forecast. The autocovariances are computed in two steps. First, the autocovariances of the fractionally integrated MA process u_t = θ(L)(1-L)^{-d}ε_t are a finite combination of those of fi_var_vals. Second, the autocovariances of x_t = φ(L)^{-1}u_t are γ_x(h) = Σ_k g(k) γ_u(h - k), where g are the autocovariances of the AR(p) process, which decay geometrically. The sum is truncated when the AR autocovariances are below 1e-15 relative to their variance, so the result is exact up to rounding error. This is equivalent to the hypergeometric representation of Sowell (1992), but it does not require the roots of the AR polynomial to be distinct. The computing time grows as the roots of the AR polynomial approach the unit circle.

References

Sowell, F. (1992). Maximum likelihood estimation of stationary univariate fractionally integrated time series models. Journal of Econometrics, 53(1-3), 165-188.

Examples

julia> arfima_var_vals(10, 0.3, [0.5], [0.2])

arfima_whittle_est(x::AbstractArray, p::Int, q::Int)

Estimates the parameters of the ARFIMA(p,d,q) model by minimizing the Whittle approximation to the Gaussian likelihood in the frequency domain.

Arguments

• x::AbstractArray: The data.
• p::Int: The order of the autoregressive polynomial.
• q::Int: The order of the moving average polynomial.

Output

• d::Real: The estimate of the fractional differencing parameter.
• AR::Array: The estimates of the AR coefficients φ, with the AR polynomial 1 - φ₁L - ... - φₚLᵖ.
• MA::Array: The estimates of the MA coefficients θ, with the MA polynomial 1 + θ₁L + ... + θ_qL^q.
• σ::Real: The estimate of the standard deviation of the innovations.

Notes

The model is φ(L)(1-L)^d x_t = θ(L)ε_t, the same parametrization used by arfima_gen. The objective is the concentrated Whittle function log(mean(I_j/g_j)) + mean(log(g_j)) over the Fourier frequencies 2πj/T, j = 1, ..., floor((T-1)/2), where g is the spectral density up to the factor σ²/(2π). Stationarity and invertibility are imposed by optimizing over d∈(-1/2,1/2) and over the partial autocorrelations of the AR and MA polynomials (Monahan, 1984). Hence, nonstationary data should be differenced before estimation. The estimator does not depend on the mean of x. It uses the Optim package, started at the log-periodogram estimate of d and zero ARMA coefficients.

References

Fox, R., and Taqqu, M. S. (1986). Large-sample properties of parameter estimates for strongly dependent stationary Gaussian time series. The Annals of Statistics, 14(2), 517-532.

Examples

julia> arfima_whittle_est(arfima_gen(500, 0, [0.5], 0.3, []), 1, 0)

arfima_whittle_llk(θ::AbstractArray, I_w::AbstractArray, w::AbstractArray, p::Int, q::Int)

Computes the concentrated Whittle objective of the ARFIMA(p,d,q) model for the unconstrained parameters θ. Internal function used by arfima_whittle_est.

coefs_to_pacf(φ::AbstractArray)

Maps the coefficients φ of a stationary AR polynomial 1 - φ₁z - ... - φₚzᵖ to its partial autocorrelations, the inverse of pacf_to_coefs. Internal function used by arfima_mle_est.

csa_cor_vals(T::Int, p::Real, q::Real)

Computes the autocorrelation function of the CSA process with parameters p and q at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• p::Real: The first parameter of the CSA process.
• q::Real: The second parameter of the CSA process.

Output

• acf::Array: The autocorrelation function of the CSA process with parameters p and q at lags 0, 1, ..., T-1.

Notes

This function uses csavarvals() to compute the autocovariance function and then normalizes by the variance (first computed value).

Examples

julia> csa_cor_vals(20, 1.4, 1.6)

csa_exact_gen(T::Int, p::Real, q::Real; μ=0, σ=1)

Generate a time series of length T with exactly the autocovariances of the cross-sectional aggregated (CSA) process with parameters p and q, using the Davies-Harte method.

Arguments

• T::Int: The length of the time series.
• p::Real: The first parameter of the CSA process.
• q::Real: The second parameter of the CSA process, q > 1.

Optional arguments

• μ::Real: The mean of the time series.
• σ::Real: The standard deviation of the innovations of the individual AR(1) processes.

Output

• x::Array: The time series.

Notes

The autocovariances are those of csa_var_vals, i.e. of the limit of the cross-sectional aggregation of AR(1) processes with squared parameters drawn from a Beta(p,q) distribution, normalized by the number of units. An error is thrown if the circulant embedding is not nonnegative definite. See davies_harte for details.

Examples

julia> csa_exact_gen(100, 1.4, 1.6)

csa_llk(p::Real, q::Real, x::AbstractArray)

Computes the log-likelihood of the CSA process with parameters p and q given the data x.

Arguments

• p::Real: The first parameter of the CSA process on the unconstrained scale; it is mapped to 1 + 2exp(p)/(1+exp(p)) ∈ (1, 3).
• q::Real: The second parameter of the CSA process on the unconstrained scale; it is mapped to 1 + 2exp(q)/(1+exp(q)) ∈ (1, 3).
• x::Array: The data.

Output

• llk::Real: The negative concentrated log-likelihood (up to constants, divided by T), suitable for minimization.

Notes

This function computes the concentrated log-likelihood function of the CSA process with parameters p and q given the data x.

Examples

julia> csa_llk(1.4, 1.8, randn(100,1))

csa_mle_est(x::AbstractArray)

Computes the maximum likelihood estimate of the parameters p and q of the CSA process and the standard deviation of the CSA process given the data x.

Arguments

• x::Array: The data.

Output

• p::Real: The maximum likelihood estimate of the first parameter of the CSA process.
• q::Real: The maximum likelihood estimate of the second parameter of the CSA process.
• σ::Real: The maximum likelihood estimate of the standard deviation of the CSA process.

Notes

This function uses the Optim package to minimize the log-likelihood function. The parameters are restricted to p∈(1,3) and q∈(1,3). The starting values are p = 2 and q = 2(1-d), with d the log-periodogram estimate, so the estimates are reproducible. The likelihood is evaluated with the Durbin-Levinson algorithm, see durbin_levinson.

Examples

julia> csa_mle_est(randn(100,1))

csa_var_matrix(T::Int, p::Real, q::Real)

Constructs the autocovariance matrix of the CSA process with parametersp and q at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• p::Real: The first parameter of the CSA process.
• q::Real: The second parameter of the CSA process.

Output

• V::Array: The autocovariance matrix of the CSA process with parameters p and q at lags 0, 1, ..., T-1.

Examples

julia> csa_var_matrix(10, 1.4, 1.8)

csa_var_vals(T::Int, p::Real, q::Real)

Computes the autocovariance function of the CSA process with parameters p and q at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• p::Real: The first parameter of the CSA process.
• q::Real: The second parameter of the CSA process.

Output

• acf::Array: The autocovariance function of the CSA process with parameters p and q at lags 0, 1, ..., T-1.

Notes

This function uses the recursive formula for the autocovariance function of the CSA process, γ(k) = B(p + k/2, q - 1) / B(p, q) for unit innovation variance, which requires q > 1.

Examples

julia> csa_var_vals(20, 1.4, 1.6)

durbin_levinson(acv::AbstractArray, x::AbstractArray)

Computes the one-step-ahead best linear predictions of a stationary time series x and their mean squared errors, given the autocovariances acv, using the Durbin-Levinson algorithm.

Arguments

• acv::AbstractArray: The autocovariances at lags 0, 1, ..., T-1 (at least length(x) values).
• x::AbstractArray: The data (zero mean).

Output

• xhat::Array: The predictions of x[t] given x[1], ..., x[t-1], for t = 1, ..., T (xhat[1] = 0).
• v::Array: The mean squared errors of the predictions.

Notes

The algorithm requires O(T^2) operations and O(T) memory, instead of the O(T^3) operations and O(T^2) memory needed to factorize the covariance matrix. The Gaussian log-likelihood follows from the prediction error decomposition: logdet(V) = sum(log.(v)) and x'V⁻¹x = sum((x - xhat).^2 ./ v), where V is the Toeplitz covariance matrix. See Brockwell and Davis (1991), Time Series: Theory and Methods, Section 5.2.

Examples

julia> durbin_levinson(fi_var_vals(100, 0.3), fi_gen(100, 0.3))

fi_cor_vals(T::Int,d::Real)

Computes the autocorrelation function of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter.

Output

• vars::Array: The autocorrelation function of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Notes

This function uses fivarvals() to compute the autocovariance function and then normalizes by the variance (first computed value).

Examples

julia> fi_cor_vals(10, 0.4)

fi_exact_gen(T::Int, d::Real; μ=0, σ=1)

Generate a stationary fractionally differenced process of length T with exactly the theoretical autocovariances, using the Davies-Harte method.

Arguments

• T::Int: The length of the time series.
• d::Real: The fractional differencing parameter, d∈(-1/2,1/2).

Optional arguments

• μ::Real: The mean of the time series.
• σ::Real: The standard deviation of the innovations.

Output

• x::Array: The time series.

Notes

In contrast to fi_gen, which applies the fractional filter truncated at the start of the sample (a nonstationary, so-called type II, process), this function draws from the stationary process. See davies_harte for details.

Examples

julia> fi_exact_gen(100, 0.4)

fi_llk(d::Real, x::AbstractArray)

Computes the log-likelihood of the fractional differenced process with parameter d given the data x.

Arguments

• d::Real: The fractional differencing parameter on the unconstrained (logit) scale; it is mapped to -1/2 + exp(d)/(1+exp(d)) ∈ (-1/2, 1/2).
• x::Array: The data.

Output

• llk::Real: The negative concentrated log-likelihood (up to constants, divided by T), suitable for minimization.

Notes

This function computes the concentrated log-likelihood function of the fractional differenced process with parameter d given the data x.

Examples

julia> fi_llk(0.4, randn(100,1))

fi_mle_est(x::AbstractArray)

Computes the maximum likelihood estimate of the fractional differencing parameter and the standard deviation of the fractional differenced process given the data x.

Arguments

• x::Array: The data.

Output

• d::Real: The maximum likelihood estimate of the fractional differencing parameter.
• σ::Real: The maximum likelihood estimate of the standard deviation of the fractional differenced process.

Notes

This function uses the Optim package to minimize the log-likelihood function. The likelihood is evaluated with the Durbin-Levinson algorithm, see durbin_levinson.

Examples

julia> fi_mle_est(randn(100,1))

fi_var_matrix(T::Int, d::Real)

Constructs the autocovariance matrix of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter.

Output

• V::Array: The autocovariance matrix of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Examples

julia> fi_var_matrix(10, 0.4)

fi_var_vals(T::Int,d::Real)

Computes the autocovariance function of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Arguments

• T::Int: The number of lags to compute.
• d::Real: The fractional differencing parameter.

Output

• vars::Array: The autocovariance function of the fractional differenced process with parameter d at lags 0, 1, ..., T-1.

Notes

This function uses the recursive formula for the autocovariance function of the fractional differenced process.

Examples

julia> fi_var_vals(10, 0.4)

har_est(x::AbstractArray; m::AbstractArray = [1 , 5 , 22])

Estimates the parameters of the Heterogenous Autoregressive (HAR) model given the data x. See Corsi (2009).

Arguments

• x::Array: The data.

Optional arguments

• m::Array: An array with the lags to use in the estimation. By default, the lags are 1, 5, and 22; as suggested by the original paper.

Output

• β::Array: The estimated parameters of the HAR model.
• σ::Real: The estimated standard deviation of the HAR model.

Examples

julia> har_est(randn(100,1))

my_toeplitz(coefs::AbstractArray)

Constructs a Toeplitz matrix from the given coefficients.

Arguments

• coefs::Array: An array of coefficients.

Output

• Toep::Array: The Toeplitz matrix constructed from the given coefficients.

Notes

This is an internal helper function; it is not exported.

Examples

julia> LongMemory.ParametricEstimators.my_toeplitz([1, 2, 3])

pacf_to_coefs(r::AbstractArray)

Maps partial autocorrelations r∈(-1,1) to the coefficients φ of a stationary AR polynomial 1 - φ₁z - ... - φₚzᵖ using the Durbin-Levinson recursion. Internal function used by arfima_whittle_est.