How to produce the minimum forecast error using R? Considering that we want to use optimize() on the interval [0,1] how can I write an R code for finding the value of β that produces the minimum forecast error without using external packages like forecast?


For simplicity assume that:

I want to use the following package:
> require(datasets)
> str(nhtemp)
 Time-Series [1:60] from 1912 to 1971: 49.9 52.3 49.4 51.1 49.4 47.9 49.8 50.9 49.3 51.9 ...

in which nhtemp is the Yearly Average Temperatures in New Haven CT.
 A: Two simple ways. Use the forecast package:
library(forecast)
fit <- ses(nhtemp,initial='simple')
beta <- fit$model$par[1]

Or use the stats package:
fit <- HoltWinters(nhtemp,gamma=FALSE,beta=FALSE)
beta <- fit$alpha

A: I hope you attempted to solve the HW problem yourself. See below for the complete code on how to do this using R and Optim function for your reference. I'm also providing you some nice references for you to try solving optimization problems in the context of exponential smoothing that would assist you to learn to solve optimization problems by hand. Based on my personal experience I have to tell that Optimization is a specialized field, you would want to take a course/reading on that before you start applying these techniques to forecasting problems.
There is an excellent article click here by Rasmussen that appeared in the Omega: The international journal of management science that would provide necessary back ground on how to optimize smoothing parameters as well as initial parameters using a nonlinear optimization tool like excel solver, you could try to replicate this yourself to learn how these parameters are optimized for exponential smoothing models.
If you  are interested in solving the problem that you stated in your post using R, as Rob pointed out you need to write a function to optimize the parameter beta and then using optim function to solve for that parameter.
You can also compare the beta values with a simple exponential smoothing in the forecast package. Both are similar (Beta value: 0.186088).
require(datasets)
str(nhtemp)

y <- as.matrix(nhtemp) ## Convert time series data to a matrix format

n = length(y)

## initialize
yh = rep(NA,60) ## Y hat
fe = rep(NA,60)
yh[[1]] = y[1] ## Initialize Yhat[1] with y [1]

## function to optimize

mind <- function (b){ for (i in 2:n)  {
                        yh[i] <- b*y[i-1]+(1-b)*yh[i-1]
                      } 
                      fe <- (yh-y)^2
                      obj <- mean(fe)
                      return(obj)
}
## Optimize to obtain parameter b

optm <- optim(par=c(0.2), fn=mind, gr = NULL, lower = 0, upper = 1,
              method = c("L-BFGS-B"))

## Print beta value

optm$par

## Check your answer with simple exponential smoothing from 
##forecast package from Rob's code above

library(forecast)
fit <- ses(nhtemp,initial='simple')
beta <- fit$model$par[1]
beta

