Time Series Analysis and Control Program Package

Description

R functions for statistical analysis and control of time series.

Details

This package provides functions for statistical analysis, prediction and control of time series. The original TIMSAC (TIMe Series Analysis and Control) or TIMSAC-72 was published in Akaike and Nakagawa (1972). After that, TIMSAC-74, TIMSAC-78 and TIMSAC-84 were published as the TIMSAC series in Computer Science Monograph.

For overview of models and information criteria for model selection, see ../doc/timsac-guide_e.pdf or ../doc/timsac-guide_j.pdf (in Japanese).

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.6, Timsac74, A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

H.Akaike, T.Ozaki, M.Ishiguro, Y.Ogata, G.Kitagawa, Y-H.Tamura, E.Arahata, K.Katsura and Y.Tamura (1985) Computer Science Monograph, No.22, Timsac84 Part 1. The Institute of Statistical Mathematics.

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Airpollution Data

Description

An airpollution data for testing perars.

Usage

  data(Airpollution)

Format

A time series of 372 observations.

Source

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Amerikamaru Data

Description

A multivariate non-stationary data for testing blomar.

Usage

  data(Amerikamaru)

Format

A 2-dimensional array with 896 observations on 2 variables.

Source

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Blsallfood Data

Description

The BLSALLFOOD data. (the Bureau of Labor Statistics, all employees in food industries, January 1967 - December 1979)

Usage

  data(Blsallfood)

Format

A time series of 156 observations.

Source

H.Akaike, T.Ozaki, M.Ishiguro, Y.Ogata, G.Kitagawa, Y-H.Tamura, E.Arahata, K.Katsura and Y.Tamura (1984) Computer Science Monographs, Timsac-84 Part 1. The Institute of Statistical Mathematics.

Time series of Canadian lynx data

Description

A time series of Canadian lynx data for testing unimar, unibar, bsubst and exsar.

Usage

  data(Canadianlynx)

Format

A time series of 114 observations.

Source

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Labor force Data

Description

Labor force U.S. unemployed 16 years or over (1972-1978) data.

Usage

  data(LaborData)

Format

A time series of 72 observations.

Source

H.Akaike, T.Ozaki, M.Ishiguro, Y.Ogata, G.Kitagawa, Y-H.Tamura, E.Arahata, K.Katsura and Y.Tamura (1985) Computer Science Monograph, No.22, Timsac84 Part 1. The Institute of Statistical Mathematics.

Power Plant Data

Description

A Power plant data for testing mulbar and mulmar.

Usage

  data(Powerplant)

Format

A 2-dimensional array with 500 observations on 3 variables.

Source

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

ARMA Model Fitting

Description

Fit an ARMA model with specified order by using DAVIDON's algorithm.

Usage

  armafit(y, model.order)

Arguments

Details

The maximum likelihood estimates of the coefficients of a scalar ARMA model

y(t) - a(1)y(t-1) -...- a(p)y(t-p) = u(t) - b(1)u(t-1) -...- b(q)u(t-q)

of a time series y(t) are obtained by using DAVIDON's algorithm. Pure autoregression is not allowed.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

# "arima.sim" is a function in "stats".
# Note that the sign of MA coefficient is opposite from that in "timsac".
y <- arima.sim(list(order=c(2,0,1), ar=c(0.64,-0.8), ma=-0.5), n = 1000)
z <- armafit(y, model.order = c(2,1))
z$arcoef
z$macoef

Power Spectrum

Description

Compute power spectrum estimates for two trigonometric windows of Blackman-Tukey type by Goertzel method.

Usage

auspec(y, lag = NULL, window = "Akaike", log = FALSE, plot = TRUE)

Arguments

Details

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

y <- arima.sim(list(order=c(2,0,0), ar=c(0.64,-0.8)), n = 200)
auspec(y, log = TRUE)

Autocorrelation

Description

Estimate autocovariances and autocorrelations.

Usage

autcor(y, lag = NULL, plot = TRUE, lag_axis = TRUE)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1 for the normal distribution 
y <- rnorm(200)
autcor(y, lag_axis = FALSE)

# Example 2 for the ARIMA model
y <- arima.sim(list(order=c(2,0,0), ar=c(0.64,-0.8)), n = 200)
autcor(y, lag = 20)

Automatic ARMA Model Fitting

Description

Provide an automatic ARMA model fitting procedure. Models with various orders are fitted and the best choice is determined with the aid of the statistics AIC.

Usage

  autoarmafit(y, max.order = NULL)

Arguments

Details

The maximum likelihood estimates of the coefficients of a scalar ARMA model

y(t) - a(1)y(t-1) -...- a(p)y(t-p) = u(t) - b(1)u(t-1) -...- b(q)u(t-q)

of a time series y(t) are obtained by using DAVIDON's variance algorithm. Where p is AR order, q is MA order and u(t) is a zero mean white noise. Pure autoregression is not allowed.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

# "arima.sim" is a function in "stats".
# Note that the sign of MA coefficient is opposite from that in "timsac".
y <- arima.sim(list(order=c(2,0,1),ar=c(0.64,-0.8),ma=-0.5), n = 1000)
autoarmafit(y)

Bayesian Seasonal Adjustment Procedure

Description

Decompose a nonstationary time series into several possible components.

Usage

  baysea(y, period = 12, span = 4, shift = 1, forecast = 0, trend.order = 2,
         seasonal.order = 1, year = 0, month = 1, out = 0, rigid = 1,
         zersum = 1, delta = 7, alpha = 0.01, beta = 0.01, gamma = 0.1,
         spec = TRUE, plot = TRUE, separate.graphics = FALSE)

Arguments

Details

This function realized a decomposition of time series y into the form

y(t) = T(t) + S(t) + I(t) + TDC(t) + OCF(t)

where T(t) is trend component, S(t) is seasonal component, I(t) is irregular, TDC(t) is trading day factor and OCF(t) is outlier correction factor. For the purpose of comparison of models the criterion ABIC is defined

ABIC = -2 \log(maximum\ likelihood\ of\ the\ model).

Smaller value of ABIC represents better fit.

Value

References

H.Akaike, T.Ozaki, M.Ishiguro, Y.Ogata, G.Kitagawa, Y-H.Tamura, E.Arahata, K.Katsura and Y.Tamura (1985) Computer Science Monograph, No.22, Timsac84 Part 1. The Institute of Statistical Mathematics.

Examples

data(LaborData)
baysea(LaborData, forecast = 12)

Bispectrum

Description

Compute bi-spectrum using the direct Fourier transform of sample third order moments.

Usage

bispec(y, lag = NULL, window = "Akaike", log = FALSE, plot = TRUE)

Arguments

Details

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.6, Timsac74, A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Examples

data(bispecData)
bispec(bispecData, lag = 30)

Univariate Test Data

Description

A univariate data for testing bispec and thirmo.

Usage

  data(bispecData)

Format

A time series of 1500 observations.

Source

H.Akaike, E.Arahata and T.Ozaki (1976) Computer Science Monograph, No.6, Timsac74 A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Bayesian Method of Locally Stationary AR Model Fitting; Scalar Case

Description

Locally fit autoregressive models to non-stationary time series by a Bayesian procedure.

Usage

blocar(y, max.order = NULL, span, plot = TRUE)

Arguments

Details

The basic AR model of scalar time series y(t) (t=1, \ldots ,n) is given by

y(t) = a(1)y(t-1) + a(2)y(t-2) + \ldots + a(p)y(t-p) + u(t),

where p is order of the model and u(t) is Gaussian white noise with mean 0 and variance v. At each stage of modeling of locally AR model, a two-step Bayesian procedure is applied

AIC of the model fitted to the new span is defined by

AIC = ns \log( sd ) + 2k,

where ns is the length of new data, sd is innovation variance and k is the equivalent number of parameters, defined as the sum of squares of the Bayesian weights. AIC of the model fitted to the preceding spans are defined by

AIC( j+1 ) = ns \log( sd(j) ) + 2,

where sd(j) is the prediction error variance by the model fitted to j periods former span.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-Stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike (1978) A Bayesian Extension of the Minimum AIC Procedure of Autoregressive Model Fitting. Research Memo. NO.126. The Institute of The Statistical Mathematics.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(locarData)
z <- blocar(locarData, max.order = 10, span = 300)
z$arcoef

Bayesian Method of Locally Stationary Multivariate AR Model Fitting

Description

Locally fit multivariate autoregressive models to non-stationary time series by a Bayesian procedure.

Usage

  blomar(y, max.order = NULL, span)

Arguments

Details

The basic AR model is given by

y(t) = A(1)y(t-1) + A(2)y(t-2) + \ldots + A(p)y(t-p) + u(t),

where p is order of the AR model and u(t) is innovation variance v.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure for the Modeling of Non-stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike (1978) A Bayesian Extension of The Minimum AIC Procedure of Autoregressive Model Fitting. Research Memo. NO.126. The institute of Statistical Mathematics.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Amerikamaru)
blomar(Amerikamaru, max.order = 10, span = 300)

Bayesian Type All Subset Analysis

Description

Produce Bayesian estimates of time series models such as pure AR models, AR models with non-linear terms, AR models with polynomial type mean value functions, etc. The goodness of fit of a model is checked by the analysis of several steps ahead prediction errors.

Usage

bsubst(y, mtype, lag = NULL, nreg, reg = NULL, term.lag = NULL, cstep = 5,
       plot = TRUE)

Arguments

Details

The AR model is given by ( mtype = 2 )

y(t) = a(1)y(t-1) + ... + a(p)y(t-p) + u(t).

The non-linear model is given by ( mtype = 2, 3 )

y(t) = a(1)z(t,1) + a(2)z(t,2) + ... + a(p)z(t,p) + u(t).

Where p is AR order and u(t) is Gaussian white noise with mean 0 and variance v(p).

Value

References

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Canadianlynx)
Regressor <- matrix(
     c( 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 1, 2, 1, 3, 1, 2, 3,
        0, 0, 0, 0, 0, 0, 0, 0, 0,  0,  0,  0, 1, 2, 2, 3, 1, 2, 3,
        0, 0, 0, 0, 0, 0, 0, 0, 0,  0,  0,  0, 0, 0, 0, 0, 1, 2, 3 ),
     nrow = 3, ncol = 19, byrow = TRUE)
z <- bsubst(Canadianlynx, mtype = 2, lag = 12, nreg = 19, Regressor)
z$arcoef.bay

Canonical Correlation Analysis of Scalar Time Series

Description

Fit an ARMA model to stationary scalar time series through the analysis of canonical correlations between the future and past sets of observations.

Usage

canarm(y, lag = NULL, max.order = NULL, plot = TRUE)

Arguments

Details

The ARMA model of stationary scalar time series y(t) (t=1,...,n) is given by

y(t) - a(1)y(t-1) - ...- a(p)y(t-p) = u(t) - b(1)u(t-1) - ... - b(q)u(t-q),

where p is AR order and q is MA order.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

# "arima.sim" is a function in "stats".
# Note that the sign of MA coefficient is opposite from that in "timsac".
y <- arima.sim(list(order=c(2,0,1), ar=c(0.64,-0.8), ma=c(-0.5)), n = 1000)
z <- canarm(y, max.order = 30)
z$arcoef
z$macoef

Canonical Correlation Analysis of Vector Time Series

Description

Analyze canonical correlation of a d-dimensional multivariate time series.

Usage

canoca(y)

Arguments

Details

First AR model is fitted by the minimum AIC procedure. The results are used to ortho-normalize the present and past variables. The present and future variables are tested successively to decide on the dependence of their predictors. When the last DIC (=chi-square - 2.0*N.D.F.) is negative the predictor of the variable is decided to be linearly dependent on the antecedents.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow= TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(1000*3), nrow = 1000, ncol = 3)
y <- mfilter(x, ar, "recursive")
z <- canoca(y)
z$arcoef

Covariance Generation

Description

Produce the Fourier transform of a power gain function in the form of an autocovariance sequence.

Usage

  covgen(lag, f, gain, plot = TRUE)

Arguments

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

spec <- raspec(h = 100, var = 1, arcoef = c(0.64,-0.8), plot = FALSE)
covgen(lag = 100, f = 0:100/200, gain = spec)

Time Series Decomposition (Seasonal Adjustment) by Square-Root Filter

Description

Decompose a nonstationary time series into several possible components by square-root filter.

Usage

  decomp(y, trend.order = 2, ar.order = 2, seasonal.order = 1, 
         period = 1, log = FALSE, trade = FALSE, diff = 1,
         miss = 0, omax = 99999.9, plot = TRUE, ...)

Arguments

Details

The Basic Model

y(t) = T(t) + AR(t) + S(t) + TD(t) + W(t)

where T(t) is trend component, AR(t) is AR process, S(t) is seasonal component, TD(t) is trading day factor and W(t) is observational noise.

Component Models

• Trend component (trend.order m1)

  m1 = 1 : T(t) = T(t-1) + v1(t)

  m1 = 2 : T(t) = 2T(t-1) - T(t-2) + v1(t)

  m1 = 3 : T(t) = 3T(t-1) - 3T(t-2) + T(t-2) + v1(t)

• AR component (ar.order m2)

  AR(t) = a(1)AR(t-1) + \ldots + a(m2)AR(t-m2) + v2(t)
  

• Seasonal component (seasonal.order k, frequency f)

  k=1 : S(t) = -S(t-1) - \ldots - S(t-f+1) + v3(t)
k=2 : S(t) = -2S(t-1) - \ldots -f\ S(t-f+1) - \ldots - S(t-2f+2) + v3(t)
  

• Trading day effect

  TD(t) = b(1) TRADE(t,1) + \ldots + b(7) TRADE(t,7)

  where TRADE(t,i) is the number of i-th days of the week in t-th data and b(1)\ +\ \ldots\ +\ b(7)\ =\ 0.

Value

An object of class "decomp", which is a list with the following components:

References

G.Kitagawa (1981) A Nonstationary Time Series Model and Its Fitting by a Recursive Filter Journal of Time Series Analysis, Vol.2, 103-116.

W.Gersch and G.Kitagawa (1983) The prediction of time series with Trends and Seasonalities Journal of Business and Economic Statistics, Vol.1, 253-264.

G.Kitagawa (1984) A smoothness priors-state space modeling of Time Series with Trend and Seasonality Journal of American Statistical Association, VOL.79, NO.386, 378-389.

Examples

data(Blsallfood)
y <- ts(Blsallfood, start=c(1967,1), frequency=12)
z <- decomp(y, trade = TRUE)
z$aic
z$lkhd
z$sigma2
z$tau1
z$tau2
z$tau3

Exact Maximum Likelihood Method of Scalar AR Model Fitting

Description

Produce exact maximum likelihood estimates of the parameters of a scalar AR model.

Usage

  exsar(y, max.order = NULL, plot = FALSE)

Arguments

Details

The AR model is given by

y(t) = a(1)y(t-1) + .... + a(p)y(t-p) + u(t)

where p is AR order and u(t) is a zero mean white noise.

Value

References

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Canadianlynx)
z <- exsar(Canadianlynx, max.order = 14)
z$arcoef.maice
z$arcoef.mle

Auto And/Or Cross Correlations via FFT

Description

Compute auto and/or cross covariances and correlations via FFT.

Usage

  fftcor(y, lag = NULL, isw = 4, plot = TRUE, lag_axis = TRUE)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1
x <- rnorm(200)
y <- rnorm(200)
xy <- array(c(x,y), dim = c(200,2))
fftcor(xy, lag_axis = FALSE)

# Example 2
xorg <- rnorm(1003)
x <- matrix(0, nrow = 1000, ncol = 2)
x[, 1] <- xorg[1:1000]
x[, 2] <- xorg[4:1003] + 0.5*rnorm(1000)
fftcor(x, lag = 20)

FPE Auto

Description

Perform FPE(Final Prediction Error) computation for one-dimensional AR model.

Usage

  fpeaut(y, max.order = NULL)

Arguments

Details

The AR model is given by

y(t) = a(1)y(t-1) + .... + a(p)y(t-p) + u(t)

where p is AR order and u(t) is a zero mean white noise.

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

y <- arima.sim(list(order=c(2,0,0), ar=c(0.64,-0.8)), n = 200)
fpeaut(y, max.order = 20)

AR model Fitting for Control

Description

Perform AR model fitting for control.

Usage

  fpec(y, max.order = NULL, control = NULL, manip = NULL)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
fpec(y, max.order = 10)

Non-stationary Test Data

Description

A non-stationary data for testing mlocar and blocar.

Usage

  data(locarData)

Format

A time series of 1000 observations.

Source

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Maximum Likelihood Computation of Markovian Model

Description

Compute maximum likelihood estimates of Markovian model.

Usage

markov(y)

Arguments

Details

This function is usually used with simcon.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.5, Timsac74, A Time Series Analysis and Control Program Package (1). The Institute of Statistical Mathematics.

Examples

x <- matrix(rnorm(1000*2), nrow = 1000, ncol = 2)
ma <- array(0, dim = c(2,2,2))
ma[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ma[, , 2] <- matrix(c( -0.2,  0.0,
                       -0.1, -0.3), nrow = 2, ncol = 2, byrow = TRUE)
y <- mfilter(x, ma, "convolution")
ar <- array(0, dim = c(2,2,3))
ar[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 2] <- matrix(c( -0.5, -0.2,
                       -0.2, -0.5), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 3] <- matrix(c( -0.3, -0.05,
                       -0.1, -0.30), nrow = 2, ncol = 2, byrow = TRUE)
z <- mfilter(y, ar, "recursive")
markov(z)

Linear Filtering on a Multivariate Time Series

Description

Applies linear filtering to a multivariate time series.

Usage

mfilter(x, filter, method = c("convolution","recursive"), init)

Arguments

Details

This is a multivariate version of "filter" function. Missing values are allowed in 'x' but not in 'filter' (where they would lead to missing values everywhere in the output). Note that there is an implied coefficient 1 at lag 0 in the recursive filter, which gives

y[i,]' =x[,i]' + f[,,1] \times y[i-1,]' + ... +f[,,p] \times y[i-p,]',

No check is made to see if recursive filter is invertible: the output may diverge if it is not. The convolution filter is

y[i,]' = f[,,1] \times x[i,]' + ... + f[,,p] \times x[i-p+1,]'.

Value

mfilter returns a time series object.

Note

'convolve(, type="filter")' uses the FFT for computations and so may be faster for long filters on univariate time series (and so the time alignment is unclear), nor does it handle missing values. 'filter' is faster for a filter of length 100 on a series 1000, for examples.

See Also

convolve, arima.sim

Examples

#AR model simulation
ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(100*3), nrow = 100, ncol = 3)
y <- mfilter(x, ar, "recursive")

#Back to white noise
ma <- array(0, dim = c(3,3,3))
ma[, , 1] <- diag(3)
ma[, , 2] <- -ar[, , 1]
ma[, , 3] <- -ar[, , 2]
z <- mfilter(y, ma, "convolution")
mulcor(z)

#AR-MA model simulation
x <- matrix(rnorm(1000*2), nrow = 1000, ncol = 2)
ma <- array(0, dim = c(2,2,2))
ma[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ma[, , 2] <- matrix(c( -0.2,  0.0,
                       -0.1, -0.3), nrow = 2, ncol = 2, byrow = TRUE)
y <- mfilter(x, ma, "convolution")

ar <- array(0, dim = c(2,2,3))
ar[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 2] <- matrix(c( -0.5, -0.2,
                       -0.2, -0.5), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 3] <- matrix(c( -0.3, -0.05,
                       -0.1, -0.30), nrow = 2, ncol = 2, byrow = TRUE)
z <- mfilter(y, ar, "recursive")

Minimum AIC Method of Locally Stationary AR Model Fitting; Scalar Case

Description

Locally fit autoregressive models to non-stationary time series by minimum AIC procedure.

Usage

  mlocar(y, max.order = NULL, span, const = 0, plot = TRUE)

Arguments

Details

The data of length n are divided into k locally stationary spans,

|<-- n_1 -->|<-- n_2 -->|<-- n_3 -->| ..... |<-- n_k -->|

where n_i (i=1,\ldots,k) denotes the number of basic spans, each of length span, which constitute the i-th locally stationary span. At each local span, the process is represented by a stationary autoregressive model.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-Stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(locarData)
z <- mlocar(locarData, max.order = 10, span = 300, const = 0)
z$arcoef

Minimum AIC Method of Locally Stationary Multivariate AR Model Fitting

Description

Locally fit multivariate autoregressive models to non-stationary time series by the minimum AIC procedure using the householder transformation.

Usage

  mlomar(y, max.order = NULL, span, const = 0)

Arguments

Details

The data of length n are divided into k locally stationary spans,

|<-- n_1 -->|<-- n_2 -->|<-- n_3 -->| ..... |<-- n_k -->|

where n_i (i=1,\ldots,k) denoted the number of basic spans, each of length span, which constitute the i-th locally stationary span. At each local span, the process is represented by a stationary autoregressive model.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-Stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Amerikamaru)
mlomar(Amerikamaru, max.order = 10, span = 300, const = 0)

Multivariate Bayesian Method of AR Model Fitting

Description

Determine multivariate autoregressive models by a Bayesian procedure. The basic least squares estimates of the parameters are obtained by the householder transformation.

Usage

  mulbar(y, max.order = NULL, plot = FALSE)

Arguments

Details

The statistic AIC is defined by

AIC = n \log(det(v)) + 2k,

where n is the number of data, v is the estimate of innovation variance matrix, det is the determinant and k is the number of free parameters.

Bayesian weight of the m-th order model is defined by

W(n) = const \times \frac{C(m)}{m+1},

where const is the normalizing constant and C(m)=\exp(-0.5 AIC(m)). The Bayesian estimates of partial autoregression coefficient matrices of forward and backward models are obtained by (m = 1,\ldots,lag)

G(m) = G(m) D(m),

H(m) = H(m) D(m),

where the original G(m) and H(m) are the (conditional) maximum likelihood estimates of the highest order coefficient matrices of forward and backward AR models of order m and D(m) is defined by

D(m) = W(m) + \ldots + W(lag).

The equivalent number of parameters for the Bayesian model is defined by

ek = \{ D(1)^2 + \ldots + D(lag)^2 \} id + \frac{id(id+1)}{2}

where id denotes dimension of the process.

Value

References

H.Akaike (1978) A Bayesian Extension of The Minimum AIC Procedure of Autoregressive Model Fitting. Research Memo. NO.126, The Institute of Statistical Mathematics.

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Powerplant)
z <- mulbar(Powerplant, max.order = 10)
z$pacoef.for
z$pacoef.back

Multiple Correlation

Description

Estimate multiple correlation.

Usage

  mulcor(y, lag = NULL, plot = TRUE, lag_axis = TRUE)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1 
y <- rnorm(1000)
dim(y) <- c(500,2)
mulcor(y, lag_axis = FALSE)

# Example 2
xorg <- rnorm(1003)
x <- matrix(0, nrow = 1000, ncol = 2)
x[, 1] <- xorg[1:1000]
x[, 2] <- xorg[4:1003] + 0.5*rnorm(1000)
mulcor(x, lag = 20)

Frequency Response Function (Multiple Channel)

Description

Compute multiple frequency response function, gain, phase, multiple coherency, partial coherency and relative error statistics.

Usage

 mulfrf(y, lag = NULL, iovar = NULL)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
mulfrf(y, lag = 20)

Multivariate Case of Minimum AIC Method of AR Model Fitting

Description

Fit a multivariate autoregressive model by the minimum AIC procedure. Only the possibilities of zero coefficients at the beginning and end of the model are considered. The least squares estimates of the parameters are obtained by the householder transformation.

Usage

mulmar(y, max.order = NULL, plot = FALSE)

Arguments

Details

Multivariate autoregressive model is defined by

y(t) = A(1)y(t-1) + A(2)y(t-2) +\ldots+ A(p)y(t-p) + u(t),

where p is order of the model and u(t) is Gaussian white noise with mean 0 and variance matrix matv. AIC is defined by

AIC = n \log(det(v)) + 2k,

where n is the number of data, v is the estimate of innovation variance matrix, det is the determinant and k is the number of free parameters.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

# Example 1
data(Powerplant)
z <- mulmar(Powerplant, max.order = 10)
z$arcoef

# Example 2
ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3,byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
z <- mulmar(y, max.order = 10)
z$arcoef

Relative Power Contribution

Description

Compute relative power contributions in differential and integrated form, assuming the orthogonality between noise sources.

Usage

  mulnos(y, max.order = NULL, control = NULL, manip = NULL, h)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
mulnos(y, max.order = 10, h = 20)

Multiple Rational Spectrum

Description

Compute rational spectrum for d-dimensional ARMA process.

Usage

mulrsp(h, d, cov, ar = NULL, ma = NULL, log = FALSE, plot = TRUE, ...)

Arguments

Details

ARMA process :

y(t) - A(1)y(t-1) -...- A(p)y(t-p) = u(t) - B(1)u(t-1) -...- B(q)u(t-q)

where u(t) is a white noise with zero mean vector and covariance matrix cov.

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1 for the normal distribution
xorg <- rnorm(1003)
x <- matrix(0, nrow = 1000, ncol = 2)
x[, 1] <- xorg[1:1000]
x[, 2] <- xorg[4:1003] + 0.5*rnorm(1000)
aaa <- ar(x)
mulrsp(h = 20, d = 2, cov = aaa$var.pred, ar = aaa$ar)

# Example 2 for the AR model
ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
z <- fpec(y, max.order = 10)
mulrsp(h = 20, d = 3, cov = z$perr, ar = z$arcoef)

Multiple Spectrum

Description

Compute multiple spectrum estimates using Akaike window or Hanning window.

Usage

  mulspe(y, lag = NULL, window = "Akaike", plot = TRUE, ...)

Arguments

Details

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

sgnl <- rnorm(1003)
x <- matrix(0, nrow = 1000, ncol = 2)
x[, 1] <- sgnl[4:1003]
# x[i,2] = 0.9*x[i-3,1] + 0.2*N(0,1)
x[, 2] <- 0.9*sgnl[1:1000] + 0.2*rnorm(1000)
mulspe(x, lag = 100, window = "Hanning")

Non-stationary Power Spectrum Analysis

Description

Locally fit autoregressive models to non-stationary time series by AIC criterion.

Usage

nonst(y, span, max.order = NULL, plot = TRUE)

Arguments

Details

The basic AR model is given by

y(t) = A(1)y(t-1) + A(2)y(t-2) +...+ A(p)y(t-p) + u(t),

where p is order of the AR model and u(t) is innovation variance. AIC is defined by

AIC = n \log(det(sd)) + 2k,

where n is the length of data, sd is the estimates of the innovation variance and k is the number of parameter.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1976) Computer Science Monograph, No.6, Timsac74 A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Examples

# Non-stationary Test Data
data(nonstData)
nonst(nonstData, span = 700, max.order = 49)

Non-stationary Test Data

Description

A non-stationary data for testing nonst.

Usage

  data(nonstData)

Format

A time series of 2100 observations.

Source

H.Akaike, E.Arahata and T.Ozaki (1976) Computer Science Monograph, No.6, Timsac74 A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Optimal Controller Design

Description

Compute optimal controller gain matrix for a quadratic criterion defined by two positive definite matrices Q and R.

Usage

  optdes(y, max.order = NULL, ns, q, r)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Multivariate Example Data
ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow= 3, ncol= 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow= 3, ncol= 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
q.mat <- matrix(c(0.16,0,0,0.09), nrow = 2, ncol = 2)
r.mat <- as.matrix(0.001)
optdes(y, ns = 20, q = q.mat, r = r.mat)

Optimal Control Simulation

Description

Perform optimal control simulation and evaluate the means and variances of the controlled and manipulated variables X and Y.

Usage

  optsim(y, max.order = NULL, ns, q, r, noise = NULL, len, plot = TRUE)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Multivariate Example Data
ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,    0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1, 0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
q.mat <- matrix(c(0.16,0,0,0.09), nrow = 2, ncol = 2)
r.mat <- as.matrix(0.001)
optsim(y, max.order = 10, ns = 20, q = q.mat, r = r.mat, len = 20)

Periodic Autoregression for a Scalar Time Series

Description

This is the program for the fitting of periodic autoregressive models by the method of least squares realized through householder transformation.

Usage

  perars(y, ni, lag = NULL, ksw = 0)

Arguments

Details

Periodic autoregressive model (i=1, \ldots, nd, j=1, \ldots, ni) is defined by

z(i,j) = y(ni(i-1)+j),

z(i,j) = c(j) + A(1,j,0)z(i,1) + \ldots + A(j-1,j,0)z(i,j-1) + A(1,j,1)z(i-1,1) + \ldots + A(ni,j,1)z(i-1,ni) + \ldots + u(i,j),

where nd is the number of periods, ni is the number of instants in one period and u(i,j) is the Gaussian white noise. When ksw is set to '0', the constant term c(j) is excluded.

The statistics AIC is defined by AIC = n \log(det(v)) + 2k, where n is the length of data, v is the estimate of the innovation variance matrix and k is the number of parameters. The outputs are the estimates of the regression coefficients and innovation variance of the periodic AR model for each instant.

Value

References

M.Pagano (1978) On Periodic and Multiple Autoregressions. Ann. Statist., 6, 1310–1317.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Airpollution)
perars(Airpollution, ni = 6, lag = 2, ksw = 1)

Plot Trend, Seasonal, AR Components and Trading Day Factor

Description

Plot trend component, seasonal component, AR component, noise and trading day factor returned by decomp.

Usage

## S3 method for class 'decomp'
plot(x, ...)

Arguments

Plot Spectrum

Description

Plot spectrum returned by mulspe and mulrsp. On and lower diagonal are real parts, and upper diagonal are imaginary parts.

Usage

## S3 method for class 'specmx'
plot(x, plot.scale = TRUE, ...)

Arguments

Prediction Program

Description

Operate on a real record of a vector process and compute predicted values.

Usage

prdctr(y, r, s, h, arcoef, macoef = NULL, impulse = NULL, v, plot = TRUE)

Arguments

Details

One step ahead Prediction starts at time R and ends at time S. Prediction is continued without new observations until time S+H. Basic model is the autoregressive moving average model of y(t) which is given by

y(t) - A(t)y(t-1) -...- A(p)y(t-p) = u(t) - B(1)u(t-1) -...- B(q)u(t-q),

where p is AR order and q is MA order.

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.6, Timsac74, A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Examples

# "arima.sim" is a function in "stats".
# Note that the sign of MA coefficient is opposite from that in "timsac".
y <- arima.sim(list(order=c(2,0,1), ar=c(0.64,-0.8), ma=c(-0.5)), n = 1000)
y1 <- y[1:900]
z <- autoarmafit(y1)
ar <- z$model[[1]]$arcoef
ma <- z$model[[1]]$macoef
var <- z$model[[1]]$v
y2 <- y[901:990]
prdctr(y2, r = 50, s = 90, h = 10, arcoef = ar, macoef = ma, v = var)

Rational Spectrum

Description

Compute power spectrum of ARMA process.

Usage

  raspec(h, var, arcoef = NULL, macoef = NULL, log = FALSE, plot = TRUE)

Arguments

Details

ARMA process :

y(t) - a(1)y(t-1) - \ldots - a(p)y(t-p) = u(t) - b(1)u(t-1) - \ldots - b(q)u(t-q)

where p is AR order, q is MA order and u(t) is a white noise with zero mean and variance equal to var.

Value

raspec gives the rational spectrum.

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1 for the AR model
raspec(h = 100, var = 1, arcoef = c(0.64,-0.8))

# Example 2 for the MA model
raspec(h = 20, var = 1, macoef = c(0.64,-0.8))

Frequency Response Function (Single Channel)

Description

Compute 1-input,1-output frequency response function, gain, phase, coherency and relative error statistics.

Usage

  sglfre(y, lag = NULL, invar, outvar)

Arguments

Value

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

ar <- array(0, dim = c(3,3,2))
ar[, , 1] <- matrix(c(0.4,  0,   0.3,
                      0.2, -0.1, -0.5,
                      0.3,  0.1,  0), nrow = 3, ncol = 3, byrow = TRUE)
ar[, , 2] <- matrix(c(0,  -0.3,  0.5,
                      0.7, -0.4,  1,
                      0,   -0.5,  0.3), nrow = 3, ncol = 3, byrow = TRUE)
x <- matrix(rnorm(200*3), nrow = 200, ncol = 3)
y <- mfilter(x, ar, "recursive")
sglfre(y, lag = 20, invar = 1, outvar = 2)

Optimal Controller Design and Simulation

Description

Produce optimal controller gain and simulate the controlled process.

Usage

simcon(span, len, r, arcoef, impulse, v, weight)

Arguments

Details

The basic state space model is obtained from the autoregressive moving average model of a vector process y(t);

y(t) - A(1)y(t-1) -\ldots- A(p)y(t-p) = u(t) - B(1)u(t-1) -\ldots- B(p-1)u(t-p+1),

where A(i) (i=1,\ldots,p) are the autoregressive coefficients of the ARMA representation of y(t).

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.6, Timsac74, A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Examples

x <- matrix(rnorm(1000*2), nrow = 1000, ncol = 2)
ma <- array(0, dim = c(2,2,2))
ma[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ma[, , 2] <- matrix(c( -0.2,  0.0,
                       -0.1, -0.3), nrow = 2, ncol = 2, byrow = TRUE)
y <- mfilter(x, ma, "convolution")

ar <- array(0, dim = c(2,2,3))
ar[, , 1] <- matrix(c( -1.0,  0.0,
                        0.0, -1.0), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 2] <- matrix(c( -0.5, -0.2,
                       -0.2, -0.5), nrow = 2, ncol = 2, byrow = TRUE)
ar[, , 3] <- matrix(c( -0.3, -0.05,
                       -0.1, -0.3), nrow = 2, ncol = 2, byrow = TRUE)
y <- mfilter(y, ar, "recursive")

z <- markov(y)
weight <-  matrix(c(0.0002,  0.0,
                    0.0,     2.9 ), nrow = 2, ncol = 2, byrow = TRUE)
simcon(span = 50, len = 700, r = 1, z$arcoef, z$impulse, z$v, weight)

Third Order Moments

Description

Compute the third order moments.

Usage

thirmo(y, lag = NULL, plot = TRUE)

Arguments

Value

References

H.Akaike, E.Arahata and T.Ozaki (1975) Computer Science Monograph, No.6, Timsac74, A Time Series Analysis and Control Program Package (2). The Institute of Statistical Mathematics.

Examples

data(bispecData)
z <- thirmo(bispecData, lag = 30)
z$tmomnt

Univariate Bayesian Method of AR Model Fitting

Description

This program fits an autoregressive model by a Bayesian procedure. The least squares estimates of the parameters are obtained by the householder transformation.

Usage

  unibar(y, ar.order = NULL, plot = TRUE)

Arguments

Details

The AR model is given by

y(t) = a(1)y(t-1) + \ldots + a(p)y(t-p) + u(t),

where p is AR order and u(t) is Gaussian white noise with mean 0 and variance v(p). The basic statistic AIC is defined by

AIC = n\log(det(v)) + 2m,

where n is the length of data, v is the estimate of innovation variance, and m is the order of the model.

Bayesian weight of the m-th order model is defined by

W(m) = CONST \times \frac{C(m)}{m+1},

where CONST is the normalizing constant and C(m)=\exp(-0.5AIC(m)). The equivalent number of free parameter for the Bayesian model is defined by

ek = D(1)^2 + \ldots + D(k)^2 +1,

where D(j) is defined by D(j)=W(j) + \ldots + W(k). m in the definition of AIC is replaced by ek to be define an equivalent AIC for a Bayesian model.

Value

References

H.Akaike (1978) A Bayesian Extension of The Minimum AIC Procedure of Autoregressive model Fitting. Research memo. No.126. The Institute of Statistical Mathematics.

G.Kitagawa and H.Akaike (1978) A Procedure for The Modeling of Non-Stationary Time Series. Ann. Inst. Statist. Math., 30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Canadianlynx)
z <- unibar(Canadianlynx, ar.order = 20)
z$arcoef

Univariate Case of Minimum AIC Method of AR Model Fitting

Description

This is the basic program for the fitting of autoregressive models of successively higher by the method of least squares realized through householder transformation.

Usage

unimar(y, max.order = NULL, plot = FALSE)

Arguments

Details

The AR model is given by

y(t) = a(1)y(t-1) + \ldots + a(p)y(t-p) + u(t),

where p is AR order and u(t) is Gaussian white noise with mean 0 and variance v. AIC is defined by

AIC = n\log(det(v)) + 2k,

where n is the length of data, v is the estimates of the innovation variance and k is the number of parameter.

Value

References

G.Kitagawa and H.Akaike (1978) A Procedure For The Modeling of Non-Stationary Time Series. Ann. Inst. Statist. Math.,30, B, 351–363.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

data(Canadianlynx)
z <- unimar(Canadianlynx, max.order = 20)
z$arcoef

White Noise Generator

Description

Generate approximately Gaussian vector white noise.

Usage

  wnoise(len, perr, plot = TRUE)

Arguments

Value

wnoise gives white noises.

References

H.Akaike and T.Nakagawa (1988) Statistical Analysis and Control of Dynamic Systems. Kluwer Academic publishers.

Examples

# Example 1
wnoise(len = 100, perr = 1)

# Example 2
v <- matrix(c(1,  0,  0,
              0,  2,  0,
              0,  0,  3), nrow = 3, ncol = 3, byrow = TRUE)
wnoise(len = 20, perr = v)

Exact Maximum Likelihood Method of Scalar ARMA Model Fitting

Description

Produce exact maximum likelihood estimates of the parameters of a scalar ARMA model.

Usage

  xsarma(y, arcoefi, macoefi)

Arguments

Details

The ARMA model is given by

y(t) - a(1)y(t-1) - \ldots - a(p)y(t-p) = u(t) - b(1)u(t-1) - ... - b(q)u(t-q),

where p is AR order, q is MA order and u(t) is a zero mean white noise.

Value

References

H.Akaike (1978) Covariance matrix computation of the state variable of a stationary Gaussian process. Research Memo. No.139. The Institute of Statistical Mathematics.

H.Akaike, G.Kitagawa, E.Arahata and F.Tada (1979) Computer Science Monograph, No.11, Timsac78. The Institute of Statistical Mathematics.

Examples

# "arima.sim" is a function in "stats".
# Note that the sign of MA coefficient is opposite from that in "timsac".
arcoef <- c(1.45, -0.9)
macoef <- c(-0.5)
y <- arima.sim(list(order=c(2,0,1), ar=arcoef, ma=macoef), n = 100)
arcoefi <- c(1.5, -0.8)
macoefi <- c(0.0)
z <- xsarma(y, arcoefi, macoefi)
z$arcoef
z$macoef