where \(\delta_i(B)=(\delta_{0i}+\delta_{1i}B^1+...+\delta_{si}B^{s_i})\) and \(\omega_i(B)=(1-\omega_{1i}B^1-...-\omega_{ri}B^{r_i})\) and \(s_i\), \(r_{i\ }\)and \(d_i\) are the orders of the transfer function for each input time series \(X_i\left(t\right)\). The numerator polynomial (order \(s_i\)) controls the lagged effect of the input variable \(X_i\left(t\right)\) on the output time series; and the denominator (order \(r_i\) ) controls the lagged effect of variable \(Y(t)\) induced by \(X_i(t)\). The backshift operator order \(d_i\) and the orders \(s_i\) and \(r_i\) are estimated by calculating the cross-correlation function between the pre-whitened input time series \(X_i(t)\) (using an ARIMA model, such that the filter \(\Phi(B)/\Theta(B)\) applied to \(X_i(t)\) produces a white noise \(v\left(t\right)\)) and the filtered output series \(\tilde{Y}(t)\ \) by applying the same ARIMA filter to the input time series. Once the orders \(s_i\), \(r_i\) and \(d_i\) are identified, the corresponding linear transfer function is fitted using the R software. In our case \(i=1,2\), since two input time series are used as covariates.
Model Validation Methods
To ensure the optimal runs over the hold back prediction (testing validation), model parameterization was achieved by minimizing together the Root Mean Squared Error (RMSE) and the Mean Absolute Scaled Error (MASE), and maximizing the correlation coefficient (R):
\(\)