1. auto.arima在判断完d,D之后,调用Arima来求解AIC最小的解, 其中最优解的形式为arima(AR, d, MA)(SAR, D, SMA)[period]。
library(forecast)
fit <- auto.arima(data)
fit$arma = c(AR, MA, SAR, SMA, period, d, D)
2.有的时间序列数据需要差分(2)消除来消除非平稳性;其中的原理是:
若数据有逐年下降或者上升的趋势,那么一阶差分后的数据很可能是平稳的,因为此时的值近似于一阶导数,而“直线方程求导后斜率是稳定的”
(1)
Arima <- function (x, order = c(0L, 0L, 0L), seasonal = list(order = c(0L,
0L, 0L), period = NA), xreg = NULL, include.mean = TRUE,
transform.pars = TRUE, fixed = NULL, init = NULL, method = c("CSS-ML",
"ML", "CSS"), n.cond, SSinit = c("Gardner1980", "Rossignol2011"),
optim.method = "BFGS", optim.control = list(), kappa = 1e+06)
{
"%+%" <- function(a, b) .Call(C_TSconv, a, b)
SSinit <- match.arg(SSinit)
SS.G <- SSinit == "Gardner1980"
upARIMA <- function(mod, phi, theta) {
p <- length(phi)
q <- length(theta)
mod$phi <- phi
mod$theta <- theta
r <- max(p, q + 1L)
if (p > 0)
mod$T[1L:p, 1L] <- phi
if (r > 1L)
mod$Pn[1L:r, 1L:r] <- if (SS.G)
.Call(C_getQ0, phi, theta)
else .Call(C_getQ0bis, phi, theta, tol = 0)
else mod$Pn[1L, 1L] <- if (p > 0)
1/(1 - phi^2)
else 1
mod$a[] <- 0
mod
}
arimaSS <- function(y, mod) {
.Call(C_ARIMA_Like, y, mod, 0L, TRUE)
}
armafn <- function(p, trans) {
par <- coef
par[mask] <- p
trarma <- .Call(C_ARIMA_transPars, par, arma, trans)
if (is.null(Z <- tryCatch(upARIMA(mod, trarma[[1L]],
trarma[[2L]]), error = function(e) NULL)))
return(.Machine$double.xmax)
if (ncxreg > 0)
x <- x - xreg %*% par[narma + (1L:ncxreg)]
res <- .Call(C_ARIMA_Like, x, Z, 0L, FALSE)
s2 <- res[1L]/res[3L]
0.5 * (log(s2) + res[2L]/res[3L])
}
armaCSS <- function(p) {
par <- as.double(fixed)
par[mask] <- p
trarma <- .Call(C_ARIMA_transPars, par, arma, FALSE)
if (ncxreg > 0)
x <- x - xreg %*% par[narma + (1L:ncxreg)]
res <- .Call(C_ARIMA_CSS, x, arma, trarma[[1L]], trarma[[2L]],
as.integer(ncond), FALSE)
0.5 * log(res)
}
arCheck <- function(ar) {
p <- max(which(c(1, -ar) != 0)) - 1
if (!p)
return(TRUE)
all(Mod(polyroot(c(1, -ar[1L:p]))) > 1)
}
maInvert <- function(ma) {
q <- length(ma)
q0 <- max(which(c(1, ma) != 0)) - 1L
if (!q0)
return(ma)
roots <- polyroot(c(1, ma[1L:q0]))
ind <- Mod(roots) < 1
if (all(!ind))
return(ma)
if (q0 == 1)
return(c(1/ma[1L], rep.int(0, q - q0)))
roots[ind] <- 1/r