longmemo (version 1.1-2)

WhittleEst: Whittle Estimator for Fractional Gaussian Noise / Fractional ARIMA

Description

Computes Whittle's approximate MLE for fractional Gaussian noise or fractional ARIMA (=: fARIMA) models, according to Beran's prescript.

Relies on minmizing Qeta() (\(= \tilde{Q}(\eta)\), which itself is based on the “true” spectrum of the corresponding process; for the spectrum, either specFGN or specARIMA is used.

Usage

WhittleEst(x,
           periodogr.x = per(if(scale) x/sd(x) else x)[2:((n+1) %/% 2)],
           n = length(x), scale = FALSE,
           model = c("fGn", "fARIMA"),
           p, q,
           start = list(H= 0.5, AR= numeric(), MA=numeric()),
           verbose = getOption("verbose"))
# S3 method for WhittleEst
print(x, digits = getOption("digits"), …)

Arguments

x

numeric vector representing a time series. Maybe omitted if periodogr.x and n are specified instead.

periodogr.x

the (raw) periodogram of x; the default, as by Beran, uses per, but tapering etc may be an alternative, see also spec.pgram.

n

length of the time series, length(x).

scale

logical indicating if x should be standardized to (sd) scale 1; originally, scale = TRUE used to be built-in; for compatibility with other methods, notably plotting spectra, scale = FALSE seems a more natural default.

model

numeric vector representing a time series.

p,q

optional integers specifying the AR and MA orders of the fARIMA model, i.e., only applicable when model is "fARIMA".

start

list of starting values; currently necessary for model = "fARIMA" and with a reasonable default for model = "fGn".

verbose

logical indicating if iteration output should be printed.

digits,…

optional arguments for print method, see print.default.

Value

An object of class WhittleEst which is basically a list with components

call

the function call.

model

= input

n

time series length length(x).

p,q

for "fARIMA": order of AR and MA parts, respectively.

coefficients

numeric 4-column matrix of coefficients with estimate of the full parameter vector \(\eta\), its standard error estimates, z- and P-values. This includes the Hurst parameter \(H\).

theta1

the scale parameter \(\hat{\theta_1}\), see Qeta.

vcov

the variance-covariance matrix for \(\eta\).

periodogr.x

= input (with default).

spec

the spectral estimate \(\hat{f}(\omega_j)\).

There is a print method, and coef, confint or vcov methods work as well for objects of class "WhittleEst".

References

Beran, Jan (1994). Statistics for Long-Memory Processes; Chapman & Hall. (Section 6.1, p.116--119; 12.1.3, p.223 ff)

See Also

Qeta is the function minimized by these Whittle estimators.

FEXPest for an alternative model with Hurst parameter, also estimated by a “Whittle” approximate MLE, i.e., a Whittle's estimator in the more general sense.

The plot method, plot.WhittleEst.

Examples

# NOT RUN {
data(NileMin)
(f.Gn.N  <- WhittleEst(NileMin))                             # H = 0.837
(f.A00.N <- WhittleEst(NileMin, model = "fARIMA", p=0, q=0)) # H = 0.899
confint(f.Gn.N)
confint(f.A00.N)

data(videoVBR)
(f.GN    <- WhittleEst(videoVBR))

## similar {but faster !}:
(f.am00  <- WhittleEst(videoVBR, model = "fARIMA", p=0, q=0))
rbind(f.am00$coef,
      f.GN  $coef)# really similar

f.am11  <- WhittleEst(videoVBR, model = "fARIMA",
                      start= list(H= .5, AR = .5, MA= .5))
f.am11
vcov(f.am11)

op <- if(require("sfsmisc"))
  mult.fig(3, main = "Whittle Estimators for videoVBR data")$old.par  else
  par(mar = c(3,1), mgp = c(1.5, 0.6, 0), mar = c(4,4,2,1)+.1)
plot(f.GN)
plot(f.am00)
plot(f.am11)

et <- as.list(coef(f.am11))
et$AR <- c(et$AR, 0, 0) # two more AR coefficients ..
f.am31  <- WhittleEst(videoVBR, model = "fARIMA", start = et)
## ... warning non nonconvergence,  but "kind of okay":
lines(f.am31, col = "red3") ## drawing on top of ARMA(1,1) above - *small* diff

f.am31 # not all three are "significant"
round(cov2cor(vcov(f.am31)), 3) # and they are highly correlated

et <- as.list(coef(f.am31))
et$AR <- unname(unlist(et[c("AR1", "AR2")]))
f.am21  <- WhittleEst(videoVBR, model = "fARIMA",
                      start = c(et[c("H","AR", "MA")]))
f.am21
lines(f.am21, col = adjustcolor("gold", .4), lwd=4)

par(op)## (reset graphic layout)

##--> ?plot.WhittleEst  for an example using  'periodogr.x'
# }
# NOT RUN {
<!-- %% ^^ "FIXME": do it here.. -->
# }