| Title: | Nonlinear Time Series Analysis |
|---|---|
| Description: | R functions for (non)linear time series analysis with an emphasis on nonparametric autoregression and order estimation, and tests for linearity / additivity. |
| Authors: | Ottar N. Bjornstad [aut, cre] |
| Maintainer: | Ottar N. Bjornstad <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0-3 |
| Built: | 2025-01-30 03:06:11 UTC |
| Source: | https://github.com/objornstad/nlts |
add.test is a function to test the permissibility of the additive autoregressive model:
add.test(x, order, n.cond = FALSE)add.test(x, order, n.cond = FALSE)
x |
A time series (vector without missing values). |
order |
a scalar representing the order to be considered. |
n.cond |
The number of observation to condition on. The default is
|
N(t) = f1(N(t-1)) + f2(N(t-2)) + ... + fd(N(t-d)) + e(t )
against the alternative:
N(t) = F(N(t-1), N(t-2), ..., N(t-d)) + e(t)
This is the Lagrange multiplier test for additivity developed by Chen et al. (1995: test II).
a vector is returned consisting of the asymptotic chi-square value, the associated d.f. and asymptotic p.val for the test of additivity.
Chen, R., Liu, J.S. & Tsay, R.S. (1995) Additivity tests for nonlinear autoregression. Biometrika, 82, 369-383. https://doi.org/10.1093/biomet/82.2.369
Bjornstad, O.N., Begon, M., Stenseth, N.C., Falck, W., Sait, S.M., & Thompson, D.J. (1998) Population dynamics of the Indian meal moth: demographic stochasticity and delayed regulatory mechanisms. Journal of Animal Ecology, 67, 110-126. https://doi.org/10.1046/j.1365-2656.1998.00168.x
data(plodia) add.test(sqrt(plodia), order = 3)data(plodia) add.test(sqrt(plodia), order = 3)
A function to estimate the contingency periodogram to test for periodicity in categorical time series.
contingency.periodogram(x, maxper = 6, exact = FALSE)contingency.periodogram(x, maxper = 6, exact = FALSE)
x |
A vector representing the categorical time series. |
maxper |
the maximum lag (period) considered. |
exact |
If TRUE the FISHER exact test is calculated |
This is the contingency periodogram of Pierre Legendre and Pierre Dutielle to test for periodicity in categorical time series. I have coded the function so as to provide both the Fisher exact test and the asymptotic chi-square test.
An object of class "contingency.periodogram" is returned consisting of a matrix with a row for each period considered. The columns are:
exact.p |
the Fisher exact test at each lag (if exact=TRUE). |
chi2 |
the asymptotic chi-square value. |
df |
the chi-square degrees-of-freedom. |
asympt.p |
the chi-squared asymptotic p-value. |
Legendre et al. (1981) The contingency periodogram: A method for identifying rhytms in series of nonmetric ecological data. Journal of Ecology, 69, 965-979. https://doi.org/10.2307/2259648
data(plodia) data<-as.factor((scale(plodia) > 0)) fit <- contingency.periodogram(data, maxper = 9) ## Not run: plot(fit)data(plodia) data<-as.factor((scale(plodia) > 0)) fit <- contingency.periodogram(data, maxper = 9) ## Not run: plot(fit)
A function to estimate the order of a time series using cross-validation of the linear autoregressive model. Coefficients are estimated using conditional least-squares. I coded this functions to estimate the order of ecological time series. Bjornstad et al. (1998, 2001)
lin.order.cls(x, order = 1:5, n.cond = 5, echo = TRUE)lin.order.cls(x, order = 1:5, n.cond = 5, echo = TRUE)
x |
A time series without missing values |
order |
The candidate orders. The default is 1:5 |
n.cond |
The number of observation to condition on. The default is 5 (must be >= max(order)) |
echo |
if TRUE a counter for the data points and the orders is produced to monitor progress. |
The time series is normalized prior to cross-validation.
Note that if the dynamics is highly nonlinear, the nonparametric
order-estimators (ll.order) may be more appropriate. (I coded
this function to use for comparison with the nonparametric methods, because
these also uses (nonlinear) conditional least-squares.)
An object of class "lin.order" is returned consisting of the following components:
order |
the grid of orders considered. |
CVd |
the cross-validation errors across the grid of orders. |
Bjornstad, O.N., Begon, M., Stenseth, N. C., Falck, W., Sait, S. M. and Thompson, D. J. 1998. Population dynamics of the Indian meal moth: demographic stochasticity and delayed regulatory mechanisms. Journal of Animal Ecology 67:110-126. https://doi.org/10.1046/j.1365-2656.1998.00168.x Bjornstad, O.N., Sait, S.M., Stenseth, N.C., Thompson, D.J. & Begon, M. 2001. Coupling and the impact of specialised enemies on the dimensionality of prey dynamics. Nature 401: 1001-1006. https://doi.org/10.1038/35059003
data(plodia) fit <- lin.order.cls(sqrt(plodia), order=1:5) ## Not run: plot(fit) summary(fit)data(plodia) fit <- lin.order.cls(sqrt(plodia), order=1:5) ## Not run: plot(fit) summary(fit)
a function to test the permissibility of the linear autoregressive model:
lin.test(x, order)lin.test(x, order)
x |
A time series (vector without missing values). |
order |
a scalar representing the order to be considered. |
N(t) = a0 + a1N(t-1) + a2N(t-2) + ... + adN(t-d) + e(t )
against the alternative:
Nt = F(N(t-1), N(t-2), ..., N(t-d)) + e(t)
This is the Tukey one-degree-of-freedom test of linearity developed by Tsay (1986). Orders up to 5 is permissible. [although the code is easily extended].
A vector is returned consisting of the asymptotic F-value, the associated numerator and denominator d.f.'s and asymptotic p.val for the test of linearity
Tsay, R.S. (1986) Nonlinearity tests for time series. Biometrika, 73, 461-466. https://doi.org/10.1093/biomet/73.2.461
data(plodia) lin.test(sqrt(plodia), order = 3)data(plodia) lin.test(sqrt(plodia), order = 3)
A function to forecast a local polynomial ‘empirical dynamic model’.
ll.edm(x, order, bandwidth, len = NA, deg = 2)ll.edm(x, order, bandwidth, len = NA, deg = 2)
x |
A time series without missing values. |
order |
The order for the nonparametric (local polynomial) autoregression. |
bandwidth |
The bandwidth for the nonparametric (local polynomial) autoregression. |
len |
The length of the predicted time-series. If NA the length of the training time series will be used. |
deg |
The degree of the local polynomial. |
The function produces a nonlinear (nonparametric) forecast using the conditional mean method of Fan et al (1996). A Gaussian kernel is used for the local polynomial autoregression.
The bandwidth and order is best estimated with the
ll.order-function.
Missing values are NOT permitted.
If deg is set to 0, the forecast uses the Nadaraya-Watson (locally
constant) estimator of the conditional expectation against lagged-abundances.
A time series with the nonlinear (nonparametric) forecast is returned
Fan, J., Yao, Q., & Tong, H. (1996) Estimation of conditional densities and sensitivity measures in nonlinear dynamical systems. Biometrika, 83, 189-206. https://doi.org/10.1093/biomet/83.1.189
Loader, C. (1999) Local Regression and Likelihood. Springer, New York. https://doi.org/10.2307/1270956
data(plodia) sim1 <- ll.edm(sqrt(plodia), order=2, bandwidth = 1.5)data(plodia) sim1 <- ll.edm(sqrt(plodia), order=2, bandwidth = 1.5)
A function to estimate the order of a time series using the nonparametric order selection method of Cheng and Tong (1992, 1994) as modified by Yao & Tong (1994; see also Fan, Yao & Tong 1996). The method uses leave-one-out cross-validation of the locally linear regression against lagged-abundances.
ll.order(x, order = 1:5, step = 1, deg = 2, bandwidth = c(seq(0.3, 1.5, by = 0.1), 2:10), cv = TRUE, echo = TRUE)ll.order(x, order = 1:5, step = 1, deg = 2, bandwidth = c(seq(0.3, 1.5, by = 0.1), 2:10), cv = TRUE, echo = TRUE)
x |
A time series without missing values. |
order |
The candidate orders. The default is 1:5. |
step |
The time step for prediction. |
deg |
The degree of the local polynomial. |
bandwidth |
The candidate bandwidths to be considered. |
cv |
if TRUE leave-one-out cross-validation will be performed. |
echo |
if TRUE a counter shows the progress |
The time series is normalized prior to cross-validation.
A Gaussian kernel is used for the locally linear regression.
The bandwidth is optimized using cross-validation. If a single bandwidth is provided, no cross validation of bandwidth will be carried out. Highly nonlinear data will require more narrow bandwidths. If NA is returned it may be because the min bandwidth considered is too small relative to the density of data.
Missing values are NOT permitted.
If deg is set to 0, the order is estimated on the basis of the
Nadaraya-Watson (locally constant) estimator of the conditional expectation
against lagged-abundances (Cheng and Tong 1992, 1994).
An object of class "ll.order" is returned consisting of the following components:
grid |
the grid of orders, bandwidths, and CV's. |
grid$order |
the orders. |
grid$CV |
the cross-validation score
across the grid of orders and bandwidths. (If |
grid$GCV |
the generalized cross-validation score. |
grid$bandwidth |
the bandwidths. |
grid$df |
the degrees of freedom of the fitted model. |
order |
the vector of orders considered. |
deg |
The degree of the local polynomial. |
Cheng, B. & Tong, H. (1992) On consistent nonparametric order determination and chaos. Journal of Royal Statistical Society B, 54, 427-449.
Cheng, B. & Tong, H. (1994) Orthogonal projection, embedding dimension and sample size in chaotic time series from a statistical perspective. Philosophical Transactions of the Royal Society London, A. , 348, 325-341. https://doi.org/10.1098/rsta.1994.0094
Fan, J., Yao, Q., & Tong, H. (1996) Estimation of conditional densities and sensitivity measures in nonlinear dynamical systems. Biometrika, 83, 189-206. ttps://doi.org/10.1093/biomet/83.1.189
Yao, Q. & Tong, H. (1994) Quantifying the influence of initial values on non-linear prediction. Journal of Royal Statistical Society B, 56, 701-725.
Bjornstad, O.N., Sait, S.M., Stenseth, N.C., Thompson, D.J., & Begon, M. (2001) Coupling and the impact of specialised enemies on the dimensionality of prey dynamics. Nature, 409, 1001-1006. https://doi.org/10.1038/35059003
Loader, C. (1999) Local Regression and Likelihood. Springer, New York. https://doi.org/10.1007/b98858
data(plodia) fit <- ll.order(sqrt(plodia), order=1:3, bandwidth = seq(0.5, 1.5, by = 0.5)) ## Not run: plot(fit) summary(fit)data(plodia) fit <- ll.order(sqrt(plodia), order=1:3, bandwidth = seq(0.5, 1.5, by = 0.5)) ## Not run: plot(fit) summary(fit)
hack to make ll.order work with locfit >1.5. not to be called by the user.
lpx(x, nn = 0, h = 0, adpen = 0, deg = 2, acri = "none", scale = FALSE, style = "none")lpx(x, nn = 0, h = 0, adpen = 0, deg = 2, acri = "none", scale = FALSE, style = "none")
x |
... |
nn |
... |
h |
... |
adpen |
... |
deg |
... |
acri |
... |
scale |
... |
style |
... |
not to be called by the user.
Catherine Loader
A function to create matrix of lagged time series. Called by various functions.
mkx(x, lags)mkx(x, lags)
x |
A univariate time series. |
lags |
The vector of time lags. |
If lags is c(1,4), say, then the function returns a matrix that
consist of columns x(t-1), x(t-4), x(t).
A matrix of lagged abundances. The last column is the current
Upmanu Lall
Lall, U. & Sharma, A. (1996) A nearest neighbor bootstrap for time series resampling. Water Resources Research, 32, 679-693. https://doi.org/10.1029/95wr02966
This is replicate 3 in Bjornstad et al. (1998).
plodiaplodia
A vector containing 55 values
Bjornstad, O. N., M. Begon, N. C. Stenseth, W. Falck, S. M. Sait, and D. J. Thompson. 1998. Population dynamics of the Indian meal moth: demographic stochasticity and delayed regulatory mechanisms. Journal of Animal Ecology 67:110-126. https://doi.org/10.1046/j.1365-2656.1998.00168.x
‘plot’ method for "contingency.periodogram" class object.
## S3 method for class 'contingency.periodogram' plot(x, ...)## S3 method for class 'contingency.periodogram' plot(x, ...)
x |
an object of class "contingency.periodogram", usually, as a result
of a call to |
... |
generic plot arguments. |
A contingency periodogram is plotted. The line represents the critical value based on the chi-squared test (95%).
‘plot’ method for class "lin.order".
## S3 method for class 'lin.order' plot(x, ...)## S3 method for class 'lin.order' plot(x, ...)
x |
an object of class "lin.order", usually, as a result of a call to
|
... |
generic plot arguments. |
A xy-plot of order against cross-validation error is produced.
‘plot’ method for class "ll.order".
## S3 method for class 'll.order' plot(x, ...)## S3 method for class 'll.order' plot(x, ...)
x |
an object of class "ll.order", usually, as a result of a call to
|
... |
generic plot arguments. |
See ll.order for details.
A xy-plot of minimum cross-validation error against order is produced.
‘plot’ method for objects of class "lomb".
## S3 method for class 'lomb' plot(x, ...)## S3 method for class 'lomb' plot(x, ...)
x |
an object of class "lomb", usually, as a result of a call to
|
... |
generic plot arguments. |
A Lomb periodogram is composed of a xy-plot of amplitude against frequency.
‘plot’ method for class "ppll".
## S3 method for class 'ppll' plot(x, ...)## S3 method for class 'ppll' plot(x, ...)
x |
an object of class "ppll", usually, as a result of a call to
|
... |
generic plot arguments. |
See prediction.profile.ll for details.
A xy-plot of one minus the cross-validation error (i.e. the prediction accuracy against prediction time step.
‘plot’ method for class "specar.ci".
## S3 method for class 'specar.ci' plot(x, period = TRUE, ...)## S3 method for class 'specar.ci' plot(x, period = TRUE, ...)
x |
an object of class "specar.ci", usually, as a result of a call to
|
period |
if TRUE x-axis is period, if FALSE frequency. |
... |
generic plot arguments. |
A xy-plot of amplitude against period (or frequency).
portman.Q uses the cumulative ACF to test for whiteness of a time series.
portman.Q(x, K)portman.Q(x, K)
x |
A time series (vector without missing values). |
K |
the maximum lag of the ACF to be used in the test. |
This is the Ljung-Box version of the the Portmanteau test for whiteness (Tong 1990). It may in particular be useful to test for whiteness in the residuals from time series models.
A vector is returned consisting of the asymptotic chi-square value, the associated d.f. and asymptotic p.val for the test of whiteness.
Tong, H. (1990) Non-linear time series : a dynamical system approach. Clarendon Press, Oxford.
data(plodia) portman.Q(sqrt(plodia), K = 10) fit <- ar(sqrt(plodia)) portman.Q(na.omit(fit$resid), K = 10)data(plodia) portman.Q(sqrt(plodia), K = 10) fit <- ar(sqrt(plodia)) portman.Q(na.omit(fit$resid), K = 10)
Calculates the leave-one-out predicted values for the optimal ll.order object
## S3 method for class 'll.order' predict(object, ...)## S3 method for class 'll.order' predict(object, ...)
object |
an object of class "ll.order", usually, as a result of a call
to |
... |
no other arguments currently allowed |
See ll.order for details.
A data frame with observed and predicted values for the optimal ll-model is returned.
A wrapper function around ll.order to calculate prediction profiles
(a la Sugihara and May 1990 and Yao and Tong 1994). The method uses
leave-one-out cross-validation of the local regression (with CV optimized
bandwidth) against lagged-abundances at various lags.
prediction.profile.ll(x, step = 1:10, order = 1:5, deg = 2, bandwidth = c(seq(0.3, 1.5, by = 0.1), 2:10))prediction.profile.ll(x, step = 1:10, order = 1:5, deg = 2, bandwidth = c(seq(0.3, 1.5, by = 0.1), 2:10))
x |
A time series without missing values. |
step |
The vector of time steps for forward prediction. |
order |
The candidate orders. The default is 1:5. |
deg |
The degree of the local polynomial. |
bandwidth |
The candidate bandwidths to be considered. |
see ll.order for details.
An object of class "ppll" consisting of a list with the following components:
step |
the prediction steps considered. |
CV |
the cross-validation error. |
order |
the optimal order for each step. |
bandwidth |
the optimal bandwidth for each step. |
df |
the degrees of freedom for each step. |
Sugihara, G., and May, R.M. (1990) Nonlinear forecasting as a way of distinguishing chaos from measurement error in time series. Nature 344, 734-741. https://doi.org/10.1038/344734a0
Yao, Q. and Tong, H. (1994) Quantifying the influence of initial values on non-linear prediction. Journal of Royal Statistical Society B, 56, 701-725.
Fan, J., Yao, Q., and Tong, H. (1996) Estimation of conditional densities and sensitivity measures in nonlinear dynamical systems. Biometrika, 83, 189-206. https://doi.org/10.1093/biomet/83.1.189
data(plodia) fit <- prediction.profile.ll(sqrt(plodia), step=1:3, order=1:3, bandwidth = seq(0.5, 1.5, by = 0.5)) ## Not run: plot(fit)data(plodia) fit <- prediction.profile.ll(sqrt(plodia), step=1:3, order=1:3, bandwidth = seq(0.5, 1.5, by = 0.5)) ## Not run: plot(fit)
‘print’ method for class "ll.order".
## S3 method for class 'll.order' print(x, verbose = FALSE, ...)## S3 method for class 'll.order' print(x, verbose = FALSE, ...)
x |
an object of class "ll.order", usually, as a result of a call to
|
verbose |
if TRUE provides a raw-printing of the object. |
... |
no other arguments currently allowed |
See ll.order for details.
A matrix summarizing the minimum cross-validation error (cv.min) and the associated Gaussian-kernel bandwidth (bandwidth.opt) and model degrees-of-freedom for each order considered.
The function to estimate the Lomb periodogram for a spectral analysis of unevenly sampled data.
spec.lomb(y = stop("no data arg"), x = stop("no time arg"), freq = NULL)spec.lomb(y = stop("no data arg"), x = stop("no time arg"), freq = NULL)
y |
vector of length n representing the unevenly sampled time series. |
x |
the a vector (of length n) representing the times of observation. |
freq |
the frequencies at which the periodogram is to be calculated. If NULL the canonical frequencies (the Fourier frequencies) are used. |
This is the Lomb periodogram to test for periodicity in time series of unevenly sampled data.
Missing values should be deleted in both x and y before execution.
An object of class "lomb" is returned consisting of the following components:
freq |
the frequencies as supplied. |
spec |
the estimated amplitudes at the different frequencies. |
f.max |
the frequency of maximum amplitude. |
per.max |
the corresponding period of maximum amplitude. |
p |
the level of significance associated with the max period. |
Lomb, N.R. (1976) Least-squares frequency-analysis of unequally spaced data. Astrophysics and Space Science 39, 447-462.
data(plodia) y <- sqrt(plodia) x <- 1:length(y) #make some missing values y[10:19] <- NA; x[10:19] <- NA #omit NAs y <- na.omit(y); x <- na.omit(x) #the lomb p'gram fit <- spec.lomb(y, x) summary(fit) ## Not run: plot(fit)data(plodia) y <- sqrt(plodia) x <- 1:length(y) #make some missing values y[10:19] <- NA; x[10:19] <- NA #omit NAs y <- na.omit(y); x <- na.omit(x) #the lomb p'gram fit <- spec.lomb(y, x) summary(fit) ## Not run: plot(fit)
A function to estimate a "confidence interval" for the power spectrum and in particular a confidence interval for the dominant period. The function uses resampling of the autoregressive parameters to attain the estimate.
specar.ci(x, order, resamp = 500, nfreq = 100, echo = TRUE)specar.ci(x, order, resamp = 500, nfreq = 100, echo = TRUE)
x |
A time series without missing values. |
order |
a scalar representing the order to be considered. If
|
resamp |
the number of resamples of the ar-coefficients from the covariance matrix. |
nfreq |
the number of points at which to save the value for the power spectrum (and confidence envelope). |
echo |
If |
A "confidence interval" for the periodogram is obtained by resampling the ar-coefficients using the variance-covariance matrix from the ar.mle object.
If a zero'th order process is chosen by using the AIC criterion, a first order process will be used.
If the dynamics is highly nonlinear, the parametric estimate of the power spectrum may be inappropriate.
An object of class "specar.ci" is returned consisting of the following components:
order |
the ar-order. |
spectrum$freq |
the spectral frequencies. |
spectrum$spec |
the estimated power-spectrum of the data. |
resamp$spectrum |
gives the quantile summary for the resampling distribution of the spectral powers. |
resamp$maxfreq |
the full vector of output for the resampled max.frequencies. |
plot.specar.ci summary.specar.ci
data(plodia) fit <- specar.ci(sqrt(plodia), order=3, resamp=10) ## Not run: plot(fit, period=FALSE) summary(fit)data(plodia) fit <- specar.ci(sqrt(plodia), order=3, resamp=10) ## Not run: plot(fit, period=FALSE) summary(fit)
‘summary’ method for class "lin.order".
## S3 method for class 'lin.order' summary(object, ...)## S3 method for class 'lin.order' summary(object, ...)
object |
an object of class "lin.order", usually, as a result of a call
to |
... |
no other arguments currently allowed |
A slightly prettified version of the object is printed.
‘summary’ method for class "ll.order".
## S3 method for class 'll.order' summary(object, GCV = FALSE, ...)## S3 method for class 'll.order' summary(object, GCV = FALSE, ...)
object |
an object of class "ll.order", usually, as a result of a call
to |
GCV |
if TRUE (or if cross-validation was not done), uses GCV values. |
... |
no other arguments currently allowed |
See ll.order for details.
A matrix summarizing the minimum cross-validation error (cv.min) and the associated Gaussian-kernel bandwidth (bandwidth.opt) and model degrees-of-freedom for each order considered.
‘summary’ method for objects of class "lomb".
## S3 method for class 'lomb' summary(object, ...)## S3 method for class 'lomb' summary(object, ...)
object |
an object of class "lomb", usually, as a result of a call to
|
... |
generic plot arguments. |
A list summarizing the analysis is printed:
period |
the dominant period. |
p.val |
the p.value. |
‘summary’ method for objects of class "specar.ci".
## S3 method for class 'specar.ci' summary(object, period = TRUE, ...)## S3 method for class 'specar.ci' summary(object, period = TRUE, ...)
object |
an object of class "specar.ci", usually, as a result of a call
to |
period |
If TRUE the summary is given in terms of the period, if false it is in terms of the frequency |
... |
generic plot arguments. |
A list summarizing the analysis is printed:
period |
the dominant period. |
p.val |
the p.value. |