I'm trying to use the `optim`

function in R to find the optimal values for $alpha$ and $l_0$ for my custom residual sum of squares (RSS) function, but the values I'm getting back aren't even close to those picked by R's `ses`

function, which were $alpha$ = 0.2971 and $l_0$ = 77260, for the `pigs`

dataset in Rob Hyndman's `fpp2`

library.

Here's my custom RSS function:

`custom_rss <- function(ts, alpha, l0) { N <- length(ts) res <- vector() for (j in 1:N-1) { y_hat <- custom_ses(ts[1:j], alpha, l0) res[j] <- (ts[j+1] - y_hat)^2 } return(sum(res)) } `

My custom function for simple exponential smoothing is:

`custom_ses <- function(ts, alpha, l0) { N <- length(ts) y_hat = 0 for (j in 1:N-1) { y_hat <- y_hat + alpha * ((1 - alpha)^j) * ts[N-j] + ((1 - alpha)^N)*l0 } return(y_hat) } `

And finally, here's how I used the `optim`

function:

`optim(par=c(0,1), fn=custom_rss, alpha=0.2971, l0=77260.0561) `

So, I'm not sure where the problem lies: my SES function, my RSS function, or the way I'm using the `optim`

function.

**Contents**hide

#### Best Answer

I think your objective and smoothing functions are not correct.

More precisely, the problems I found are:

- The way you pass the parameters to the optimizer
- The absence of constraints on the parameters ($alpha in [0, 1])$
- The formula for y_hat. I used the Weighted Average Form in https://otexts.org/fpp2/ses.html#eq:7-ses

The parameters to optim() should be passed through the vector par. Then, from this vector, you can extract alpha and l0 in your objective function definition. Second, I think it is appropriate to estimate those parameter in such a way you are sure they will respect the necessary constraints. This can be achieved by transforming their values within the optimization routine.

Below you will find a modified example inspired by your R codes above

`rm(list = ls()); cat("14"); graphics.off() library(fpp2) # Objective function custom_rss = function(par, ts) { # Be seure to have alpha in [0,1] alpha = plogis(par[1]) l0 = par[2] N = length(ts) fit = custom_ses(ts, alpha, l0) res = ts - fit out = sum(res^2)/(N-1) return(out) } # Smoother function custom_ses = function(ts, alpha, l0) { N = length(ts) y_hat = c(l0, 0 * (2:N)) for (j in 1:(N-1)) { new = alpha * ts[j] + (1-alpha) * y_hat[j] y_hat[j+1] = new # print(y_hat) } return(y_hat) } # The data data("pigs") ts = pigs # optim-based estimates est = optim(par = c(0.1, 5000), fn = custom_rss, ts = ts) alpha = plogis(est$par[1]) l0 = est$par[2] fit = custom_ses(ts, alpha = alpha, l0 = l0) # Ses estimates sesMod = ses(ts) # Check cbind(ses = sesMod$model$par, custom = c(alpha, l0)) # Plots plot(1:length(ts), ts, type = "l", lwd = 1, col = 8) points(1:length(ts), ts, pch = 16, col = 8, cex = 0.5) lines(1:length(ts), sesMod$fitted, col = 2, type = "l", lwd = 2) lines(1:length(ts), fit, col = 3, type = "l", lty = 2, lwd = 2) legend("bottomright", legend = c("Data", "ses", "custom"), col = c(8, 2, 3), lty = c(1,1,2), lwd = 2) `

### Similar Posts:

- Solved – How to write log-likelihood for beta-geometric with optim() in R
- Solved – How to write log-likelihood for beta-geometric with optim() in R
- Solved – Using Nelder-Mead algorithm for minimizing the sum of absolute deviations from the residuals
- Solved – Calculating gradient of a function for optimization
- Solved – Calculating gradient of a function for optimization