Heating and Cooling Degree-Days Climate Change Projections for Portugal

Climate change is expected to influence cooling and heating energy demand of residential buildings and affect overall thermal comfort. Towards this end, the heating degree-day (HDD), the cooling degree-day (CDD) and the HDD+CDD were computed from an ensemble of 7 high-resolution bias-corrected simulations attained from EURO-CORDEX under RCP4.5 and RCP8.5. These three indicators were analyzed for 1971−2000 (from E-OBS) and 2011−2040 and 2041−2070, under both RCPs. Results show that the overall spatial distribution of HDD trends for the 3 time-periods points out an increase of energy demand to heat internal environments in Portugal's northern-eastern regions, most significant under RCP8.5. It is projected an increase of CDD values for both scenarios; however, statistically significant linear trends were only found for 2041−2070 under RCP4.5. The need for cooling is almost negligible for the remaining periods, though linear trend values are still considerably higher for 2041−2070 under RCP8.5. By the end of 2070, higher amplitudes for all indicators are depicted for southern Algarve and Alentejo regions, mainly under RCP8.5. For 2041−2070 the Centre and Alentejo (North and Centre) regions present major positive differences for HDD(CDD) under RCP4.5(RCP8.5), within the 5 NUTS II regions predicting higher heating(cooling) requirements for some locations.

The projected rise in temperatures [8] is expected to pose greater risks to urban areas.
The extent of the risk depends on human vulnerability and adaptation effectiveness, namely in the construction sector. Indoor environment conditions contribute greatly to human wellbeing, as most people spend around 90% of their time indoors, mainly at home or in the workplace [9]. The fluctuations in outdoor air temperatures [10] will have not only a substantial impact on human comfort, but also on building energy use [11] mainly in the existing residential buildings. Therefore estimate air temperature fluctuation projections have relevant implications for estimating its future impacts on residential heating and cooling related energy demand.
Several studies used multiple methods to estimate future residential heating and cooling energy demand in buildings. While some authors choose simple approaches such as using current climate, discarding climate variability [12] or choosing a warm past year to represent a warming climate [13], others opt to use climate models using several datasets, namely global climate simulation models (GCMs) [11,[14][15][16][17][18][19]. The most common methods used to determine residential demand in the future use parametric energy balance and degree-day methods. The degree-day method is a simple and widely used approach to relate outdoor temperature with the heating/cooling energy requirements.
In this study, we have employed the degree-day method following the methodology used by [20] and, later on, by [11]. This methodology defines a base temperature (Tb) for the heating and cooling season and allows the computation of the respective outdoor air temperature deviations from maximum and minimum temperatures. The base temperature is a point at which internal gains equals the heat loss, acting as a threshold below (or above) which heating (or cooling) appliances are needed or not to operate to maintain indoor thermal comfort. Under the Portuguese Regulation on the Energy Performance of Residential Buildings (REPRS) [21,22], these temperatures are 18C related to the degrees-day of heating (HDD) and 25C to the degrees-day of cooling (CDD).
Further details will be provided in the Materials and Methods section; however, it is worth mention that the REPRS is in line with the European Directive 2010/31/EC [23], which aims at reducing the greenhouse gas emissions by 20% by 2020 and in 80% until 2050, in relation to the 1990 emissions levels. Therefore, this objective includes the adoption of standard methodologies for calculating energy consumption, quality requirements for new and existing building envelopes, periodic inspection of boilers and air conditioning central systems, as well as building energy certification.
Three key energy performance indicators were computed in this work: the HDD, the CDD and the global indicator HDD+CDD, obtained from an ensemble mean of seven biased corrected regional climate models (RCMs) for mainland Portugal. kriging (OK), and ordinary cokriging (OCK) are the most frequently used techniques in environmental studies for spatial interpolation of data [24][25][26][27][28][29]. Several geostatistical techniques were performed in this study to attain the most accurate spatial representation of the different indicators.
This study's main goal is to analyse the impacts of climate change on heating or cooling related energy demand for residential buildings thermal comfort by computing

Materials and Methods
The overall methodology framework can be depicted in Figure 2 but will be detailed in the following subsections.  [39,40,41].
In this study, seven RCMs were retrieved from EURO-CORDEX (Table 1). in Figure 3. The E-OBS datasets and the respective GCMs that have a coarser spatial resolution on 0.11° regular grid, overlap thus allowing this bias correction. In this study, we used the quantile-quantile bias correction. This method assumes that the distribution function of a variable may change in the future. However, this methodology allows the correction of the complete distribution, tails included. Further details regarding this methodology can be found in [5]. Bias correction was applied to daily mean, minimum  (Table 1) was used to compute the HDD and the CDD.

Heating (HDD) and cooling (CDD) degree-day
Energy consumption linked to the thermal comfort of buildings is related to the HDD and the CDD. The HDD translates the amount of energy needed (i.e., to a building with a heating system) on a given day or period to heat the indoor environment in a climate considered cold to a specific base temperature (18°C). The CDD reflects the amount of energy required (i.e., for a building with a cooling system) on a given day or period to cool the indoor environment in a climate considered warm to a specific base temperature (25°C). The theoretical formulation for calculating the HDD and CDD can be carried out in several ways. Calculations can be performed using monthly or annual data or with more sophisticated models. Although the base temperature values may differ, depending on the country under analysis, in this work, the daily values for the HDD should be determined using a base temperature (Tb in Table 2) of 18°C, while the daily CDD values using a base temperature (Tb in Table 3) of 25°C [21,22]. Daily HDD and CDD values are then calculated following the cases in Tables 2 and 3 [11,20], respectively, in which Ta is the mean calculated from the Tx and Tn temperature values.  [11,20]. the annual values for HDD should be calculated as the cumulative sum of the daily HDD values for the 'cold season' in which there is now a need to 'heat up' the internal environment of the buildings. This heating station is considered to start on the first 10day mean after 1 October when the average daily temperature is below 15°C and ends in the last 10-day mean before 31 May in which that temperature is still below 15°C. [42,43] proposed a combined degree-day index by summing HDD and CDD The spatial representation of these indicators will be presented after careful consideration of the best interpolation techniques that will be explored herein.

Geostatistical techniques
Geostatistical methods have been shown superior to the conventional and deterministic methods for spatial interpolation of rainfall [24]. Kriging and cokriging are two spatial interpolation methods that have been widely used to create spatially continuous climate-related data [25]. They estimate the value of a variable or indicator of interest at an unmonitored location based on the values at neighbouring monitored locations by fitting a semi-variogram model, which is a function of spatial distance. Simple kriging (SK) and Ordinary kriging (OK) differ by the methods used to model the means of primary and secondary variables. SK assumes that local means are relatively constant and equal to the population mean, which is well known. The population mean is used as a factor in each local estimate, along with the samples in the local neighbourhood [44].
Estimated primary and secondary local means could differ from the means calculated on the whole dataset. Consequently, OK did not require knowledge of the primary and secondary local means [44].
Cokriging allows additional predictor variables that exhibit inter-correlations with the variable of interest, possibly producing better prediction performance than the kriging method. This can help to minimise the error variance of the estimation [45]. The standard form of cokriging is the OCK method. This usually reduces the prediction error variance and specifically outperforms the kriging method if the secondary variable, the digital terrain model (DTM), is highly correlated (correlation coefficient higher than 0.75 ) with the primary variable and many more points are known [26].
Kriging was used to interpolate temperature and precipitation in the Mediterranean by [46]. [27] used OK to interpolate monthly temperature anomalies but preferred IDW for precipitation. In [47], IDW was chosen since it captures well local variations and captures exact values at co-located grid points for several climate variables. [29] compared IDW, OK, and OCK to predict air temperature at unmeasured Turkey sites. The OCK with elevation as an auxiliary variable proved to be the best technique to predict temperature against the criteria of model efficiency and relative root mean squared error (RMSE). Covariables derived from DTM are widely used to adjust topographic conditions [28] in interpolation techniques. However, the best technique's choice must be carefully evaluated since the temperature is not solely determined by elevation and land cover but also by atmospheric circulation patterns in the northern hemisphere [48]. Moreover, it has been reported that in some areas, precipitation was not related to elevation [49].
[50] studied the spatial pattern of CDD on a typical normal and extremely hot summer day using OCK geospatial mapping technique. Results revealed reasonable predictability of city-wide CDD with the OCK method, which uses two co-variables: "elevation of the weather station" and "building volume density within the 1,000 m radius neighbouring area".
[20] used OK to project future HDD, CDD and HDD+CDD in the USA.
In this investigation, the ArcGIS Geostatistical Analyst was used, and three techniques were tested: IDW, OK and OCK. The input datasets, in this case, HDD, CDD and HDD+CDD were evaluated regarding 1) data distribution, 2) global trends, and 3) directional influences.
First, all datasets were tested regarding their normality (frequency histograms for the attributes) being subject to a transformation when skewed since the normal distribution datasets generate better results. Trend analyses identify the presence or absence of trends in the input dataset and identify which polynomials order best fit the trend. Local variation can be added to the surface by modelling the trend using one of the smooth functions, removing it from the data and allowing the subsequent analysis. Therefore, this evaluation was performed for all variables. Lastly, since a directional influence will affect the semi-variogram and the fit of the model, the -semivariogram model's anisotropy must also be evaluated. The directional influence can be statistically quantified and accounted for when making the map.
Following the methodology previously presented, the HDD and HDD+CDD datasets histograms showed that their distributions were not normal, so a logarithmic transformation was performed; conversely, since CDD showed a normal distribution, no transformation was done. Regarding the trend analysis ( Figure 4), an upward trend in the West-East direction was detected for all HDD datasets. Due to mainland Portugal location, for CDD and HDD+CDD, this trend is expected since the energy requirements for heating or cooling increase from west to east (oceanic influence). The HDD and CDD datasets trends in the North-South direction are also predictable since the heating(cooling) requirements decrease(increase) towards the south.
The HDD+CDD dataset trends are similar to the HDD since the HDD values are relatively higher than the CDD, strongly influencing the sum. Consequently, these results substantiate the need to test the semi-variogram models with trend-removing functions.
A first-order trend removal function was thus used since the trends proved to be almost linear.
In this study, IDW and eleven semi-variograms were tested for both OK and OCK:  [24]. Results for the cross-validation statistics can be observed in Table 4. Statistically significant trends (at a 5% significance level) were also assessed by using the rank-based non-parametric Spearman's rho (SR) statistical test [53,54]. This nonparametric test can be used to detect monotonic trends in time series and is widely used in hydro-meteorological studies. The magnitude of the slope of the trend was estimated using Theil and Sen's approach [55,56]. The slope was estimated by Lastly, the areolar mean (for mainland Portugal) for each indicator was computed, and statistically significant linear trends were obtained for 30-years-time periods between 1971 and 2070 under RCP4.5 and RCP8.5. Only the statistically significant linear regression models will be presented for each period (p-value < 5%, e.g., at a 5% significance level) with the associated indicator time-series.      Given the results previously attained (Figures 8 and 9

Case Study: NUTS II
An analysis for a case study within the NUTS II region (Figure 1 (Figures 1 and 11). Overall, results show that regions with higher projected cooling or heating demands present higher increases under both RCPs until 2070. However, higher needs are predicted for the cooling needs since for the heating requirements, like previously stated, non-significant increases were detected. Therefore, it can be concluded that for degree-day values, future spatial distribution for 2011−2040 no significant changes are projected on a national scale, although on a regional scale, that might not be the case (Figure 11).  (Figures 11b and 12). global climate projections [34] provide these regional simulations in the RCPs [35,36]. In this case, calculations and subsequence analysis were made under RCP4.5 and RCP8.5

Preprints
scenarios. An observational dataset of corresponding temperatures E-OBS was used to bias correct the simulations. In this study, we used the quantile-quantile bias correction, which assumes that the distribution function may change in the future.
From a seven-member bias-corrected ensemble of maximum and minimum daily temperatures, the HDD, CDD and HDD+CDD indicators were computed. Having in mind that, although the base temperature values may differ, depending on the country under analysis, in this work, the daily values for the HDD were determined using a base temperature of 18°C, while the daily CDD values using a base temperature of 25°C following the Portuguese legislation [21]. Daily HDD and CDD values were then calculated following the [11,20] methodology. As a result of these methodological changes due to the specifications of the Portuguese Law, the magnitude of the indicators and trends attained in this work and other studies that encompasses Portugal within Europe cannot be directly compared, that is the case of [11,57]. Proposed by [42,43] Geostatistical analysis of the three indicators was performed following the methodology previously presented. As such, for HDD and HDD+CDD datasets histograms showed that their distributions did not follow the normal distribution, therefore a logarithmic transformation was performed; conversely, since CDD followed a normal distribution, no transformation was done. The trend analysis (Figure 4) showed an upward trend in the West-East direction also detected for all HDD datasets. Due to mainland Portugal location, for CDD and HDD+CDD, this trend is expected since the energy requirements for heating or cooling increase from west to east due to the oceanity influence on climate. The HDD and CDD datasets trends in the North-South direction are also predictable since the heating(cooling) requirements decrease(increase) towards the south. The HDD+CDD dataset trends were similar to the HDD since the HDD values are relatively higher than the CDD, consequently, strongly influence the sum. Subsequently, these results substantiate the need to test the semi-variogram models with trendremoving functions. A first-order trend removal function was thus used since the trends proved to be almost linear.
In this study, IDW and eleven semi-variograms were tested for both OK and OCK: The statistically significant anomalies were assessed by the Mann-Whitney-Wilcoxon test (MWW) at a 5% significance level [51,52]. Statistically significant trends (at a 5% significance level) were also assessed by using the rank-based non-parametric Spearman's rho (SR) statistical test [53,54]. Lastly, the areolar mean (for mainland Portugal) for each indicator was computed, and statistically significant linear trends were obtained for 30years-time periods between 1971 and 2070 under RCP4.5 and RCP8.5 scenarios. Only the statistically significant linear regression models were presented for each period (p-value < 5%, e.g., at a 5% significance level) with the associated indicator time-series.
The main outcomes of this study will be summarised herein:   The Portuguese Regulation on the Energy Performance of Residential Buildings [21], as aforementioned, is in line with the European Directive [23], which aims at reducing the greenhouse gas emissions by 20% by 2020 and in 80% until 2050, in relation to the 1990 emissions levels. This study allowed to conclude that major differences in heating and cooling energy demand can be expected for mainland Portugal under both RCPs and until 2070. The predicted regional differences in residential buildings stock heating and cooling