Table 7 summarises the MSPE for the influenza season. The best result for each set of covariates is shaded in light gray and the overall winner for each horizon is in bold font. In the models with the predictor Raten,s(w) (rows (a), (e), (f) and (g)) the gain, in terms of MSPE, of the functional GLS models (GLS–AR(1), iGLS–AR(1), iGLS–AR(2) and iGLS–AR(p)) is relatively small with respect to functional LM models because the Raten,s(w) partly accounts for the temporal dependence. Furthermore, in some sense, the inclusion of the predictor Raten,s(w) in the model is akin to the estimation of the dependence structure. The models without influenza rate (rows (b), (c), (d) and (h)) begin with a worse result in the LM setting, but their results become competitive (or even become the best ones) with the inclusion of the serial dependence. The difference between the GLS or iGLS setting is that the latter allows more flexibility, not only defining a different dependence structure in each county, but also in the estimation of that dependence. This is particularly useful when the forecast horizon increases. The GLS setting must fix the order of the AR in advance and, when the number of regions is high, it is a tough assumption to consider the order of the serial dependence model fixed for all of them. For n + 1 the best models are (b) and (c) with GLS–AR(1) and iGLS–AR(1) specifications, using the curve of temperature of last 14 days as the predictor and a simple AR(1) structure for the adjustment of the residuals. The best autoregressive model estimated by the iGLS–AR(p) model has been, in most cases, of order 1. For n + 2, in some regions, an AR(1) or AR(2) model may be insufficient; the best result is achieved with the iGLS–AR(p) procedure, which presents greater flexibility in estimating the different p order for each county.