//:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** @author Hao Peng, John Miller * @version 1.6 * @date Sat Jun 13 01:27:00 EST 2017 * @see LICENSE (MIT style license file). * * @see http://en.wikipedia.org/wiki/Autoregressive%E2%80%93moving-average_model * @see http://www.emu.edu.tr/mbalcilar/teaching2007/econ604/lecture_notes.htm * @see http://www.stat.berkeley.edu/~bartlett/courses/153-fall2010 * @see http://www.princeton.edu/~apapanic/ORFE_405,_Financial_Time_Series_%28Fall_2011%29_files/slides12-13.pdf * * @title Seasonal Auto-Regressive, Integrated, Moving Average (SARIMA) Model */ package scalation.analytics package forecaster import scala.math.{max, sqrt} import scalation.linalgebra._ import scalation.math.double_exp import scalation.minima.QuasiNewton import scalation.plot.Plot import scalation.random.{Normal, Uniform} import scalation.util.banner import ForecasterVec.testInSamp //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Companion object for class `SARIMA`. Includes features related to differencing * and automated order selection. * @see www.jstatsoft.org/article/view/v027i03/v27i03.pdf */ object SARIMA { //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Difference the time series. Return both the completely differenced time series * and the intermediate (differenced once) time series (needed to scale back results later). * @param y the time series to be differenced * @param d the order of simple differencing * @param dd the order of seasonal differencing * @param period the seasonal period */ def difference (y: VectoD, d: Int, dd: Int, period: Int): (VectoD, VectoD) = { val y_ = ARIMA.difference (y, d) (differenceSeason (y_, dd, period), y_) } // difference //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Difference for seasonality. * @param y_ the time series to be seasonally differenced, * this is usually the intermediate result after simple differencing. * @param dd the order of seasonal differencing * @param period the seasonal period */ def differenceSeason (y_ : VectoD, dd: Int, period: Int): VectoD = { dd match { case 0 => y_ case 1 => VectorD (for (i <- 0 until y_.dim-period) yield y_(i+period) - y_(i)) case 2 => VectorD (for (i <- 0 until y_.dim-2*period) yield y_(i+2*period) - 2*y_(i+period) + y_(i)) case 3 => VectorD (for (i <- 0 until y_.dim-3*period) yield y_(i+3*period) - 3*y_(i+2*period) + 3*y_(i+period) - y_(i)) case _ => flaw ("differenceSeason", "SARIMA does not support seasonal differencing higher than order 3"); null } // match } // differenceSeason //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Transform the fitted values on the training data of a differenced time series back * to the original scale. * @see stats.stackexchange.com/questions/32634/difference-time-series-before-arima-or-within-arima * @param xp the vector of predictions/fitted values of a differenced time series * @param x_ the intermediate result after differencing for trend, but before differencing for seasonality * @param y the original time series * @param d the order of simple differencing * @param dd the order of seasonal differencing * @param period the seasonal period */ def transformBack (xp: VectoD, x_ : VectoD, y: VectoD, d: Int, dd: Int, period: Int): VectoD = { val xp_ = transformBackSeason (xp, x_, dd, period) ARIMA.transformBack (xp_, y, d) } // transformBack //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Transform the fitted values on the training data of a differenced time series back * to the original scale. Undo seasonal differencing only. * @see stats.stackexchange.com/questions/32634/difference-time-series-before-arima-or-within-arima * @param xp the vector of predictions/fitted values of a differenced time series * @param x_ the intermediate result after differencing for trend, but before differencing for seasonality * @param dd the order of seasonal differencing * @param period the seasonal period */ def transformBackSeason (xp: VectoD, x_ : VectoD, dd: Int, period: Int): VectoD = { dd match { case 0 => xp case 1 => val tb = new VectorD (xp.dim + period) for (i <- 0 until period) tb(i) = x_(i) for (i <- 0 until x_.dim-period) tb(i+period) = xp(i) + x_(i) tb case 2 => val tb = new VectorD (xp.dim + 2*period) for (i <- 0 until 2*period) tb(i) = x_(i) for (i <- 0 until x_.dim-2*period) tb(i+2*period) = xp(i) + 2*x_(i+period) - x_(i) tb case 3 => val tb = new VectorD (xp.dim + 3*period) for (i <- 0 until 3*period) tb(i) = x_(i) for (i <- 0 until x_.dim-3*period) tb(i+3*period) = xp(i) + 3*x_(i+2*period) - 3*x_(i+period) + x_(i) tb case _ => flaw ("transformBackSeason", "SARIMA does not support seasonal differencing higher than order 3"); null } // match } // transformBackSeason //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Transform the forecast values of a differenced time series back to the * original scale. * @see stats.stackexchange.com/questions/32634/difference-time-series-before-arima-or-within-arima * @param xf the vector of forecasted values of a differenced time series * @param x_ the intermediate result after differencing for trend, but before differencing for seasonality * @param xx the original zero-mean time series * @param d the order of simple differencing * @param dd the order of seasonal differencing * @param period the seasonal period * @param t the time being forecasted (@see the 'forecast' method) */ def transformBackF (xf: VectoD, x_ : VectoD, xx: VectoD, d: Int, dd: Int, period: Int, t: Int): VectoD = { val xf_ = transformBackFSeason (xf, x_, dd, period, t-d) ARIMA.transformBackF (xf_, xx, d, t) } // transformBackF //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Transform the forecast values of a differenced time series back to the * original scale. Undo seasonal differencing only. * @see stats.stackexchange.com/questions/32634/difference-time-series-before-arima-or-within-arima * @param xf the vector of forecasted values of a differenced time series * @param x_ the intermediate result after differencing for trend, but before differencing for seasonality * @param dd the order of seasonal differencing * @param period the seasonal period * @param t the time point being forecasted (@see the 'forecast' method) */ def transformBackFSeason (xf: VectoD, x_ : VectoD, dd: Int, period: Int, t: Int): VectoD = { dd match { case 0 => xf case 1 => val tb = x_.slice (t - period, t) ++ xf for (i <- period until tb.dim) tb(i) += tb(i-period) tb.slice (period) case 2 => val tb = x_.slice (t - 2*period, t) ++ xf for (i <- 2*period until tb.dim) tb(i) += (2*tb(i-period) - tb(i-2*period)) tb.slice (2*period) case 3 => val tb = x_.slice (t - 3*period, t) ++ xf for (i <- 3*period until tb.dim) tb(i) += (3*tb(i-period) - 3*tb(i-2*period) + tb(i-3*period)) tb.slice (3*period) case _ => flaw ("transformBackFSeason", "SARIMA does not support seasonal differencing higher than order 3"); null } // match } // transformBackFSeason } // SARIMA object import SARIMA._ //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMA` class provides basic time series analysis capabilities for Auto- * Regressive 'AR' Integrated 'I' Moving-Average 'MA' models. In a * 'SARIMA(p, d, q)' model, 'p' and 'q' refer to the order of the Auto-Regressive * and Moving-Average components of the model; 'd' refers to the order of * differencing. Given time series data stored in vector 'y', its next value 'y_t = y(t)' * may be predicted based on prior values of 'y' and its noise: *

* y_t = δ + Σ(φ_i y_t-i) + Σ(θ_i e_t-i) + e_t *

* where 'δ' is a constant, 'φ' is the auto-regressive coefficient vector, * 'θ' is the moving-average coefficient vector, and 'e' is the noise vector. *------------------------------------------------------------------------------ * If 'd > 0', then the time series must be differenced first before applying * the above model. Seasonal differencing, auto-regressive and moving average * factors can be incorporated into the model by applying seasonal differencing * (possibly in addition to simple differencing) first, then add the seasonal * auto-regressive and moving average terms, that rely on lagged values and errors, * respectively, from one or more seasonal periods in the past, on the right hand * side of the equation. *------------------------------------------------------------------------------ * @param y the original input vector (time series data) * @param d the order of Integration/simple differencing * @param dd the order of seasonal differencing * @param period the seasonal period (at least 2) */ class SARIMA (y: VectoD, d: Int = 0, dd: Int = 0, period: Int = 2) extends ARIMA (y, d) { if (period < 2) flaw ("SARIMA", "the seasonal period must be at least 2") private val DEBUG = true // debug flag private var pp = 0 // seasonal AR order private var qq = 0 // seasonal MA order private var z_ : VectoD = null // intermediate results after simple differencing private var φφ: VectoD = null // seasonal AR(pp) coefficients private var θθ: VectoD = null // seasonal MA(qq) coefficients differenced = d > 0 || dd > 0 // flag indicating whether differencing will be applied init (y) // initialize vectors and parameters //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Initialize variables based on the working time-series 'v'. * Set/change the working time series. May be used to set the time series * to a different time window in order to produce newer forecast. * @param v the working vector/time-series */ override def init (v: VectoD): Unit = { mu = v.mean // sample mean val zz = difference (v, d, dd, period) // difference (simple/seasonal) the time series z = zz._1 // time series after differencing z_ = zz._2 // intermediate results after simple differencing zp = new VectorD (z.dim) // predicted values prior to undifferencing/uncentering e = new VectorD (z.dim) // vector of errors/residuals sig2 = z.variance // sample variance } // init //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Return the model name including its current hyper-parameter. */ override def modelName: String = { if (dd > 0) s"SARIMA ($p, $d, $q) x ($pp, $dd, $qq)_${period}" else s"SARIMA ($p, $d, $q)" } // modelName //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Return the parameter vector (concatenation of φ, θ, φφ and θθ). */ override def parameter: VectoD = φ ++ θ ++ φφ ++ θθ //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Set values for the models orders 'p', 'q', 'pp' and 'qq'. * @param pq the vector of model orders */ override def setPQ (pq: VectoI): SARIMA = { val n = pq.dim if (n > 0) p = pq(0) φ = new VectorD (p) // AR coefficients if (n > 1) q = pq(1) θ = new VectorD (q) // MA coefficients if (n > 2) pp = pq(2) φφ = new VectorD (pp) // seasonal AR coefficients if (n > 3) qq = pq(3) θθ = new VectorD (qq) // seasonal MA coefficients cap = Array (p, q, pp*period, qq*period).max // greatest lag params = p + q + pp + qq + (if (differenced) 0 else 1) // number of parameters this } // setPQ //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Train/fit an `SARIMA` model to the times series data. Must call 'setPQ' first. * It uses BFGS, a Quasi-Newton optimizer, to minimize the negative log-likelihood. */ override def train (): SARIMA = { val solver = new QuasiNewton (nll) // nonlinear optimizer val b = new VectorD (params + 1) // parameter values if (! differenced) b(b.size-2) = mu // sample mean, initial est. for μ parameter b(b.size-1) = sqrt (sig2) // sample standard deviation, initial est. for σ parameter solver.solve (b) // find b that maximizes likelihood δ = μ * (1 - φ.sum) // update drift value // δ = stats.mu * (1 - φ.sum) if (DEBUG) { showParameterEstimates () // show ARMA parameters println (s"θθ = $θθ") // seasonal moving average coefficients println (s"φφ = $φφ") // seasonal auto-regressive coefficients } // if this } // train //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The negative log-likelihood function to be minimized. * @see spia.uga.edu/faculty_pages/monogan/teaching/ts/Barima.pdf * @see stats.stackexchange.com/questions/77663/arima-estimation-by-hand * @param input the input parameter vector */ override def nll (b: VectoD): Double = { if (b.size != params + 1) flaw ("nll", "input parameter vector size incorrect") for (i <- 0 until p) φ(i) = b(i) for (i <- p until p+pp) φφ(i-p) = b(i) for (i <- p+pp until p+pp+q) θ(i-p-pp) = b(i) for (i <- p+pp+q until p+pp+q+qq) θθ(i-p-pp-q) = b(i) if (! differenced) μ = b(b.size-2) σ2 = b.last~^2 updateFittedValues () } // nll //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Update the vector of fitted values 'zp', the vector of errors 'e', and * return the negative log-likelihood '-ll'. * @see `Fit` for definition of 'll'. */ override def updateFittedValues (): Double = { if (! differenced) for (i <- z.range) z(i) = y(i) - μ // for undifferenced time series, center using est. μ zp(0) = z(0) // no past values or errors => copy actual for (t <- 1 until zp.dim) { e(t-1) = z(t-1) - zp(t-1) // error in previous forecast var sum = 0.0 for (j <- 0 until p if t-j > 0) sum += φ(j) * z(t-1-j) for (j <- 0 until pp if t-(1+j)*period >= 0) sum += φφ(j) * z(t-(1+j)*period) for (j <- 0 until q if t-j > 0) sum += θ(j) * e(t-1-j) for (j <- 0 until qq if t-(1+j)*period >= 0) sum += θθ(j) * e(t-(1+j)*period) zp(t) = sum } // for -ll (e.normSq / m, σ2, m) // return negative log likelihood } // updateFittedValues //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Return the vector of predicted/fitted values on the training/full data. * Based on 'zp' calculated in the 'updateFittedValues' method. */ override def predictAll (): VectoD = { if (differenced) transformBack (zp, z_, y, d, dd, period) else zp + μ } //predictAll //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** Produce h-steps-ahead forecast for SARIMA models. * @see ams.sunysb.edu/~zhu/ams586/Forecasting.pdf * @param t the time point from which to make forecasts (in the original scale) * @param h the number of steps to forecast, must be at least one */ override def forecast (t: Int = y.dim, h: Int = 1): VectoD = { if (t > y.dim) flaw ("forecast", s"t ($t) cannot be greater than y.dim (${y.dim})") val tz = t - d - dd * period // scale t to match vector z and e if (tz < cap) flaw ("forecast", s"tz ($tz) must be at least cap ($cap)") val zf = new VectorD (cap + h) // forecasted centered values val e_ = new VectorD (cap + h) // available observed errors for (i <- 0 until cap if tz-cap+i >= 0) { // seed with first cap values zf(i) = z(tz-cap+i) // copy first cap values e_(i) = e(tz-cap+i) // unveil first cap errors (observed in training) } // for for (i <- cap until zf.dim) { // start at t = cap (enough for first value to forecast) var sum = 0.0 for (j <- 0 until p) sum += φ(j) * zf(i-1-j) for (j <- 0 until pp) sum += φφ(j) * zf(i-(1+j)*period) for (j <- 0 until q) sum += θ(j) * e_(i-1-j) for (j <- 0 until qq) sum += θθ(j) * e_(i-(1+j)*period) zf(i) = sum } // for val f = zf.slice (cap) // dump first cap values if (differenced) transformBackF (f, z_, y, d, dd, period, t) else f + μ // return the vector of forecasts } // forecast } // SARIMA class //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest` object is used to test the `SARIMA` class. * > runMain scalation.analytics.forecaster.SARIMATest */ object SARIMATest extends App { println ("SARIMA!") val m = 100 val noise = Uniform (-5, 5) val y = VectorD (for (i <- 0 until m) yield i + noise.gen) val ts = new SARIMA (y, 1) // time series data: y vs. t, apply 1st order differencing // ts.plotFunc2 (ts.acf, "ACF") // must turn on DEBUG so that ACF is actually computed // Build SARIMA(1, 0, 0), SARIMA(2, 0, 0) and SARIMA(0, 0, 3) models ts.setPQ (VectorI (1, 0)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll (), "Plot of y, SARIMA(1, 0, 0) vs. t", true) ts.setPQ (VectorI (2, 0)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll (), "Plot of y, SARIMA(2, 0, 0) vs. t", true) ts.setPQ (VectorI (0, 3)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll () , "Plot of y, SARIMA(0, 0, 3) vs. t", true) println (ts.forecastAll(3)) } // SARIMATest object //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest2` object is used to test the `SARIMA` class. * > runMain scalation.analytics.forecaster.SARIMATest2 */ object SARIMATest2 extends App { val m = 30 val noise = Normal (0, 2) val y = VectorD (for (i <- 0 until m) yield i + noise.gen) println (s"y = $y") val p = 1 val d = 1 // order of differencing val q = 1 val steps = 2 // number of steps for the forecasts val ts = new SARIMA (y, d) // time series data: y vs. t // Build ARIMA(p, 0), ARIMA(0, q) and ARIMA(p, q) models for the (differenced) time series data ts.setPQ (VectorI (p, 0)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll (), s"Plot of y, SARIMA($p, 0, 0) vs. t", true) val ar_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA($p, 0, 0) model = $ar_f") ts.setPQ (VectorI (0, q)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll (), s"Plot of y, SARIMA(0, 0, $q) vs. t", true) val ma_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA(0, 0, $q) model = $ma_f") ts.setPQ (VectorI (p, q)).train ().eval_ () println (ts.report) new Plot (null, y, ts.predictAll (), s"Plot of y, SARIMA($p, 0, $q) vs. t", true) val arma_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA($p, 0, $q) model = $arma_f") } // SARIMATest2 object //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest3` object is used to test the `SARIMA` class. * Forecasting lake levels. * @see cran.r-project.org/web/packages/fpp/fpp.pdf * > runMain scalation.analytics.forecaster.SARIMATest3 */ object SARIMATest3 extends App { import ForecasterVec.y // val (pp, dd, qq) = (1, 1, 1) // seasonal hyper-parameters // val period = 4 // seasonal period // val d = 0 // apply d-th order differencing - no differencing val d = 1 // apply d-th order differencing - first differences for (h <- 1 to 2) { // forecasting horizon for (p <- 0 to 4) { // auto-regressive hyper-parameter settings for (q <- 0 to 4 if p+q > 0) { // moving-average hyper-parameter settings banner (s"Test3: SARIMA ($p, $d, $q) with h = $h") val mod = new SARIMA (y, d) // create an SARIMA model mod.setPQ (VectorI (p, q)).train ().eval_ () // set p and q, train and evaluate the SARIMA model println (mod.report) val yp = mod.predictAll () // predicted/fitted values val yfa = mod.forecastAll (h) val yfb = mod.forecastAll2 (h) assert (yfa == yfb) val yf = yfa.col(h) // forecasted values - h steps ahead val yf2 = yfb.col(h) // forecasted 2 values - h steps ahead new Plot (null, y, yp, s"Plot of y & yp, predicted SARIMA ($p, $d, $q) vs. t", true) new Plot (null, y, yf, s"Plot of y & yf, forecasted (h = $h) SARIMA ($p, $d, $q) vs. t", true) new Plot (null, y, yf2, s"Plot of y & yf2, forecasted2 (h = $h) SARIMA ($p, $d, $q) vs. t", true) if (h == 1) ForecasterVec.differ (yp, yf, allow = true) ForecasterVec.differ (yf, yf2, allow = true) val skip = max (p, q) // skip the cap start-up banner (s"SARIMATest3: QoF (@h = $h) for yf = mod.forecastAll") println (s"rSq (yf) for h = $h is ${rSqF (y, yf)}") println (s"rSq (yf) for h = $h is ${rSqF (y, yf, max (p, q))}, skip = $skip") println (s"sMAPE (yf) for h = $h is ${smapeF (y, yf)}") println (s"sMAPE (yf) for h = $h is ${smapeF (y, yf, max (p, q))}, skip = $skip") banner (s"SARIMATest3: QoF (@h = $h) for yf2 = mod.forecastAll2") println (s"rSq (yf2) for h = $h is ${rSqF (y, yf2)}") println (s"rSq (yf2) for h = $h is ${rSqF (y, yf2, max (p, q))}, skip = $skip") println (s"sMAPE (yf2) for h = $h is ${smapeF (y, yf2)}") println (s"sMAPE (yf2) for h = $h is ${smapeF (y, yf2, max (p, q))}, skip = $skip") } // for } // for } // for } // SARIMATest3 object //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest4` object is used to test the `SARIMA` class. * > runMain scalation.analytics.forecaster.SARIMATest4 */ object SARIMATest4 extends App { val path = BASE_DIR + "travelTime.csv" val data = MatrixD (path) val t = data.col(0) val y = data.col(1) val p = 1 val q = 1 val d = 1 // apply 1st order differencing val steps = 1 // number of steps for the forecasts val ts = new SARIMA (y, d) // time series data: y vs. t println (s"y = $y") // Build SARIMA(p, d, q) models ts.setPQ (VectorI (p, 0)).train ().eval_ () println (ts.report) new Plot (t, y, ts.predictAll (), s"Plot of y, SARIMA($p, $d, 0) vs. t", true) val ar_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA($p, $d, $q) model = $ar_f") ts.setPQ (VectorI (0, q)).train ().eval_ () println (ts.report) new Plot (t, y, ts.predictAll (), s"Plot of y, SARIMA(0, $d, $q) vs. t", true) val ma_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA($p, $d, $q) model = $ma_f") ts.setPQ (VectorI (p, q)).train ().eval_ () println (ts.report) new Plot (t, y, ts.predictAll (), s"Plot of y, SARIMA($p, $d, $q) vs. t", true) val arma_f = ts.forecast (h = steps) println (s"$steps-step ahead forecasts using SARIMA($p, $d, $q) model = $arma_f") } // SARIMATest4 object //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest5` object is used to test the `SARIMA` class. * > runMain scalation.analytics.forecaster.SARIMATest5 */ object SARIMATest5 extends App { val m = 50 val d = 0 // levels of differencing val steps = 2 // number of steps for the forecasts val sig2 = 10000.0 val noise = Normal (0.0, sig2) val y = VectorD (for (i <- 0 until m) yield 40 * (i-1) - (i-2) * (i-2) + noise.gen) val mod = new SARIMA (y, d) // time series model SARIMA banner (s"Build SARIMA(1, $d, 0) model") mod.setPQ (VectorI (1, 0)).train ().eval_ () // train for SARIMA(1, d, 0) model testInSamp (y, 2, mod, 1) for (p <- 1 to 3) { banner (s"Build SARIMA($p, $d, $p) model") mod.setPQ (VectorI (p, p)).train ().eval_ () // retrain for SARIMA(p, d, p) model testInSamp (y, p+p, mod, 1) banner ("Make Forecasts") val yf = mod.forecast (steps) println (s"$steps-step ahead forecasts using SARIMA($p, $d, $p) model = $yf") } // for } // SARIMATest5 object //:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::: /** The `SARIMATest6` object is used to test the `SARIMA` class on real data: * Forecasting lake levels. Select the best number of lags. * @see cran.r-project.org/web/packages/fpp/fpp.pdf * > runMain scalation.analytics.forecaster.SARIMATest6 */ object SARIMATest6 extends App { import ForecasterVec.y val d = 0 // level of differencing val mod = new SARIMA (y, d) // create model for time series data mod.setPQ (VectorI (1, 1)).train ().eval_ () // train the model and evaluate val res = mod.forwardSel () println (s"forwardSel: $res") for ((sp, sq) <- Array ((1, 0), (2, 1), (1, 1), (1, 2), (0, 1))) { val res = mod.forwardSel2 (VectorI (sp, sq)) println (s"forwardSel2: $res") } // for } // SARIMATest6 object