A 1-km global dataset of historical (1979-2013) and future (2020-2100) Köppen-Geiger climate classification and bioclimatic variables

The Köppen-Geiger classification scheme provides an effective and ecologically meaningful way to characterize climatic conditions and has been widely applied in climate change studies. Significant changes in the Köppen climates have been observed and projected in the recent two centuries. Current accuracy, temporal coverage, spatial and temporal resolution of historical and future climate classification maps cannot sufficiently fulfil the current needs of climate change research. 10 Comprehensive assessment of climate change impacts requires a more accurate depiction of fine-grained climatic conditions and continuous long-term time coverage. Here, we present a series of improved 1-km Köppen-Geiger climate classification maps for six historical periods in 1979-2013 and four future periods in 2020-2099 under RCP2.6, 4.5, 6.0, and 8.5. The historical maps are derived from multiple downscaled observational datasets and the future maps are derived from an ensemble of bias-corrected downscaled CMIP5 projections. In addition to climate classification maps, we calculate 12 bioclimatic 15 variables at 1-km resolution, providing detailed descriptions of annual averages, seasonality, and stressful conditions of climates. The new maps offer higher classification accuracy than existing climate map products and demonstrate the ability to capture recent and future projected changes in spatial distributions of climate zones. On regional and continental scales, the new maps show accurate depictions of topographic features and correspond closely with vegetation distributions. We also provide a heuristic application example to detect long-term global-scale area changes of climate zones. This high-resolution 20 dataset of the Köppen-Geiger climate classification and bioclimatic variables can be used in conjunction with species distribution models to promote biodiversity conservation and to analyze and identify recent and future interannual or interdecadal changes in climate zones on a global or regional scale. The dataset referred to as KGClim, is publicly available at http://doi.org/10.5281/zenodo.5347837 for historical climate and http://doi.org/10.5281/zenodo.4542076 for future climate.

inconsistent temporal spans of climatology datasets to the period 1980-2016, by adding interpolated temperature change offsets or multiplying precipitation factors, which may lead to biased coverage of the historical period. Current classification accuracy, temporal coverage, spatial and temporal resolution of historical and future climate classification maps cannot sufficiently fulfil the current needs of climate change research. Significant changes in the Köppen climates have been observed and projected in 65 the recent two centuries (Belda et al., 2014;Chan & Wu, 2015;Chen & Chen, 2013;Rohli, Andrew, Reynolds, Shaw, & Vázquez, 2015;Yoo & Rohli, 2016). Previous studies found that large-scale shifts in climate zones have been observed over more than 5% of the total land area since the 1980s, and approximately 20.0% of the total land area is projected to experience climate zone changes under RCP8.5 by 2100 (Cui, Liang, & Wang, 2021). Detection of recent and future changes in climate zones with the application of the Köppen climate maps needs more accurate depiction of fine-grained climatic conditions, 70 continuous and longer temporal coverage.
This creates the urgent need for global maps of the Köppen climate classification with increased accuracy, finer spatial and temporal resolutions. Currently available global observational datasets of temperature and precipitation collected during the recent centuries, and the global climate simulations under alternative future climate scenarios have offered the possibility to create a comprehensive dataset for past and future climates. In this study, we presented an improved long-term  climate classification map series for 1) six historical 30-yr periods of the observational record (1979-2008,1980-2009, 1981-2010, 1982-2011, 1983-2012, 1984-2013) and four future 30-yr periods (2020-2049, 2040-2069, 2060-2089, 2070-2099) under four RCPs (RCP2.6, 4.5, 6.0 and 8.5). To improve the classification accuracy and achieve a resolution as fine as 1-km (30 arcsecond), we combined multiple datasets, including WorldClim V2 (Fick & Hijmans, 2017), CHELSA V1.2 (Booth et al., 2014), CRU TS v4.03 (New, Hulme, & Jones, 2000), UDEL (Willmott & Matsuura, 2001), GPCC datasets (Beck, Grieser, & 80 Rudolf, 2005) and bias-corrected downscaled Coupled Model Intercomparison Project Phase 5 (CMIP5) model simulations (Navarro-Racines, Tarapues, Thornton, Jarvis, & Ramirez-Villegas, 2020) (Table 1). We used the WorldClim Historical Climate Data V2 (Fick & Hijmans, 2017) to downscale the 0.5° climatology datasets including CRU, UDEL and GPCC, and derive high resolution climate data for the historical periods. To determine the final climate class, we used the climate class with the highest agreement level from an ensemble of climate maps derived from different combinations of surface air 85 temperature and precipitation products, as implemented in Beck et al. (2018). In addition to the Köppen-Geiger climate maps, we also calculated 12 bioclimatic variables at the same 1-km resolution using these climate datasets for the same historical and future periods. This dataset can be used to in conjunction with SDMs to promote biodiversity conservation, or to map plant functional type distributions for earth system model simulations, or to analyse and identify recent and future changes in climate zones on a global or regional scale. 90 To validate the Köppen-Geiger climate classification maps, we used the station observations from Global Historical Climatology Network-Daily (GHCN-D) (Menne, Durre, Vose, Gleason, & Houston, 2012), and Global Summary Of the Day (GSOD) (National Climatic Data Center, NESDIS, NOAA, & U.S. Department of Commerce, 2015) database. At the regional and continental scale, we compared our Köppen-Geiger climate classification maps with previous map products, associated maps of forest cover, and elevation distribution, for 1) regions with large spatial gradients in climates, including central and 95 eastern Africa, Europe, North America, and 2) regions with sharp elevation gradients, including Tibetan Plateau, central Rocky Mountains, central Andes. Further, we conducted sensitivity analysis with respect to classification temporal scale, dataset input, and data integration methods. We also provided a heuristic example which used climate classification map series to detect the long-term area changes of climate zones, showing how the Köppen-Geiger climate classification map series can be applied in climate change studies. 100  (New et al., 2000), University of Delaware Precipitation and Air Temperature (UDEL) (Willmott & Matsuura, 2001) and Global Precipitation Climatology Centre (GPCC) (Beck et al., 2005) datasets. To decide the datasets to use, we conducted a sensitivity analysis on the input climatology datasets and utilized monthly air temperature datasets from CRU, UDEL, GHCN_CAMS Gridded 2m Temperature (Fan & Dool, 2008) and monthly precipitation datasets from GPCC, UDEL, NOAA's PRECipitation REConstruction over Land (PREC/L) (Chen, Xie, Janowiak, & Arkin, 2002). Evaluation results 110 indicated that incorporating only CRU, UDEL temperature datasets and UDEL, GPCC precipitation datasets and excluding GHCN_CAMS and PREC/L datasets led to higher accuracy in the classification results. Therefore, we chose CRU, UDEL, and GPCC datasets as the classification system input to boost the final accuracy.

Datasets
To explicitly correct topographic effect, we used 1-km CHELSA V1.2 and WorldClim V2 datasets in addition to the 0.5° resolution datasets. The CHELSEA dataset statistically downscaled temperature data from the ERA-Interim climatic reanalysis. 115 For precipitation data, it incorporated multiple orographic predictors and performed bias correction (Karger et al., 2017). With major topo-climatic drivers considered, the CHELSA dataset demonstrated good performance in ecological science studies.
CHELSA data exhibited comparable accuracy for temperatures and better predicted precipitation patterns based on the validation results (Karger et al., 2017).
We produced the future Köppen classification map series using the CCAFS climate statistically bias-corrected and downscaled 120 CMIP5 projections (Navarro-Racines et al., 2020). The CCAFS presented a global database of future climates developed by a climate model bias correction method based on the CMIP5 GCM simulations (Taylor, Stouffer, & Meehl, 2012) archive, coordinated by the World Climate Research Programme in support of the IPCC Fifth Assessment Report (AR5) (Hartmann et al., 2013). The total is 35 GCMs, and all RCPs, RCP 2.6, 4.5, 6.0 and 8.5 (Table S1). Projections are available at varied coarse scales (70-400km). To achieve high-resolution (1km) climate representations, downscaling method has been applied with the 125 use of the WorldClim data (Fick & Hijmans, 2017) . Technical evaluation showed that the bias-correction method that CCAFS data applied reduced climate model bias by 50-70%, which could potentially address the bias issue in model simulations for the threshold-based Köppen classification scheme (Navarro-Racines et al., 2020).
The Köppen climate classification scheme was first introduced by Wladimir Köppen (1884). It is one of the earliest quantitative classification systems of Earth's climates. Its modification, the Köppen-Geiger classification (KGC) was first published in 1936 (Köppen, 1936), developed by Wladimir Köppen and Rudolf Geiger. KGC identifies climates based on their effects on plant growth from the aspects of warmth and aridity, and classifies climate into five main climate classes and 30 subtypes (Rubel & Kottek, 2011). The five main climate zones distinguish between plants of the tropical climate zone (A), the arid 135 climate zone (B), the temperate climate zone (C), the boreal climate zone (D) and the polar climate zone (E), referring to the five major climate zones (Sanderson, 1999). All these main climate zones are thermal zones except the arid (B) climate zone, which is defined based on precipitation threshold. This research followed the Köppen-Geiger climate classification as described in Kottek et al. (2006), and Rubel & Kottek (2010). This latest version of the KGC scheme was first presented by Geiger (1961) (Table 2). Several existing  climate map products, including Peel et al. (2007), Kriticos et al. (2012), and Beck et al. (2018) applied the KGC scheme modified following Russell (1931). Russell (1931) adjusted the definition of the boundary of temperate (C) and boreal (D) climate zones using the coldest monthly temperature > 0 °C instead of >-3 °C. This threshold was proposed because the 0°C line fits the distribution of the topographical features and vegetation in western United States, where at that time meteorological stations were sparsely distributed (Jones, 1932). However, the application of 0°C boundary to the global climates has not been 145 validated. Therefore, this research didn't utilize the Russell's modification (1931) and followed the latest version KGC proposed by Geiger (1961).  and (b) delta downscaling method with January temperature from CRU dataset. Baseline  and present-day climate data (e.g. 1979-2008) are from CRU, UDEL, or GPCC datasets, which have a coarse spatial resolution of 0.5 o . Precipitation anomaly is change factor of monthly precipitation from baseline to present-day climates. Temperature delta is change in monthly air temperature from baseline to present-day climates. WorldClim  climate data is adjusted by multiplying 30 arc-second interpolated anomaly (for precipitation) or adding 30 arc-second interpolated delta (for temperature) to generate the downscaled climate surfaces with 30 arc-second resolution.

Statistical downscaling
Due to limited number of available observational datasets with high resolution and long-term continuous temporal coverage, the research implemented the delta method by applying a delta change or change factor (Hay, Wilby, & Leavesley, 2000;Wilby & Wigley, 1997) onto the WorldClim historical observations (Fick & Hijmans, 2017) to achieve 30-yr average climatology data with a 1-km resolution based on the CRU, UDEL and GPCC datasets. The delta method is a statistical 160 downscaling method that assumes that the relationship between climatic variables remain relatively constant at local scale (Wilby & Wigley, 1997). We applied delta method to downscale the long-term (30-yr) mean climates using coarse-resolution monthly climatology datasets. The delta changes or change factors are calculated as the differences between the 30-yr longterm means of temperature or precipitation of baseline  and present-day climates. The delta method comprises the following four steps: 1) calculate 30-yr averages for baseline  and present day of monthly temperature and 165 precipitation; 2) calculate anomaly for precipitation and delta for temperature; 3) apply thin-plate splines interpolation (TPS) to create 1km surface of precipitation anomaly and temperature delta; 4) multiply anomaly or add delta to historical climates based on WorldClim dataset (Fig. 1).
First, using monthly time series from CRU, UDEL and GPCC datasets, we calculated 30-yr means as a baseline , for each climatology dataset and each variable. We used 1970-2000 as baseline period, for consistency with WorldClim 170 Historical Climate Data V2. Next, we calculated 30-yr means for each month and each 30-yr present-day period in 1979-2013 We then calculated anomalies as proportional differences between present-day and baseline in total precipitation and delta as difference in temperature. To derive 30 arc-second (1-km) anomaly or delta surfaces, we applied thin-plate splines (TPS) interpolation (Craven & Wahba, 1978;Franke, 1982;Schempp, Zeller, & Duchon, 1977) on precipitation anomaly and temperature delta. TPS has been widely used in climate science (Hijmans, Cameron, Parra, Jones, & Jarvis, 2005;Navarro-175 Racines et al., 2020) as it produced a smooth and continuous surface, which is infinitely differentiable. Last, we multiplied the change factor or added the delta to the WorldClim  data to get downscaled present-day monthly climate data.
Our future Köppen-Geiger map series are based on an ensemble of maps derived from the CCAFS bias-corrected and downscaled climate projections, which include 35 CMIP5 GCMs, and 4 RCPs (Navarro-Racines et al., 2020). Large misclassifications exist within the GCMs as detected in previous assessment of large areas ranging between 20-50% of the 180 total land area (Cui, Liang, & Wang, 2021). Deficiencies in model physics are also more likely to contribute to uncertainties in the maps than grid size or reference dataset limitations (Tapiador, Moreno, & Navarro, 2019). Multi-model mean and deltachange method can mitigate the bias effects from the threshold-based classification scheme and have been utilized to simulate better results of climate classification (Hanf, Körper, Spangehl, & Cubasch, 2012). Therefore, we chose the CCAFS biascorrected and downscaled CMIP5 projections (Navarro-Racines et al., 2020) to reduce the amplified errors due to uncertainty 185 of climate projections. Navarro-Racines et al. (2020) interpolated anomalies of original GCM outputs using thin plate spline spatial interpolation to achieve a baseline climate with a 1km surface. Then they applied delta method to the interpolated baseline climates to correct the model biases (Hay et al., 2000;Ho, Stephenson, Collins, Ferro, & Brown, 2012). The historical Köppen-Geiger climate classification map series was generated using the highest confidence class from an ensemble of maps using all combinations of surface air temperature and precipitation products (Fig. 2), as described in Beck et al. (2018). The highest confidence was given to the most common climate class for each grid cell. The final historical climate map series were derived using the climate class with the highest level of confidence in an ensemble of 3 × 3 = 9 classification 195 maps based on combinations of the 3 precipitation datasets (CRU, UDEL, and CHELSA) and 3 surface air temperature datasets (GPCC, UDEL, and CHELSA). To further test the sensitivity of the method using the climate with the highest level of agreement, we incorporated another data integration method using the mean of multiple datasets. We quantified the degree of confidence placed in the Köppen-Geiger climate map series using the degree of confidence at the grid cell level calculated by dividing the occurrence frequency of the climate class with the highest level of agreement by the ensemble size. The calculated 200 confidence level can be viewed as the agreement degree in classification resulted derived from different climatology datasets.

Data Integration
The future Köppen-Geiger climate classification map series under 4 RCPs, were derived based on the most common climate class from an ensemble of future climate maps. We generated a future Köppen-Geiger climate classification map for each climate model projection, using the CCAFS bias-corrected and downscaled CMIP5 GCM dataset. For example, the future Köppen-Geiger climate classification map series under RCP8.5 was derived from an ensemble of 30 maps based on 30 CMIP5 205 models. The level of confidence was estimated using the ratio between the frequency of the climate class with the highest level of agreement in the future map results, and the ensemble size.

Validation
We validated the historical climate maps using the station observations from Global Historical Climatology Network-Daily (GHCN-D) (Menne et al., 2012) and Global Summary Of the Day (GSOD) database (National Climatic Data Center et al., 210 2015) as reference data. GHCN-D dataset provides daily climate data over global land areas and contains records from over 80,000 weather stations worldwide, about one third of which have both temperature and precipitation data available (Menne et al., 2012). GSOD dataset includes global daily summary data over 9,000 stations, of which the historical data from 1973 being the most complete (National Climatic Data Center et al., 2015). For each station, time series of monthly temperature and precipitation were calculated from the daily observations with months with <15 daily values discarded. Then if ≥6 months 215 are present, monthly climatology were generated subsequently by averaging the monthly means for the given 30-yr period.
We removed duplicate stations in the two datasets and discarded stations with gap years or missing data in the given 30 years.
For each station and each 30-yr period, we applied the Köppen-Geiger climate classification, and then evaluated overall classification performance for each climate map using total accuracy, which is defined as the percentage of correct classes, and average precision, which is averaged fraction of correct classification for all climate classes. 220 Using the same validation datasets and station selection process, we also evaluated the previous climate maps from Beck et al., (2018) Kriticos et al., (2012), Peel et al., (2007), and Kottek et al., (2006). We applied the same Köppen-Geiger climate classification criteria described in the preivous studies to assess the overall accuracy of the map products. To further validate the climate classification results, we performed sensitivity analysis on the data integration method, the climate classificaiton time scale, and climatology dataset input, using the same validation datasets from GHCN-D and GSOD. In addition, we 225 compared the climate classification results with forest cover and elevation maps, and with the two high-resolution comparable climate map products, Beck et al., (2018) (1-km) and Kriticos et al., (2012) (0.167 o ), at regional and continental scale. The forest cover map we used is the 2000 30m Landsat-based forest cover map (Hansen et al., 2013). The elevation data is from the NASA SRTM Digital Elevation 30m data (Farr et al., 2007).   Regional distributions of climatic conditions are largely created by local variation in topography in rugged terrain (Dobrowski et al., 2013;Franklin et al., 2013). The climate classification and confidence level maps of mountainous areas of Central Rocky Mountains and Tibetan Plateau are shown in Figure 4 and 5 respectively. For each combination of precipitation and surface 255 air temperature datasets, we generated a Köppen-Geiger climate classification map (see Fig. 4a and 5a for 1979-2008 maps for the central Rocky Mountains and Tibetan Plateau). The final Köppen-Geiger classification map is derived based on the most common climate type among all the climate maps ( Fig. 4b and 5b). We then calculated corresponding confidence levels to quantify the uncertainty in the classification maps ( Fig. 4c and 5c). The uncertainty in climate classification in mountainous areas is attributed to the uncertainty existing in climate data, espectially precipitation data. In rugged terrain, CHELSA 260 precipitation data shows more detailed precipitation patterns, causing disagreement in classificaion results of the 3 rd level climate classes which depict precipitation seasonality.  (1979-2008, 1980-2009, 1981-2010, 1982-2011, 1983-2012, 265 1984-2013). (a) Small-scale accuracy of historical Köppen-Geiger climate maps. (b) Small-scale precision of historical Köppen-Geiger climate maps. Climate classificaiton has been applied for each station. The small-scale accuracy and precision are calculated based on the classification results of all the stations within the given region, with a minimum of 3 stations in the 5° search radius.

Results and Discussion
We validated the historical climate maps using the station observations from Global Historical Climatology Network-Daily 270 (GHCN-D) (Menne et al., 2012) and Global Summary Of the Day (GSOD) database (National Climatic Data Center et al., 2015). Figure 6 shows the small-scale distributions of total accuracy and average precision for historical Köppen-Geiger climate map series with 10° grid cells. Due to uneven distributions of weather stations, remote areas in the Pacific islands, Central Africa, and Amazon Forest suffer from a lack of station observations or an underrepresented validation results. Table 3 Continental and global overall accuracy, average precision, and confidence level of the historical Köppen-Geiger climate map series (1979-2008, 1980-2009, 1981-2010, 1982-2011, 1983-2012, 1984-2013). The overall accuracy is calculated as the percentage of correct climate classes using ground observations, and average precision is averaged fraction of correct classification for all climate classes. Confidence level values shows the 95% confidence interval of the confidence level for each continent, and the whole globe. All the values are presented in percentage.

Region
Africa We summarized the overall accuracy, average precision, and confidence levels for each continent and the whole globe ( and precision values, the continental average confidence levels range from 91.55% to 94.93%, and the global level is 92.90% (Table S2). Overall, the spatial patterms of total accuracy and average precision show good correspondence with classification confidence levels (Figure 3), indicating a potential of confidence level to represent classification uncertainty. Using the same validation datasets from GHCN-D and GSOD, we tested sensitivity of the climate map series using different combinations of temperature and precipitation dataset, and different method of data integration (Table 3). Results indicated 285 an average total accuracy of the 1km Köppen-Geiger classification maps generated with all the CHELSA, donwscaled CRU, GPCC and UDEL datasets and with only downscaled CRU, GPCC, UDEL datasets as 82.39% and 81.17% respecively.
Using the mean of multiple datasets which can potentially reduce the data bias, led to better classification results. We estiomtated the total accuracy of the previous high resolution Köppen-Geiger climate map products using the same validation datasets. We applied the same classification system described in the previous studies and the same time period of 290 the previous climate map product to process the station observation data and estimate their overall accuracy. Compared with the previous high resolution Köppen-Geiger climate map products, Beck et al. (2018) and Kriticos et al., (2012), the newly generated Köppen-Geiger climate map series showed greater accuracy in total.

295
We conducted sensitivity analysis of the Köppen classification scheme and tested multiple time scales, 10-yr, 20-yr, and 30yr. The selection criteria of station observations were adjusted accordingly based on the time scale utilized. Accuracy results exhibited decreasing accuracy for shorter time scale (Fig. 7). Further, we estimated the total accuracy for the Köppen-Geiger climate classification maps from previous studies, Beck et al., (2018) Kriticos et al., (2012), Peel et al., (2007), and Kottek et al., (2006), using the same validation dataset and consistent Köppen-Geiger climate classification scheme the corresponding 300 study applied. The validation results demonstrate that the new Köppen-Geiger maps have comparatively higher overall accuracy than all the previous studies.  (Hansen et al., 2013). The elevation data is the NASA SRTM Digital Elevation 30m data (Farr et al., 2007). The representative period of each map is listed in parentheses. At the regional and continental scale, we compared our Köppen-Geiger climate classification maps with previous map products for regions with large spatial gradients in climates, including central and eastern Africa, Europe, North America, and regions 320 with sharp elevation gradients, including Tibetan Plateau, central Rocky Mountains, central Andes (Fig. 8). We compared the new 1-km Köppen-Geiger climate classification maps from our study for time periods of 1980-2009, and 1984-2013 with the high-resolution Köppen-Geiger maps from two previous studies, Beck et al., (2018), which has a resolution of 1-km and temporal coverage of 1980, and Kriticos et al., (2012, which has a resolution of 0.0167o and covers . The Köppen classifications demonstrate good correlation with natural landscape distributions (Belda et al., 2014;Köppen, 1936;325 Trewartha, 1954). To show the agreement between the improved Köppen-Geiger climate classification maps and regional lanscape distributions, we also showed maps of forest cover, and elevation distribution for these regions. Figure 8 illustrate the enhanced regional details of the maps.

Regional and continental scale comparison
Compared with the Köppen-Geiger climate maps from previous studies with only one time period, the series of the Köppen-Geiger climate maps from our study demonstrate the ability to capture recent changes in spatial distributions of climate zones. 330 For example, our maps can detect the significant changes in the climate zones specifically driven by the accelerated global warming since the 1980s, for example, the poleward movements of boreal (D) and polar (E) climates in high latitudes in North America shown in the comparison between the 1980-2009 and 1984-2013 Köppen-Geiger climate maps (Fig. 8d). Another example is the expansion of savanna (Aw) climate into temperature (Cw) climate zone, witnessed in Central Africa (Fig. 8e).
Another improvement of the new series of the Köppen-Geiger climate maps is the application of threshold of -3 oC as the 335 boundary of temperate (C) and boreal (D) climate zones, which show better agreement with global boreal forest distributions at regional scale compared with Russell's modification of 0 oC (1931), which Beck et al., (2018), and Kriticos et al., (2012) utilized. Based on the comparison results of the Köppen climate zones and the biome classifications from the World Wildlife Federation , the boreal (D) climate zone largely corresponds to the distribution of boreal forest (Cui, Liang, & Wang, 2021). For example, evidenced in Figure (Fig. 8b). The abrupt changes in climate along the edges of the Andes mountains are also well described in the new maps (Fig.   8f).
In addition, the distribution of tropical (A), temperate (C) and boreal(D) climate zones in the new Köppen-Geiger maps correspond closely with tree lines in the forest cover maps. The temperate (C) and boreal(D) climate distributions based on the 350 Köppen-Geiger maps show a better agreemenet with the forest distributions of the Middle and Southern Rocky Mountains than Beck et al., (2018), and Kriticos et al., (2012) (Fig. 8a). For another example, the boundaries of the tropocal rainforest in Central Africa and South America are clearly delineated in the in the new Köppen-Geiger maps (Fig. 8e and 8f). Beyond the Köppen-Geiger climate classification maps, we calculated a set of bioclimatic variables from the monthly climate 355 data (see full list in Table 5). The bioclimatic variables at 1-km sptial resolution can capture regional environmental variations expecially in mountainous areas and areas with strong climate variations. These bioclimatic variables can be used in studies of environmental, agricultural and biological sciences, for example, development of species distribution modeling and assessment of biological impacts induced by climate change. The variables provide descriptions of annual averages, and seasonality of climates. The warmest half year or the coldest half year is defined as the period of the warmest six months or 360 the coldest six months. We validated the bioclimatic variables from different datasets with station data from GHCN-D (Menne et al., 2012) andGSOD database (National Climatic Data Center et al., 2015) (Fig. 9). We calculated a linear regression model for the 12 bioclimatic variables for each 10° grid cell (Fig. 10). The 30-yr average mean annual temperature (MAT) from CHELSA dataset shows 370 overall highest fit with station data, with CRU, and UDEL datasets showing smaller, but still strong correlation with station data. The 30-yr average mean annual precipitation (MAP) estimates from GPCC, UDEL, and CHELSA datasets have considerable uncertainties, indicated by relatively low correlation with station observations. In current precipitation datasets, there exist a varied degree of discrepancy in annual estimates over multiple time scales (Sun et al., 2018). 375 Figure 10. Small-scale comparison of annual temperature (MAT) and mean annual precipitation (MAP) variables derived from different datasets with station data. Small-scale correlation between the 30-yr average mean annual temperature (MAT) and mean annual precipitation (MAP) data and ground observations for three historical periods (1979-2008, 1981-2010, 1983-2012). The station data is from GHCN-D and GSOD database. The figure shows the R 2 value for 10° grid cells.    Future climate classifications derived from the diverse GCM projections for four RCPs, which are inherently uncertain (Gleckler, Taylor, & Doutriaux, 2008;Winsberg, 2012), provide a proxy of global distributions of climatic conditions and can represent potential spatial changes in climate zones under global warming. The large uncertainty and strong disagreement in projected climate classification maps at high latitudes and in regions with rugged terrain can be indicated by relatively low 405 confidence levels. Figure 12 and 13 show the future Köppen-Geiger climate classification maps based on GCM projections under RCP8.5 and associated confidence levels for the central Rocky Mountains and Tibetan Plateau. We generated a future Köppen-Geiger climate classification map for each bias-corrected and downscaled CMIP5 GCM projection (see Fig. 12a and 13a for 2070-2099 maps for the central Rocky Mountains and Tibetan Plateau). Noticeable regional changes in climate zones have been projected by comparing the 2070-2099 and 1979-2008 climate classification maps (see Fig. 12b and 12c for the 410 central Rocky Mountains, and Fig. 13b and 13c for Tibetan Plateau). Changes in climatic conditions under global warming have significant impacts on biodiversity and ecological systems. Area changes of climate zones can indicate spatial shrinkage or expansion of analogous climatic conditions, potentially implying threats for species range contraction or opportunities for range expansion (Cui, Liang, & Wang, 2021). To examine the area 420 changes of climate zones, we calculated the total area covered by each climate type for each historical and future periods under high-emission RCP8.5 scenario (Fig. 14). Our results of changes in area occupied by different climate zones demonstrate good agreement with results from previous studies (Chan & Wu, 2015). Results show that accelerated anthropogenic global warming since the 1980s has caused large-scale changes in climate zones and the shifts into warmer and drier climates are projected in this century. The tropical and arid climates are expanding into large areas in mid latitudes 425 whereas the high-latitude climates will experience significant area shrinkage.

Conclusion
Changes in broad-scale climatic conditions, driven by anthropogenic global warming, lead to the redistribution of species diversity and the reorganization of ecosystems. Distributions of the Earth's climatic conditions have been widely characterized based on the Köppen climate classification system. The Köppen climate classification maps require fine resolutions of at least 430 1-km to detect relevant microrefugia and promote effective conservation. Studies examining recent and future interannual or interdecadal changes in climate zones at regional scale needs more accurate depiction of fine-grained climatic conditions, continuous and longer temporal coverage.
We presented an improved long-term Köppen-Geiger climate classification map series for six historical 30-yr periods in 1979-2013 and four future 30-yr periods in 2020-2099 under RCP2.6, 4.5, 6.0 and 8.5. To improve the classification accuracy and 435 achieve a resolution as fine as 1-km, we combined multiple datasets, including WorldClim V2, CHELSA V1.2, CRU TS v4.03, UDEL, GPCC datasets and bias-corrected downscaled CMIP5 model simulations from CCAFS. The historical climate maps are based on the most common climate type from an ensemble of climate maps derived from combinations of observational climatology datasets. The future climate maps are based on an ensemble of climate maps derived from 35 GCMs. We estimated the corresponding confidence levels to quantify the uncertainty in climate maps. We also calculated 12 bioclimatic variables 440 at the same 1-km resolution using these climate datasets for the same historical and future periods to provide data of annual averages, seasonality, and stressful conditions of climates.
To validate the Köppen-Geiger climate classification maps, we used the station observations from GHCN-D and GSOD database. Our validation results show that the new Köppen-Geiger maps have comparatively higher overall accuracy than all the previous studies. Although the new maps exhibit improved overall accuracy, relatively lower confidence level and larger discrepancy in classification results are found especially in mountainous regions and major climate transitional zones located in mid and high latitudes. The confidence levels can provide a useful quantification of classification uncertainty.
Compared with climate maps from previous studies with a single present-day period, the series of the Köppen-Geiger climate maps from our study demonstrate the ability to capture recent and future projected changes in spatial distributions of climate zones. On regional and continental scale, the new maps show accurate depictions of topographic features and correspond 450 closely with vegetation distributions. Our Köppen-Geiger climate classificaion maps can offer a descriptive and ecological relevant way to provide insights into changes in spatial distributions of climate zones.
One of the limitations is that the future Köppen-Geiger climate maps built on dowscaled climate model projections exsit unaviodable uncertainties. The classification agreement levels of GCMs are relatively low at high latitudes and in regions with rugged terrain. The main sources of model discrepancies and uncertainties are deficiencies in model physics and varied model 455 resolution. The climate model outputs have coarse spatial resolution varying from 70-400 km and cannot well represent future climate change at the same scale of 1-km as our baseline climatology. Through bias-corretion and donwscaling methods, we made assumptions that local relationships between climatic variables remain constant across different scales, leading to a compromise between spatial scale and climate model physics.
We also tested the sensitivity of classification results to different time scale, dataset input, and data integration methods. Results 460 show that 30-yr time scale exhibited the highest accuacy results. Moreover, using the mean of multiple datasets from CHELSA, CRU, UDEL, and GPCC could lead to better classification results. Last, we provided a heuristic example which used climate classification map series to detect the long-term area changes of climate zones, showing how the new Köppen-Geiger climate classification map series can be applied in climate change studies. With improved accuracy, high spatial resolution, long-term continuous time coverage, this global dataset of the Köppen-Geiger climate classification and bioclimatic variables can be 465 used to in conjunction with species distribution models to promote biodiversity conservation, and to analyse and identify recent and future interannual or interdecadal changes in climate zones on a global or regional scale.