# R/exposmooth.R In aTSA: Alternative Time Series Analysis

#### Documented in expsmooth

#' Simple Exponential Smoothing
#' @description Performs a simple exponential smoothing for univariate time series
#' with no trend or seasonal pattern.
#' @param x a numeric vector or univariate time series.
#' @param trend the type of trend. See details.
#' @param alpha the smoothing parameter for constant component. The default is \code{0.2}.
#' @param beta the smoothing parameter for linear component. The default is \code{0.10557}.
#' @param gamma the smoothing parameter for quadratic component. The default is \code{0.07168}.
#' @param lead the number of steps ahead for which prediction is required.
#' The default is \code{0}.
#' @param plot a logical value indicating to print the plot of original data v.s smoothed
#' data. The default is \code{TRUE}.
#' @details Simple exponential smoothing is a weighted average between the most recent
#' observation and the most recent forecasting, with weights \eqn{\alpha} and
#' \eqn{1 - \alpha}, respectively. To be precise, the smoothing equation of single exponential
#' smoothing (constant model, \code{trend = 1}) is given by
#' \deqn{level[t] = \alpha *x[t] + (1 - \alpha)*level[t-1],}
#' and the forecasting equation is
#' \deqn{hat{x}[t+1|t] = level[t],}
#' for \eqn{t = 1,...,n}.
#' The initial value \eqn{level[0] = x[1]}. For example, \eqn{hat{x}[1|0] = level[0]},
#' \eqn{hat{x}[2|1] = level[1]},..., etc.
#'
#' Let \eqn{x1[t]} be the smoothed values of single exponential smoothing. The double
#' exponential smoothing (\code{trend = 2}, a linear model) is to apply a single
#' exponential smoothing again to the smoothed sequence \eqn{x1[t]}, with a new smoothing
#' parameter \code{beta}. Similarly, we denote the smoothed values of double
#' exponential smoothing to be \eqn{x2[t]}. The triple exponential smoothing
#' (\code{trend = 3}, a quadratic model) is to apply the single exponential smoothing
#' to the smoothed sequence \eqn{x2[t]} with a new smoothing parameter \code{gamma}. The
#' default smoothing parameters (weights) \code{alpha}, \code{beta}, \code{gamma} are
#' taken from the equation \code{1 - 0.8^{1/trend}} respectively, which is similar
#' to the FORECAST procedure in SAS.
#'
#' @note Missing values are removed before the analysis.
#'
#' @return A list with class \code{"es"} containing the following components:
#' \item{estimate}{the smoothed values.}
#' \item{pred}{the predicted values when \code{lead} > 0.}
#' \item{accurate}{the accurate measurements.}
#' @author Debin Qiu
#' @examples x <- rnorm(100)
#' es <- expsmooth(x) # trend = 1: a constant model
#' plot(x,type = "l")
#' lines(es\$estimate,col = 2)
#' expsmooth(x,trend = 2) # trend = 2: a linear model
#' expsmooth(x,trend = 3) # trend = 3: a quadratic model
#' @importFrom stats ts
#' @importFrom stats frequency
#' @importFrom stats is.ts
#' @export
expsmooth <- function(x,trend = 1,alpha = 0.2,beta = 0.10557,
gamma = 0.07168,lead = 0,plot = TRUE)
{
if (NCOL(x) > 1L)
stop("'x' must be a univariate time series")
if (any(c(alpha,beta,gamma) > 1) || any(c(alpha,beta,gamma) < 0))
stop("'alpha' or 'beta' or 'gamma' must be between 0 and 1")
if (is.ts(x))
x <- ts(x[is.finite(x)],start = time(x)[1],frequency = frequency(x))
else
x <- x[is.finite(x)]
if (length(x) < 1L)
stop("'x' must have at least one observation")
single.expo <- function(x,parameter) {
m <- length(x)
x.hat <- c(x[1],numeric(m - 1))
for (i in 2:m)
x.hat[i] <- parameter*x[i-1] + (1 - parameter)*x.hat[i - 1]
return(x.hat)
}
if (trend == 1)
x.hat <- single.expo(x,alpha)
else if (trend == 2)
x.hat <- single.expo(single.expo(x,alpha),beta)
else if (trend == 3)
x.hat <- single.expo(single.expo(single.expo(x,alpha),beta),gamma)
else stop("'trend' must be one of 1,2,3")
if (is.ts(x))
x.hat <- ts(x.hat,start = time(x),frequency = frequency(x))
result <- list(estimate = x.hat)
stop("'lead' must be a positive number")
}
if (plot) {
plot(x,main = "original v.s smoothed data",type = "l")
lines(x.hat,col = 2)
}
result <- c(result,list(accurate = accurate(x,x.hat,trend,output = FALSE)))
class(result) <- "es"
return(result)
}


## Try the aTSA package in your browser

Any scripts or data that you put into this service are public.

aTSA documentation built on May 29, 2017, 11:44 a.m.