| Title: | Butterworth-Induced Autoregressive Model |
|---|---|
| Description: | Implements the Butterworth-Induced Autoregressive ('BTWAR') model, where autoregressive coefficients are obtained from analog Butterworth filter prototypes mapped into the discrete-time domain using the Matched Z-Transform. The framework establishes a structured connection between frequency-domain filter design and time-domain autoregressive modeling. Model order selection is performed via nested rolling-origin cross-validation. Method described in Bras-Geraldes, Rocha and Martins (2026) <doi:10.3390/math14030479>. |
| Authors: | Carlos Bras-Geraldes [aut, cre, cph], J. Leonel Rocha [aut, cph] |
| Maintainer: | Carlos Bras-Geraldes <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0.1 |
| Built: | 2026-05-20 10:04:03 UTC |
| Source: | https://github.com/cgeraldes/btwar |
Determines the differencing order required to stationarise the training
series using make_stationary() and applies the same order to the
test series, ensuring both sets are on a comparable scale.
apply_stationarity(y_tr_raw, y_te_raw, max_d = 3L, alpha = 0.05)apply_stationarity(y_tr_raw, y_te_raw, max_d = 3L, alpha = 0.05)
y_tr_raw |
Numeric vector. Raw training time series. Must contain at least 4 observations. |
y_te_raw |
Numeric vector. Raw test (hold-out) time series. Must contain at least 1 observation. |
max_d |
Non-negative integer. Maximum differencing order passed to
|
alpha |
Numeric in |
This function is a convenience wrapper intended for use before fitting a
btwar_fit model when the user wants to inspect or control
the stationarity step explicitly. Internally, btwar_fit
calls make_stationary() directly.
A named list with four elements:
y_trNumeric vector. Stationarised training series
of length .
y_teNumeric vector. Differenced test series of
length .
dInteger. Number of differences applied.
stationaryLogical. TRUE if the ADF test
declared the training series stationary.
set.seed(1) y <- cumsum(rnorm(200)) # random walk (non-stationary) n_tr <- 160 split <- apply_stationarity( y_tr_raw = y[seq_len(n_tr)], y_te_raw = y[seq.int(n_tr + 1L, length(y))] ) split$d # differencing order applied split$stationary # was stationarity achieved?set.seed(1) y <- cumsum(rnorm(200)) # random walk (non-stationary) n_tr <- 160 split <- apply_stationarity( y_tr_raw = y[seq_len(n_tr)], y_te_raw = y[seq.int(n_tr + 1L, length(y))] ) split$d # differencing order applied split$stationary # was stationarity achieved?
Fits a Butterworth-Induced Autoregressive (BTWAR) model using nested rolling-origin cross-validation for order and attenuation parameter selection.
btwar_fit( y_tr_raw, y_te_raw, fs = 2, method = "ls", As_vec = seq(20, 60, by = 5), N_vec = 2:20, max_d = 3L, alpha_stationarity = 0.05, min_train = NULL, verbose = TRUE )btwar_fit( y_tr_raw, y_te_raw, fs = 2, method = "ls", As_vec = seq(20, 60, by = 5), N_vec = 2:20, max_d = 3L, alpha_stationarity = 0.05, min_train = NULL, verbose = TRUE )
y_tr_raw |
Numeric vector. Training time series. |
y_te_raw |
Numeric vector. Test (hold-out) time series. |
fs |
Positive numeric. Sampling frequency. Default |
method |
Character string. Scaling estimation method, one of
|
As_vec |
Numeric vector. Grid of stopband attenuation values (dB)
to search over. Default |
N_vec |
Integer vector. Grid of AR orders to search over.
Default |
max_d |
Non-negative integer. Maximum differencing order allowed
during the stationarity transformation step. Default |
alpha_stationarity |
Numeric in |
min_train |
Non-negative integer or |
verbose |
Logical. If |
An object of class "btwar", which is a named list
containing:
parametersList with d_used, N_opt,
A_opt, fc_opt, alpha, phi, and
poles_z.
performanceList with rmse_train and
rmse_test.
trainList with y_real and yhat for
the training set.
testList with y_real and yhat for
the test set.
muTraining mean used for centering.
callThe matched call.
set.seed(42) y <- cumsum(rnorm(900)) train_series <- y[1:600] test_series <- y[601:900] result <- btwar_fit( y_tr_raw = train_series, y_te_raw = test_series, fs = 2, method = "ls", N_vec = 2:5, As_vec = c(20, 40), verbose = FALSE ) summary(result)set.seed(42) y <- cumsum(rnorm(900)) train_series <- y[1:600] test_series <- y[601:900] result <- btwar_fit( y_tr_raw = train_series, y_te_raw = test_series, fs = 2, method = "ls", N_vec = 2:5, As_vec = c(20, 40), verbose = FALSE ) summary(result)
Extract AR Coefficients from a BTWAR Model
## S3 method for class 'btwar' coef(object, ...)## S3 method for class 'btwar' coef(object, ...)
object |
Object of class |
... |
Additional arguments (currently unused). |
Named numeric vector of AR coefficients.
Estimates the one-sided power spectrum of a time series using the Fast Fourier Transform (FFT). The series is mean-centred prior to transformation to remove the DC component.
compute_spectrum(x, fs = 1)compute_spectrum(x, fs = 1)
x |
Numeric vector. The input time series. Must contain at least 2 observations. |
fs |
Positive numeric. Sampling frequency in Hz (or consistent
units). Default |
Let denote the -th coefficient of the DFT of
. The raw periodogram is defined as
and the corresponding frequencies are
Only the non-redundant (one-sided) portion of the spectrum is returned,
i.e. frequencies from up to (but not including) the Nyquist
frequency .
A named list with two elements:
fNumeric vector of frequencies
( to , exclusive), length
.
PNumeric vector of power spectral estimates,
same length as f.
set.seed(1) x <- arima.sim(n = 256, model = list(ar = 0.8)) sp <- compute_spectrum(x, fs = 1) plot(sp$f, sp$P, type = "l", xlab = "Frequency (Hz)", ylab = "Power", main = "One-sided power spectrum")set.seed(1) x <- arima.sim(n = 256, model = list(ar = 0.8)) sp <- compute_spectrum(x, fs = 1) plot(sp$f, sp$P, type = "l", xlab = "Frequency (Hz)", ylab = "Power", main = "One-sided power spectrum")
Fitted Values from a BTWAR Model
## S3 method for class 'btwar' fitted(object, ...)## S3 method for class 'btwar' fitted(object, ...)
object |
Object of class |
... |
Additional arguments (currently unused). |
Numeric vector of fitted values on the training set.
Computes the mean squared error (MSE) between observed and predicted
values. MSE is commonly used for optimisation and theoretical analysis,
and is the square of rmse.
mse(y_true, y_pred)mse(y_true, y_pred)
y_true |
Numeric vector of observed (true) values. |
y_pred |
Numeric vector of predicted values. Must have the same
length as |
A single non-negative numeric value.
y <- c(1, 2, 3) yhat <- c(1.1, 1.9, 3.2) mse(y, yhat)y <- c(1, 2, 3) yhat <- c(1.1, 1.9, 3.2) mse(y, yhat)
Draws the Butterworth filter magnitude response (in dB) associated with
the selected BTW-AR model, using the optimal order and
cutoff frequency stored in the fitted object. Reference lines
are drawn at the cutoff frequency, the Nyquist frequency, and their
corresponding magnitude levels.
plot_bode( x, fs = 2, n_freq = 2000L, colour_magnitude = "blue", colour_cutoff = "orange", colour_nyquist = "red", lwd_magnitude = 1.2, lwd_ref = 1, base_size = 11, title = NULL )plot_bode( x, fs = 2, n_freq = 2000L, colour_magnitude = "blue", colour_cutoff = "orange", colour_nyquist = "red", lwd_magnitude = 1.2, lwd_ref = 1, base_size = 11, title = NULL )
x |
Object of class |
fs |
Positive numeric. Sampling frequency (Hz) used to derive the
Nyquist frequency as |
n_freq |
Positive integer. Number of frequency points used to
evaluate the magnitude response. Default |
colour_magnitude |
Character string. Colour for the magnitude
response curve. Default |
colour_cutoff |
Character string. Colour for the cutoff frequency
reference lines (vertical dashed and horizontal dotted).
Default |
colour_nyquist |
Character string. Colour for the Nyquist frequency
reference lines (vertical dashed and horizontal dotted).
Default |
lwd_magnitude |
Positive numeric. Line width for the magnitude
response curve. Default |
lwd_ref |
Positive numeric. Line width for all reference lines.
Default |
base_size |
Positive numeric. Base font size passed to
|
title |
Character string. Plot title. If |
The Butterworth filter magnitude response is
evaluated on a log-spaced frequency axis from Hz to the
Nyquist frequency. The magnitude at the cutoff frequency equals
dB (), and the magnitude at the Nyquist
frequency reflects the actual stopband attenuation achieved by the
selected filter.
A ggplot object. The plot is also printed
as a side effect when called interactively.
btwar_fit, plot.btwar,
plot_zpoles
sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_bode(fit, fs = 12000)sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_bode(fit, fs = 12000)
Displays the one-sided amplitude spectrum of a numeric series, computed via the Fast Fourier Transform (FFT), with reference lines at the cutoff frequency and the Nyquist frequency.
plot_freq_amplitude( x, fs, fc, colour_spectrum = "blue", colour_cutoff = "red", colour_nyquist = "black", lwd_spectrum = 0.6, lwd_ref = 1, base_size = 11, title = NULL )plot_freq_amplitude( x, fs, fc, colour_spectrum = "blue", colour_cutoff = "red", colour_nyquist = "black", lwd_spectrum = 0.6, lwd_ref = 1, base_size = 11, title = NULL )
x |
Numeric vector. The time series to analyze. The mean is removed internally before computing the FFT. |
fs |
Positive numeric. Sampling frequency (Hz). Must match the
value passed to |
fc |
Positive numeric. Cutoff frequency (Hz) to display as a
reference line. Typically read from |
colour_spectrum |
Character string. Colour for the amplitude
spectrum segments. Default |
colour_cutoff |
Character string. Colour for the cutoff frequency
reference line. Default |
colour_nyquist |
Character string. Colour for the Nyquist frequency
reference line. Default |
lwd_spectrum |
Positive numeric. Line width for the amplitude
spectrum segments. Default |
lwd_ref |
Positive numeric. Line width for the reference lines.
Default |
base_size |
Positive numeric. Base font size passed to
|
title |
Character string. Plot title. Default |
A ggplot object. The plot is also printed
as a side effect when called interactively.
btwar_fit, plot.btwar,
plot_bode
sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_freq_amplitude(as.numeric(sim$train), fs = 12000, fc = fit$parameters$fc_opt)sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_freq_amplitude(as.numeric(sim$train), fs = 12000, fc = fit$parameters$fc_opt)
Displays the poles of the selected BTW-AR model in the complex Z-plane, together with the unit circle. Optionally overlays poles from one or more external sources (e.g., a fitted ARIMA model or the true data-generating poles).
plot_zpoles( x, external_list = NULL, colour_selected = "green", external_colours = c("red", "darkgreen", "purple", "orange", "brown", "navy"), external_shapes = c(17L, 15L, 18L, 8L, 10L, 12L), lim = 1.5, base_size = 10, title = NULL )plot_zpoles( x, external_list = NULL, colour_selected = "green", external_colours = c("red", "darkgreen", "purple", "orange", "brown", "navy"), external_shapes = c(17L, 15L, 18L, 8L, 10L, 12L), lim = 1.5, base_size = 10, title = NULL )
x |
Object of class |
external_list |
Named list of external pole sets to overlay.
Each element must be a |
colour_selected |
Character string. Colour for the BTW-AR
(selected) poles. Default |
external_colours |
Character vector. Colours assigned to external
pole sets in order. Recycled if fewer colours than sets are provided.
Default |
external_shapes |
Integer vector. |
lim |
Positive numeric. Half-width of the plot window on both axes.
Default |
base_size |
Positive numeric. Base font size passed to
|
title |
Character string. Plot title. If |
Poles of the selected BTW-AR model are extracted from
x$parameters$poles_z, which is populated by btwar_fit
via the internal cross-validation routine. The unit circle is drawn as a
reference: poles inside the circle correspond to stable filters.
Colours and shapes for external pole sets are assigned sequentially from
external_colours and external_shapes, with circular
recycling when the number of sets exceeds the palette length.
A ggplot object. The plot is also printed
as a side effect when called interactively.
sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_zpoles(fit) df_true <- data.frame(Re = Re(sim$poles), Im = Im(sim$poles)) plot_zpoles(fit, external_list = list("True Signal" = df_true))sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot_zpoles(fit) df_true <- data.frame(Re = Re(sim$poles), Im = Im(sim$poles)) plot_zpoles(fit, external_list = list("True Signal" = df_true))
Generates a ggplot2 line chart of fitted or test-set predictions
from a "btwar" object. Optionally overlays the observed series
and/or an external comparison series (e.g., a benchmark model).
## S3 method for class 'btwar' plot( x, dataset = c("train", "test"), fs = 1, show_observed = TRUE, external = NULL, external_name = "External", title = NULL, colour_btwar = "#2166AC", colour_observed = "#333333", colour_external = NULL, lwd_btwar = 0.9, lwd_observed = 0.55, lwd_external = 0.7, alpha = 0.9, base_size = 11, palette = "Dark 2", layer_order = NULL, ... )## S3 method for class 'btwar' plot( x, dataset = c("train", "test"), fs = 1, show_observed = TRUE, external = NULL, external_name = "External", title = NULL, colour_btwar = "#2166AC", colour_observed = "#333333", colour_external = NULL, lwd_btwar = 0.9, lwd_observed = 0.55, lwd_external = 0.7, alpha = 0.9, base_size = 11, palette = "Dark 2", layer_order = NULL, ... )
x |
Object of class |
dataset |
Character string. Which dataset to display: |
fs |
Positive numeric. Sampling frequency used to construct the time
axis as |
show_observed |
Logical. If |
external |
Optional numeric vector of the same length as the selected
dataset. If supplied, it is plotted as an additional comparison series.
Default |
external_name |
Character string. Legend label for the external
series. Default |
title |
Character string. Plot title. If |
colour_btwar |
Character string. Colour for the BTWAR prediction
line. Default |
colour_observed |
Character string. Colour for the observed series
line. Default |
colour_external |
Character string. Colour for the external series
line. If |
lwd_btwar |
Positive numeric. Line width for the BTWAR prediction
line. Default |
lwd_observed |
Positive numeric. Line width for the observed series
line. Default |
lwd_external |
Positive numeric. Line width for any external series
line. Default |
alpha |
Numeric in |
base_size |
Positive numeric. Base font size passed to
|
palette |
Character string. Name of the |
layer_order |
Character vector. Controls the drawing order (back to
front) and legend order of the series. Must contain the names of all
active series: |
... |
Additional arguments (currently unused). |
Series are drawn in the following fixed order (back to front): Observed, BTWAR, external.
A ggplot object. The plot is also printed
as a side effect when called interactively.
btwar_fit, plot_zpoles,
fitted.btwar, residuals.btwar
sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot(fit, dataset = "train", fs = 12000) plot(fit, dataset = "test", fs = 12000, show_observed = TRUE)sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) plot(fit, dataset = "train", fs = 12000) plot(fit, dataset = "test", fs = 12000, show_observed = TRUE)
Returns the poles of an autoregressive (AR) model in the complex Z-plane
by finding the roots of the characteristic polynomial
.
poles_AR(phi)poles_AR(phi)
phi |
Numeric vector of AR coefficients
|
The characteristic polynomial is
which, after multiplication by , becomes the polynomial whose
roots are computed by polyroot. Coefficients are
reversed internally so that polyroot receives them in ascending
degree order.
A complex vector of length containing the Z-plane poles.
A model is stable if and only if all poles lie strictly inside the unit
circle, i.e., all(Mod(poles_AR(phi)) < 1).
yhat_ar, yhat_arma,
plot_zpoles, polyroot
# AR(2) with phi = c(0.6, -0.3) phi <- c(0.6, -0.3) poles <- poles_AR(phi) print(poles) all(Mod(poles) < 1) # TRUE => stable # Use with plot_zpoles sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) df <- data.frame(Re = Re(poles), Im = Im(poles)) plot_zpoles(fit, external_list = list("AR(2) poles" = df))# AR(2) with phi = c(0.6, -0.3) phi <- c(0.6, -0.3) poles <- poles_AR(phi) print(poles) all(Mod(poles) < 1) # TRUE => stable # Use with plot_zpoles sim <- simulate_ar_split(phi_real = c(0.34, -0.22, 0.16), n = 1000, burn = 200, prop_train = 0.7) fit <- btwar_fit(y_tr_raw = as.numeric(sim$train), y_te_raw = as.numeric(sim$test), fs = 12000, method = "ls", N_vec = 3:3) df <- data.frame(Re = Re(poles), Im = Im(poles)) plot_zpoles(fit, external_list = list("AR(2) poles" = df))
Generates one-step-ahead predictions for a new time series using a
fitted "btwar" model.
## S3 method for class 'btwar' predict(object, newdata, ...)## S3 method for class 'btwar' predict(object, newdata, ...)
object |
Object of class |
newdata |
Numeric vector. New time series to predict from. |
... |
Additional arguments (currently unused). |
Numeric vector of predicted values.
Compact console display of a fitted "btwar" object.
## S3 method for class 'btwar' print(x, ...)## S3 method for class 'btwar' print(x, ...)
x |
Object of class |
... |
Additional arguments (currently unused). |
Invisibly returns x.
Print a BTWAR Model Summary
## S3 method for class 'summary.btwar' print(x, ...)## S3 method for class 'summary.btwar' print(x, ...)
x |
Object of class |
... |
Additional arguments (currently unused). |
Invisibly returns x.
Residuals from a BTWAR Model
## S3 method for class 'btwar' residuals(object, ...)## S3 method for class 'btwar' residuals(object, ...)
object |
Object of class |
... |
Additional arguments (currently unused). |
Numeric vector of residuals from the training set.
Computes the root mean squared error (RMSE) between observed and predicted values. RMSE penalises large deviations more strongly than MAE due to squaring, and is expressed in the same units as the response variable.
rmse(y_true, y_pred)rmse(y_true, y_pred)
y_true |
Numeric vector of observed (true) values. |
y_pred |
Numeric vector of predicted values. Must have the same
length as |
A single non-negative numeric value.
y <- c(1, 2, 3) yhat <- c(1.1, 1.9, 3.2) rmse(y, yhat)y <- c(1, 2, 3) yhat <- c(1.1, 1.9, 3.2) rmse(y, yhat)
Generates a stationary AR(p) time series via arima.sim,
discards a burn-in period, and partitions the result into a training set and
a hold-out test set.
simulate_ar_split( phi_real, n = 2000L, burn = 200L, prop_train = 0.8, seasonal = FALSE, freq = NA )simulate_ar_split( phi_real, n = 2000L, burn = 200L, prop_train = 0.8, seasonal = FALSE, freq = NA )
phi_real |
Numeric vector of true AR coefficients
|
n |
Integer. Number of observations to retain after burn-in.
Default |
burn |
Integer. Number of initial observations discarded as burn-in.
Default |
prop_train |
Numeric in |
seasonal |
Logical. If |
freq |
Integer or |
Random number generation is controlled externally by the user via
set.seed when reproducibility is required.
The function does not modify the global random number generator state.
For reproducible results, users should call set.seed() prior
to invoking this function.
A named list with four elements:
seriesThe full simulated series (length n).
trainTraining observations (first
floor(n * prop_train) values).
testHold-out test observations (remaining values).
polesComplex vector of true AR poles in the Z-plane.
# Reproducible AR(3) with phi = (0.6, -0.3, 0.2) set.seed(123) result <- simulate_ar_split(phi_real = c(0.6, -0.3, 0.2)) length(result$train) # 1600 length(result$test) # 400 # AR(3) as a monthly time series, 90/10 split set.seed(123) result2 <- simulate_ar_split( phi_real = c(0.6, -0.3, 0.2), n = 1200, prop_train = 0.9, seasonal = TRUE, freq = 12 )# Reproducible AR(3) with phi = (0.6, -0.3, 0.2) set.seed(123) result <- simulate_ar_split(phi_real = c(0.6, -0.3, 0.2)) length(result$train) # 1600 length(result$test) # 400 # AR(3) as a monthly time series, 90/10 split set.seed(123) result2 <- simulate_ar_split( phi_real = c(0.6, -0.3, 0.2), n = 1200, prop_train = 0.9, seasonal = TRUE, freq = 12 )
Returns a "summary.btwar" object with model parameters and
performance metrics. The associated print method displays a
formatted summary including the AR coefficients.
## S3 method for class 'btwar' summary(object, ...)## S3 method for class 'btwar' summary(object, ...)
object |
Object of class |
... |
Additional arguments (currently unused). |
An object of class "summary.btwar" (invisibly).
Computes one-step-ahead predictions for a time series under a pure
autoregressive model of order .
yhat_ar(x, phi)yhat_ar(x, phi)
x |
Numeric vector. The observed time series. |
phi |
Numeric vector of AR coefficients
|
No intercept is included. Center x before calling this function
if the series has a non-zero mean.
A numeric vector of the same length as x. The first
elements are NA (insufficient history); element
() contains
.
poles_AR, yhat_arma,
btwar_fit
x <- as.numeric(arima.sim(list(ar = c(0.6, -0.3)), n = 200)) phi <- c(0.6, -0.3) yh <- yhat_ar(x, phi) # First p values are NA head(yh, 5) # RMSE on the predictable portion obs <- x[(length(phi) + 1):length(x)] hat <- yh[(length(phi) + 1):length(x)] sqrt(mean((obs - hat)^2))x <- as.numeric(arima.sim(list(ar = c(0.6, -0.3)), n = 200)) phi <- c(0.6, -0.3) yh <- yhat_ar(x, phi) # First p values are NA head(yh, 5) # RMSE on the predictable portion obs <- x[(length(phi) + 1):length(x)] hat <- yh[(length(phi) + 1):length(x)] sqrt(mean((obs - hat)^2))
Computes one-step-ahead predictions for a time series under an
autoregressive moving-average (ARMA) model of orders .
yhat_arma(x, ar, ma)yhat_arma(x, ar, ma)
x |
Numeric vector. The observed time series. |
ar |
Numeric vector of AR coefficients
|
ma |
Numeric vector of MA coefficients
|
Residuals prior to the start index are initialised to zero. No intercept
is included; center x before calling if the series has a non-zero
mean. For a pure AR model, prefer yhat_ar, which is
slightly more efficient.
A numeric vector of the same length as x. Elements
are NA; element contains
where residuals are accumulated
recursively and initialised to zero.
x <- as.numeric(arima.sim(list(ar = 0.6, ma = 0.4), n = 200)) yh <- yhat_arma(x, ar = 0.6, ma = 0.4) # First max(p, q) values are NA head(yh, 5) # RMSE on the predictable portion start <- max(length(0.6), length(0.4)) + 1L sqrt(mean((x[start:length(x)] - yh[start:length(x)])^2))x <- as.numeric(arima.sim(list(ar = 0.6, ma = 0.4), n = 200)) yh <- yhat_arma(x, ar = 0.6, ma = 0.4) # First max(p, q) values are NA head(yh, 5) # RMSE on the predictable portion start <- max(length(0.6), length(0.4)) + 1L sqrt(mean((x[start:length(x)] - yh[start:length(x)])^2))