This fits an exponential smoothing state space model with model = 'ANN' to y,
having first performed classic multiplicative seasonal adjustment. A drift value is also calculated
by lsfit(0:(length(y) - 1), y)$coef[2] / 2. In combination with forecast.thetam(), this provides
identical results to forecast::thetaf(...). The purpose of splitting it into a `model` and
`forecast` functions is to make the approach consistent with other modeling / forecasting approaches
used in hybridModel().