# A Köppen-Geiger classification derivative tailored for numerical modeling of ecosystems within watershed dynamics

Benoit Chauveau, Arnaud Pujol, Ahmed N Mama

► To cite this version:

Benoit Chauveau, Arnaud Pujol, Ahmed N Mama. A Köppen-Geiger classification derivative tailored for numerical modeling of ecosystems within watershed dynamics. Science of the Total Environment, 2026, 1015, pp.181445. 10.1016/j.scitotenv.2026.181445. hal-05491403

HAL Id: hal-05491403

https://ifp.hal.science/hal-05491403v1

Submitted on 3 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Research Paper

# A Köppen-Geiger classification derivative tailored for numerical modeling of ecosystems within watershed dynamics

Benoit Chauveau<sup>D*</sup>, Arnaud Pujol<sup>D1</sup>, Ahmed N. Mama

IFP Energies Nouvelles, 1 et 4 avenue de Bois-Preau, Rueil-Malmaison, 92852, France

<figure><img src="image_2.png" /></figure>

## HIGHLIGHTS

* A new framework is proposed to generate climate scenarios based on Köppen-Geiger classes.

* Clustering reveals distinct climatic patterns (or subclasses) within each Köppen-Geiger class.

* Pattern differences are linked to regional and climatic signatures.

* Subclasses are merged into a limited set of representative climate groups.

* A lake model simulation demonstrates the applicability of the approach.

## GRAPHICAL ABSTRACT

### Illustration of the general workflow

Generating ready-to-use climate scenarios for different climatic contexts and designed for various watershed dynamics modeling

<figure><img src="image_3.png" /></figure>

## ARTICLE INFO

Dataset link: http://dx.doi.org/10.5281/zenodo.016419447

**Keywords:**

Köppen-Geiger

Climatic classification

Environment

Climate change

Hydrodynamics model

## ABSTRACT

The dynamics of watersheds result from the interplay of multiple, coupled hydrological and ecological processes, which are strongly modulated by climatic conditions. Addressing these processes through numerical modeling often requires the generation of climate scenarios that accurately reflect short-term climatic variability in key parameters such as air temperature, precipitation, solar radiation, wind speed, relative humidity and cloud cover. In this study, we introduce a new framework to generate climate scenarios that capture the specific climatic features of the regions represented by the Köppen-Geiger classification. To this end, clustering methods were applied to identify climatic curve patterns, leading to the subdivision of each Köppen-Geiger class into subclasses defined by distinct temporal patterns, resulting in a total of 108 subclasses. A complementary analysis was then conducted to explore whether these pattern differences could be related to partial climatic signatures or specific regional contexts. This led to the reorganization of these subclasses into 10 groups, for which such characteristic signatures were proposed.

Climatic scenarios can be directly derived from the newly generated subclasses, providing ready-to-use climate forcings for hydro-environmental numerical models. The approach was illustrated through a case study on Lake Mendota using one of the generated scenarios. Result quality is evaluated by comparing the errors between vertical profiles simulated using the climatic scenario and those simulated using observed climatic data. Finally, the simulated primary productivity remains close to the reference value, with an

<footnote>* Corresponding author.</footnote>

<footnote>E-mail addresses: benoit.chauveau@ifpen.fr (B. Chauveau), arnaud.pujol@ifpen.fr (A. Pujol).</footnote>

<footnote n="1">¹ This author contributed equally to this work.</footnote>
approximate difference of 6 %. This proof-of-concept experiment, which calls for further applications, shows that the framework reproduces the key environmental drivers required for hydro-ecological and environmental modeling. Its extension across diverse climatic contexts therefore holds strong potential for predictive purposes.

# 1. Introduction

Over the past decades, a range of tools dedicated to the numerical simulation of hydro-environmental systems has emerged, with the aim of improving our understanding of how these systems respond to projected climatic and anthropogenic forcings. Among the available models, some require input climate data with a high temporal resolution, typically on a daily basis, or even hourly. This is the case, for example, for watershed dynamics models such as SWAT (Arnold et al., 1998; Neitsch et al., 2011), MIKE SHE (Abbott et al., 1986; DHI, 2012), VIC (Liang et al., 1994; Hamman et al., 2018), SHETRAN (Ewen et al., 2000; Boorman et al., 2007), or lake system dynamics models such as GLM (Hipsey et al., 2019) or FLAKE (Mironov et al., 2010). In these studies, a relatively straightforward approach is employed: climate data are either sourced directly from recent measurements (when aiming to reproduce observed dynamics) or extrapolated over one or two years based on climatic records from the years preceding the simulated period.

It becomes inherently more challenging to adapt climatic parameters when the study period lies in the distant past without direct climate records, or when aiming to develop predictive approaches at multi-decadal or longer timescales. In this context, it is common to represent climatic conditions (estimated either through numerical approaches or indirect field observations) using maps of climatic classes. This raises the central question of whether a consistent climate scenario for modeling ecosystems within watershed dynamics can be derived from a single climatic class representation.

Climatic classifications provide a coherent framework for representing a limited set of climate types. Several systems exist, each defined by characteristics tailored to its specific objectives. The first quantitative classification was introduced by Köppen (Köppen, 1900) and later refined into a widely accepted version (Köppen, 1936), commonly known as the Köppen-Geiger classification. This approach defines climatic classes based on air temperature and precipitation regimes. Since then, several modifications have emerged, giving rise to new classifications such as the Köppen-Trewartha system (Trewartha and Horn, 1980), which redefined the criteria while preserving the original philosophy of Köppen-Geiger (Belda et al., 2014). More differentiated approaches have also been proposed, incorporating additional climatic parameters to define classes—for example, the classifications by Thornthwaite (1948), Holdridge (1967), which consider potential evapotranspiration. Finally, it is worth noting that numerous classifications based on atmospheric circulation patterns have also been developed (see the review by Huth et al. (2008)).

The Köppen-Geiger classification remains the reference standard for several reasons. It has underpinned numerous climatic studies for over a century, providing extremely valuable spatial and temporal coverage. Its relatively simple applicability is also a key factor in its success, as it allows the reconstruction of class distribution maps from climate records, such as temperature and precipitation, both readily available from various open-access databases over the past fifty years. Finally, its well-established link with vegetation and ecosystems makes it highly effective for interdisciplinary approaches.

Using a climatic database for the 1951–2000 period, Kottek et al. (2006) produced the first digital Köppen-Geiger world map. Beck et al. (2018, 2023) later produced high-resolution Köppen-Geiger maps for various 30-year periods (1991–2020, 2041–2070, 2071–2099). The first period relies on a compilation of an extensive database integrating multiple sources, whereas the later periods are derived from climate model simulations designed to represent different IPCC climate scenarios.

Although climate models are useful to generate Köppen-Geiger maps for future periods, they can also be applied to reconstruct climatic conditions for past geological periods, thereby providing the necessary inputs for producing Köppen-Geiger paleomaps. This approach has already been applied in Chaboureau (2012), Chumakov et al. (1995), Scotese (2000), and more recently in Burgener et al. (2023).

Building a climatic scenario representative of a given climatic region and suitable for hydro-environmental modeling requires linking this region to the characteristic patterns of several environmental variables. In this paper, a pattern is defined as set of parameters represented by a time series over a relevant period. The climatic parameters identified include air temperature and precipitation, but also solar radiation, cloud cover, relative humidity, and wind speed. To our knowledge, no such link exists with climatic classifications such as Köppen-Geiger or its derivatives. The aim of this study is precisely to address this gap. To this end, the exploitation of a climatic database forms the cornerstone of this study.

The first part introduces the methodological approach and outlines the key objectives. The second part develops the database sampling strategy. The third section focuses on the initial stage of data processing, which leads to the construction of climatic patterns associated with the Köppen-Geiger classification. A key message emerging from this stage is that multiple climatic patterns may coexist within a single Köppen-Geiger class. Building on this observation, we propose a subdivision of the existing classes into climatic subclasses, each characterized by a specific climatic pattern. These emerging subclasses define a derivative of the Köppen-Geiger classification, which we consider to be particularly suited for the numerical modeling of ecosystems within watershed dynamics. In the fourth section, a complementary data-processing strategy is introduced in order to examine if specific climatic signatures can emerge from groups of subclasses, independently to the initial Köppen-Geiger class division. The fifth section is dedicated to illustrating an application of the derived Köppen-Geiger classes in the context of hydro-environmental modeling. To this end, we explore the modeling of the dynamics of a modern lacustrine system. Finally, a discussion of the strengths and limitations of the methodology is provided.

# 2. Objectives and methodological approach

The primary objective of this study is to develop a reproducible framework for generating realistic climatic scenarios suitable for hydro-environmental and ecosystem modeling, particularly in contexts where direct climate observations are sparse or unavailable (e.g., past or future periods). The proposed approach relies on present-day climate data to construct plausible climatic forcings, under the assumption that meaningful climatic analogs can be identified and exploited in a scientifically consistent manner. To this end, the framework is anchored in the Köppen-Geiger climate classification, which provides a widely used and structured representation of global climate regimes.

The methodological strategy follows three main stages:

* the selection and sampling of multi-variable climate datasets,

* the identification of representative climatic patterns within Köppen-Geiger classes,

* the generation of synthetic climate scenarios consistent with these patterns.

The selection of climate datasets is guided by two key principles: the availability of time series for a broad set of relevant climatic variables, and reproducibility through the use of open-access data sources. Based on these criteria, the NASA POWER database and MODIS satellite
products are used as the primary data sources. Details regarding dataset selection and sampling strategies are provided in Section 3.

Climate data are first organized according to the Köppen-Geiger classification, which distinguishes 30 climate classes grouped into five major categories (A-E). Permanently ice-covered regions (subgroup EF) are excluded, as they are poorly suited to conventional hydro-environmental modeling. The analysis therefore focuses on the remaining classes represented in the sampled dataset.

Within each Köppen-Geiger class, climatic information is analyzed at an annual scale, consistent with the temporal basis of the Köppen-Geiger criteria. Daily climate data are thus aggregated into twelve-month representations that describe the typical seasonal evolution of climatic variables. Preliminary inspection of these annualized curves suggests that multiple distinct climatic structures may coexist within a single class, motivating the need for a systematic, data-driven characterization of intra-class variability.

To address this issue, the framework incorporates multivariate statistical methods to identify coherent climatic patterns within each class. These patterns form the basis for defining subclasses, which serve as internally consistent climatic units for scenario construction. Additional classification steps are subsequently introduced to explore similarities between subclasses across different Köppen-Geiger classes, with the objective of providing a broader climatic context to support scenario selection.

The practical relevance of the proposed framework is illustrated through an application to the numerical modeling of a lacustrine system. A case study focusing on Lake Mendota (USA) is used to assess the ability of subclass-based climatic scenarios to reproduce lake dynamics simulated with the General Lake Model (GLM). This application serves as a proof of concept for the use of derivative Köppen-Geiger subclasses in hydro-environmental modeling.

Overall, the goal of this work is not to redefine the Köppen-Geiger classification, but to extend its applicability by associating each class with a set of representative climatic patterns and scenario-generation tools. The derived subclasses and the scripts used for climate scenario generation are made openly available via a Zenodo repository (DOI: 10.5281/zenodo.16419447).

## 3. Data sampling strategy

In this section, the data sampling strategy used to collect data for reconstructing annual variations of key climatic parameters for each Köppen-Geiger class is described.

Both temporal and spatial collections are based on the assumption that the most realistic distribution of Köppen-Geiger classes is provided by Beck et al. (2018, 2023).

### 3.1. Data sources

NASA provides an open-source dataset, the NASA POWER Project (2025) database, which records daily (and even hourly for some variables) variations of numerous climatic parameters over the past fifty years. This dataset is readily accessible via NASA's APIs, ensuring automated data retrieval. An additional advantage is that it offers a consistent set of variables (air temperature and precipitation, but also solar radiation, relative humidity, and wind speed) within a single framework. Using a single coherent database ensures internal consistency between parameters and avoids biases or incompatibilities that may arise when combining multiple sources.

The only exception concerns cloud cover data which is not directly available in the POWER project database but that can be deduced from Infrared (longwave) radiation through the Stefan-Boltzmann law. However, longwave radiation is strongly modulated by air temperature and humidity, as atmospheric emissivity largely depends on water vapor. Thus, it integrates the combined influence of cloud cover, temperature, and relative humidity—two of which are already included in

the dataset. For this reason, it provides less independent information than direct cloud cover observations, supporting our choice to source cloud data from a different database. Parameter independence is a critical aspect of the data-processing methods applied in this work and is further discussed in Section 4.

Satellite-based measurements from MODIS (Platnick et al., 2015) provide independent estimates of cloud cover. This source offers one of the most reliable datasets for satellite-based cloud retrieval. However, its temporal coverage is more limited than that of the NASA POWER database, which may raise several issues discussed in Section 3.4. To validate its use, additional checks were performed, including comparisons with the ERA5 reanalysis dataset (Hersbach et al., 2023).

The following subsections provide detailed descriptions of the specific use of these datasets.

### 3.2. Temporal data collection

Temporal sampling aligns with the 1991–2020 period of the dataset used by Beck et al. (2018, 2023) to create their reference present-day Köppen-Geiger map. This ensures that the extracted temporal variations of climatic parameters can be consistently associated with the Köppen-Geiger classes, which are themselves derived from climatological variations over the same 30-year period. This aspect is essential, as shorter sampling windows may alter the spatial distribution of Köppen-Geiger classes (Beck et al., 2005; Spinoni et al., 2015). Applied to our methodology, this means that maintaining consistency with the reference period is critical to ensure that the climatic patterns extracted from a given region are truly representative of the Köppen-Geiger class assigned to that region in the map of Beck et al. (2023).

Data were collected at a daily time step, corresponding to the upper range of time resolutions commonly used in hydro-environmental models, as mentioned in the introduction. This temporal resolution allows for the capture of fine sub-annual and even sub-monthly variations, which may contain characteristic climatic signatures specific to a given region.

As previously mentioned, NASA data provide hourly resolution for many climatic variables. However, several reasons justify our decision to retain a daily time step. The first is that adopting an hourly resolution would not improve the identification of climatic patterns using the method described in Section 4, since this identification relies exclusively on annual mean values. The second reason is the intention to preserve an annual scheme representative of a given climatic class, without it being distorted by high-frequency fluctuations that can vary substantially from one region to another within the same class. This choice does not preclude the possibility of superimposing in subsequent step high-frequency signals on the patterns, either on specific climatic parameters or across all of them, depending on the needs of the targeted hydro-environmental simulation. We specifically explore this aspect in Section 7.3, which focuses on the modeling of a lacustrine system. Such high-frequency fluctuations may reflect local geographic and vegetative conditions, but they are also strongly controlled by latitude—particularly in the case of solar radiation.

### 3.3. Spatial data collection

The spatial resolution we adopted aligns closely with that of the NASA POWER database. Specifically, precipitation data are provided at a resolution of 0.1°, while most meteorological variables are available at $1/2^\circ \times 5/8^\circ$ latitude/longitude. Primary solar radiation is an exception, offered at a coarser resolution of 1°. We assume that variations in solar radiation are smooth and continuous, and that a 1° resolution is sufficient to accurately capture them. Consequently, we consider that the NASA POWER database provides a spatial resolution of $1/2^\circ \times 5/8^\circ$. Accordingly, we sampled one point — defined by a latitude, a longitude, and a Köppen-Geiger index — for each $1/2^\circ \times 5/8^\circ$ grid cell, with
<figure><img src="image_4.png" /><figcaption>**Fig. 1.** Comparison between MODIS and ERA5 cloud cover datasets. Each point represents a measurement pairing of annual mean cloud cover, with MODIS values on the x-axis and ERA5 values on the y-axis. The dashed line corresponds to the 1:1 line, indicating perfect agreement between the two datasets. The solid black line represents the linear regression fit. The color scale used for the points follows the convention commonly applied in Köppen-Geiger climate maps.</figcaption></figure>

the only constraint being that it corresponds to the most representative Köppen index within that area. A grid cell is considered to be located in a continental domain if more than 50% of its pixels are classified as continental. Following this strategy, a total of 48,910 locations were ultimately sampled. Table 1 reports the distribution of the sampled Köppen indices alongside their distribution in the 1991–2020 Köppen–Geiger map. The overall distribution of indices is well preserved with this sampling approach. However, some classes — already sparsely represented in the 1991–2020 Köppen–Geiger map — are entirely absent from our dataset. This outcome simply reflects the fact that the spatial resolution of the NASA POWER database, along with its intrinsic data quality, is insufficient to adequately capture these less common climatic classes.

## 3.4. Cloud cover data and temporal sampling

As mentioned in the introduction, the NASA POWER database does not provide direct cloud cover data, whereas MODIS offers one of the most reliable datasets for cloud cover retrieval. The MODIS data discussed here are derived from the Level-3 MODIS Atmosphere Eight-Day Global Product, from which cloud cover was extracted for 48,910 sampled locations over the period 2001–2020. This freely available dataset provides a spatial resolution of 1° and a temporal resolution of 8-day intervals.

The absence of cloud cover data between 1991 and 2000 raises an important question regarding data consistency. We considered whether it would be relevant to extend the time series using another dataset. In this context, we explored the ERA5 reanalysis, which provides cloud cover data from 1980 to the present. However, these data are intrinsically less reliable than those derived from the MODIS program, as they are derived from a numerical weather model (IFS) that assimilates a wide range of observations, including satellite radiances, but does not come from direct satellite observations. A direct comparison between MODIS and ERA5 cloud cover revealed a systematic discrepancy, with the best-fit regression line being $CC_{ERA5} = 0.88 CC_{MODIS} + 0.07$ (Fig. 1). This indicates that ERA5 tends to underestimate cloud cover at higher values while slightly overestimating it at lower values, reflecting structural differences between the two datasets. Such differences likely

stem from the contrasting nature of the products: MODIS provides direct satellite-based retrievals at high spatial resolution, whereas ERA5 cloud cover is a model-diagnosed variable influenced by radiative transfer parameterizations and cloud microphysics schemes.

Merging these two datasets would therefore introduce artificial discontinuities in the time series, undermining the robustness of our analyses. In other words, complementing the dataset with ERA5 cloud cover to span the target period (1991–2020) would not provide any added value.

Another way to retrieve cloud cover data is to estimate it from longwave radiation. It is therefore legitimate to examine the correlation between MODIS-derived cloud cover and cloud cover estimated from longwave radiation. To do so, we apply a standard empirical formula derived from Crawford and Duchon (1999) and Brutsaert (1975). Without detailing the formula, it is worth noting that the relationship also involves air temperature and relative humidity to link cloud cover with longwave radiation. Fig. 2 highlights the strong correlation between these two parameters when considering, for each spatial sampling point, the average values over the studied time period. This first indicates that the shorter time period covered by the MODIS data compared to the full 1991–2020 period does not appear to introduce bias in the data. This observation reinforces the relevance of using the MODIS dataset despite its shorter sampling period. We emphasize here that, although an excellent correlation was found between cloud cover estimated from longwave radiation in the NASA POWER database and MODIS cloud cover, MODIS remains the preferred choice. Indeed, while this result supports the use of longwave-derived cloud cover estimates in the absence of direct observations, such a correlation does not capture the potentially complex, nonlinear relationships between cloud cover and other temporal variables (temperature, relative humidity). Since the approach introduced in this paper relies on clustering algorithms (current Section 4) based on multivariate distance metrics rather than on simple linear relationships between variables, the use of MODIS cloud cover data is therefore more appropriate.

## 3.5. Data consistency

Reconstructing temperature and precipitation curves from a different database and at a different spatial resolution scale inevitably
<caption>**Table 1**</caption>

Sampling of the reference Köppen-Geiger map and percentage of indices modified if re-estimated from the NASA POWER database.

<table><thead><tr><th>Index</th><th>Reference<br>(1991–2020)</th><th>Initial</th><th># Points</th><th>% Change</th></tr></thead><tbody><tr><td>1</td><td>3.73%</td><td>3.72%</td><td>1820</td><td>20.77%</td></tr><tr><td>2</td><td>2.69%</td><td>2.67%</td><td>1307</td><td>42.39%</td></tr><tr><td>3</td><td>9.77%</td><td>9.87%</td><td>4828</td><td>11.33%</td></tr><tr><td>4</td><td>13.29%</td><td>13.31%</td><td>6511</td><td>12.27%</td></tr><tr><td>5</td><td>5.03%</td><td>4.96%</td><td>2426</td><td>40.93%</td></tr><tr><td>6</td><td>4.95%</td><td>4.98%</td><td>2434</td><td>22.68%</td></tr><tr><td>7</td><td>6.04%</td><td>6.00%</td><td>2934</td><td>51.91%</td></tr><tr><td>8</td><td>0.95%</td><td>0.99%</td><td>485</td><td>22.47%</td></tr><tr><td>9</td><td>0.63%</td><td>0.59%</td><td>287</td><td>40.07%</td></tr><tr><td>10</td><td>—</td><td>—</td><td>—</td><td>—</td></tr><tr><td>11</td><td>2.43%</td><td>2.40%</td><td>1174</td><td>37.48%</td></tr><tr><td>12</td><td>0.93%</td><td>0.88%</td><td>428</td><td>50.00%</td></tr><tr><td>13</td><td>—</td><td>—</td><td>—</td><td>—</td></tr><tr><td>14</td><td>3.83%</td><td>3.82%</td><td>1868</td><td>14.40%</td></tr><tr><td>15</td><td>2.24%</td><td>2.24%</td><td>1098</td><td>28.32%</td></tr><tr><td>16</td><td>0.08%</td><td>0.03%</td><td>15</td><td>46.67%</td></tr><tr><td>17</td><td>0.17%</td><td>0.13%</td><td>63</td><td>30.16%</td></tr><tr><td>18</td><td>0.46%</td><td>0.44%</td><td>213</td><td>34.74%</td></tr><tr><td>19</td><td>1.90%</td><td>1.76%</td><td>863</td><td>97.68%</td></tr><tr><td>20</td><td>0.01%</td><td>0.01%</td><td>4</td><td>100.00%</td></tr><tr><td>21</td><td>0.78%</td><td>0.79%</td><td>385</td><td>12.21%</td></tr><tr><td>22</td><td>1.34%</td><td>1.28%</td><td>625</td><td>39.04%</td></tr><tr><td>23</td><td>2.54%</td><td>2.48%</td><td>1211</td><td>28.16%</td></tr><tr><td>24</td><td>0.31%</td><td>0.31%</td><td>154</td><td>93.51%</td></tr><tr><td>25</td><td>1.33%</td><td>1.31%</td><td>639</td><td>7.67%</td></tr><tr><td>26</td><td>8.23%</td><td>8.17%</td><td>3994</td><td>17.73%</td></tr><tr><td>27</td><td>17.52%</td><td>17.73%</td><td>8673</td><td>7.28%</td></tr><tr><td>28</td><td>0.16%</td><td>0.15%</td><td>75</td><td>49.33%</td></tr><tr><td>29</td><td>8.65%</td><td>8.99%</td><td>4396</td><td>11.67%</td></tr><tr><td>30</td><td>—</td><td>—</td><td>—</td><td>—</td></tr></tbody></table>

leads to situations where the Köppen-Geiger classification obtained from Beck's dataset does not match the classification derived from an alternative source. In our case, we observe a change in the Köppen classification for 21% of the sampled locations (see Table 1 and Fig. 3). To assess the impact of spatial resolution, we evaluated the inconsistency rate between Köppen-Geiger classifications from Beck's maps at different spatial resolutions (0.01°, 0.1°, and 0.5°). We found discrepancies in classification of around 10%. We also examined which reference map best aligned with the classifications derived from NASA POWER data, and found that the 1 km resolution map provided the most consistent match (following the strategy which consists in associating to a sampling window the most representative Köppen-index found in the map).

A closer examination of these inconsistencies reveals that certain classes (particularly those defined by strong seasonal contrasts in precipitation between summer and winter) are especially difficult to reproduce using the NASA POWER database.

Several options could have been considered to handle these inconsistencies. A first approach would be to restrict the analysis to

areas where the Köppen-Geiger indices are correctly reproduced using NASA POWER data. However, this would lead to the exclusion of a significant portion of the dataset. A second option would have been to reconstruct a Köppen-Geiger map entirely from NASA data in order to maintain full consistency between the collected climate variables and the Köppen-Geiger indices provided in Beck et al. (2018, 2023). However, the potential benefit of such consistency does not justify the loss of accuracy in describing the distribution of climate classes. We therefore consider maintaining the Köppen-Geiger classification as provided in Beck et al. (2018, 2023) as the most relevant choice. Indeed, it appears more important to focus on accurately capturing the variability of climate parameters — beyond temperature and precipitation — within each class, rather than modifying the classification itself.

Finally, it is worth noting that many of the locations where the Köppen-Geiger indices derived from NASA data differ from the reference map fall within climatic transition zones, where temperature and precipitation conditions place them at the boundary between multiple classes (Fig. 3-c)
<figure><img src="image_5.png" /><figcaption>Fig. 2. Comparison between computed longwave radiation (based on annual mean cloud cover from MODIS) and longwave radiation from the NASA POWER database. Each point represents a pair of values for a given location, with longwave radiation calculated from cloud cover on the x-axis and NASA POWER values on the y-axis. The dashed line corresponds to the 1:1 line, indicating perfect agreement. The solid black line represents the linear regression fit. The color scale for the points follows the convention commonly used in Köppen-Geiger climate classification maps.</figcaption></figure>

## 4. Derivated Köppen-Geiger classes and climatic scenarios

Although Köppen-Geiger classes are commonly treated as inter-nally homogeneous, the analysis of annualized climatic curves suggests that multiple distinct patterns may coexist within a single class. Fig. 4 illustrates this point for the BSh class, chosen here as a repre-sentative example. For each climatic variable, raw daily data (blue) aggregated over multiple years are shown together with polynomial fits (red) describing the mean annual seasonal cycle. While the ensemble of curves delineates a recognizable class-level pattern, the dispersion observed around these mean trajectories indicates the presence of struc-turally different sub-patterns. This motivates a systematic, data-driven identification of intra-class climatic variability.

To formalize this observation, we apply a classification framework combining Principal Component Analysis (PCA) and k-means clustering to annual mean values of six independent climatic variables: precip-itation (P), air temperature (T), wind speed (WS), cloud cover (CC), relative humidity (RH), and shortwave radiation (SR). Variable inde-pendence is preserved to avoid artificial correlations that could bias subgroup identification. To improve discrimination between patterns, the climatic dataset is augmented with geographic and topographic descriptors (latitude (LAT), longitude (LON), distance to coastline (DC), and elevation (Z)).

The number of subgroups within each Köppen-Geiger class is first estimated using the elbow method, which identifies diminishing re-turns in within-cluster variance reduction. This heuristic estimate is then constrained by a variance-explained criterion, requiring that the selected clustering configuration account for 80%-90% of the total variance. When necessary, the number of clusters is adjusted to satisfy this condition, ensuring a balance between parsimony and statistical representativeness.

Applying this procedure to all sampled Köppen-Geiger classes yields 108 distinct clusters, each representing a coherent family of climatic curves. Internal consistency within clusters is supported by correlation matrices showing stable co-variations among climatic variables. In the specific case of the BSh class, five subgroups are identified (Fig. 4). While the unclassified data already display recognizable seasonal structures, clustering reveals that different subsets of locations are more

accurately represented by distinct climatic patterns. This subdivision also retains a degree of geographical organization, as illustrated in Fig. 5. Equivalent maps for all subclasses are provided in the Supplementary Material.

These families of curves define a set of data-driven climatic sub-classes, referred to hereafter as *derivative Köppen-Geiger classes*. For each subclass, polynomial representations are fitted to the annualized curves of the six climatic variables (air temperature, wind speed, short-wave radiation, cloud cover, relative humidity, and precipitation). To preserve the periodic structure of the seasonal cycle, continuity con-straints are imposed on both value and slope at the December–January transition.

Each fitted variable *p* yields a mean seasonal curve

$$ \bar{p}_p(t) \in \mathbb{R}^T, $$

together with an empirical uncertainty envelope defined as

$$ \bar{p}_{p,\text{polymax}}(t) = \bar{p}_p(t) + c \sigma_{p,\text{res}}, $$

$$ \bar{p}_{p,\text{polmin}}(t) = \bar{p}_p(t) - c \sigma_{p,\text{res}}, $$

with confidence factor

$$ c = 1, $$

and residual standard deviation

$$ \sigma_{p,\text{res}} = \sqrt{\frac{1}{N} \sum_{i=1}^{N} \frac{1}{T} \sum_{t=1}^{T} (p_{p}^{(i)}(t) - \bar{p}_{p}(t))^{2}}. $$

Here, $p_{p}^{(i)}(t)$ denotes the observed value of variable *p* on day *t* for sample *i*, $\bar{p}_{p}(t)$ the fitted subclass-mean polynomial, *N* the number of samples in the subclass, and *T* = 365 the number of days per year.

Together with the inter-variable correlation matrix, these poly-nomial representations and uncertainty envelopes (confidence band represented by the green area in the derivative Köppen-Geiger classes (Fig. 4)) define the statistical structure of each derivative class and constitute the basis for the climate scenario generation framework introduced in the following section.
<figure><img src="image_6.png" /><figcaption>Fig. 3. (a). Sampling map displaying Köppen–Geiger classifications from (Beck et al., 2018). (b). Sampling map based on Köppen–Geiger classifications estimated using the POWER project database. (c). Locations where classifications differ between the two sources; colors correspond to the Köppen–Geiger classes derived from the POWER project data.</figcaption></figure>

## 5. Multi-constraint scenario generation

Let $\mathcal{P}$ denote the set of variables for which the user specifies an annual-mean constraint. The user may impose one or more annual-mean constraints $\{m_p\}$ on a subset of variables $p \in \mathcal{P}$. For each constrained parameter $p$, the generated time series $x_p(t)$ must satisfy the following annual mean condition:

$$ \frac{1}{T} \sum_{t=1}^{T} x_p(t) = m_p. $$

To generate time series that fulfill these constraints while remaining physically and statistically consistent, we model the full set of variables as a single concatenated vector:

$$ X = \begin{bmatrix} \mathbf{X}_1 \\ \vdots \\ \mathbf{X}_p \end{bmatrix} \in \mathbb{R}^{P \times T}, $$

where $\mathbf{X}_p \in \mathbb{R}^T$ is the daily values (e.g., 365 days) for variable $p$. where $\mathbf{X}_p \in \mathbb{R}^T$ denotes the daily values (typically 365 days) of variable $p$. We represent $X$ with a multivariate normal distribution to

account for cross-variable correlations and to allow conditioning on annual means. Although climate time series are not strictly Gaussian, this approximation remains useful for generating internally consistent scenarios:

$$ X \sim \mathcal{N}(\mu, \Sigma), $$

where the mean vector $\mu$ corresponds to the stacked class-average seasonal curves:

$$ \mu = \begin{bmatrix} \bar{p}_1 \\ \vdots \\ \bar{p}_p \end{bmatrix}, $$

and the covariance matrix $\Sigma$ is defined as:

$$ \Sigma = \Sigma_{\text{cross}} \otimes I_T. $$

which assumes temporal independence across days while preserving cross-variable dependence through $\Sigma_{\text{cross}}$. This formulation does not reproduce daily autocorrelation, but it ensures that generated scenarios remain consistent across variables and compatible with prescribed annual means.

Although day-to-day correlations are not directly encoded in the covariance structure (due to the use of $I_T$), the temporal evolution of each variable is still realistically captured via the polynomial fits used to define $\mu$. These smooth polynomial trends represent the typical seasonal dynamics observed historically for each climate variable. As a result, the generated scenarios exhibit coherent intra-annual behavior even in the absence of explicitly modeled temporal dependence.

Each element of the cross-variable covariance matrix $\Sigma_{\text{cross}}$ captures how two variables co-vary, based on their empirical correlation and respective residual variability:

$$ [\Sigma_{\text{cross}}]_{ij} = \rho_{ij} \sigma_{i,\text{res}} \sigma_{j,\text{res}}, $$

where $\rho_{ij}$ is the Pearson correlation coefficient between variables $i$ and $j$, and $\sigma_{i,\text{res}}$ is the residual standard deviation of variable $i$.

The Pearson correlation coefficient $\rho_{ij}$ is computed as:

$$ \rho_{ij} = \frac{\sum_{i=1}^{T} (x_i(t) - \bar{x}_i) (x_j(t) - \bar{x}_j)}{\sqrt{\sum_{i=1}^{T} (x_i(t) - \bar{x}_i)^2} \sqrt{\sum_{i=1}^{T} (x_j(t) - \bar{x}_j)^2}}, $$

where $\bar{x}_i$ and $\bar{x}_j$ are the mean values of the time series $x_i(t)$ and $x_j(t)$ over the $T$ time steps.

This formulation separates the dependency structure ($\rho_{ij}$) from the scale of variability of each variable. The correlation $\rho_{ij}$ measures the linear relation between two variables, and the product $\sigma_{i,\text{res}}$, $\sigma_{j,\text{res}}$ rescales the covariance according to the variability of each variable. Its modular structure allows individual variances or correlations to be modified without altering the overall structure.

In practice, this modeling choice allows us to encode the realistic cross-dependencies between climatic variables (e.g., temperature and relative humidity), while ignoring temporal autocorrelation. Temporal independence is enforced by assuming a block-diagonal covariance structure, $\Sigma = \Sigma_{\text{cross}} \otimes I_T$, where each time step is treated as statistically identical but multivariately coupled to other variables.

Let us partition $X$ into two blocks: $X_c$, the variables with imposed annual-mean constraints, and $X_u$, the unconstrained variables. Given a target realization $\mathbf{x}_c^*$ that satisfies the user-defined annual means, we compute the conditional expectation of the unconstrained variables as:

$$ \mathbb{E}[\mathbf{X}_u | \mathbf{X}_c = \mathbf{x}_c^*] = \mu_u + \sum_{uc} \Sigma_{cc}^{-1} (\mathbf{x}_c^* - \mu_c), $$

where $\Sigma_{uc}$ is the cross-covariance between $X_u$ and $X_c$.

This expression adjusts the unconstrained variables based on the deviation between the constrained realization and the class-mean curves, modulated by the correlations between variables.

By merging the constrained realization $\mathbf{x}_c^*$ and the corresponding conditional expectation for $X_u$, we obtain a full multivariate scenario $\hat{X}$ that satisfies all imposed annual-mean constraints exactly, preserves
<figure><img src="image_7.png" /><figcaption>Fig. 4. BSh Köppen-Geiger class. Inside the black rectangle: in blue, thirty years (twenty years for cloud cover) of raw data retrieved from the POWER project and the MODIS database, covering six climatic parameters and aggregated annually (SW: shortwave radiation; T: air temperature; RH: relative humidity; P: precipitation; WS: wind speed; CC: cloud cover); in red, the corresponding polynomial regression. Outside the black rectangle: derived classes, each exhibiting a distinct set of curves and associated confidence bands (green shaded areas), extracted from the BSh dataset.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 5. Geographical distribution of the various points belonging to the BSh class, categorized by their respective subclasses.</figcaption></figure>

observed inter-variable relationships and remains seasonally consistent with the class it belongs to.

To ensure physical plausibility, each time series is finally clipped to stay within its empirical envelope:

$$ \hat{x}_p(t) \in [\bar{p}_{p,polymin}(t), \bar{p}_{p,polymax}(t)]. $$

To illustrate the scenario generation process under user-defined constraints, we consider two cases based on sub-climate 3 of the BSh class (denoted BSh.3). In both examples, synthetic daily series are generated using the polynomial seasonal patterns fitted for this subclimate (see above). The reference annual mean values for BSh.3 are 22.58 °C for temperature and 1.46 mm d⁻¹ for precipitation.

* **Example 1: Temperature constraint only.** A single constraint is imposed on temperature, corresponding to a residual standard deviation of -1 (i.e., an annual mean temperature of 20.27 °C). The generated series (Fig. 6a) for all variables remain consistent with the correlation structure associated with this sub-climate (Fig. 6b).

* **Example 2: Temperature and precipitation constraints.** We impose two constraints this time: the same annual mean temperature as in Example 1, and a specified annual precipitation corresponding to a residual standard deviation of 1 (i.e., an annual mean precipitation of 2.55 mm d⁻¹). The resulting precipitation curve respects the prescribed mean, and the temperature curve remains consistent with the first case (Fig. 6c). Importantly,
<figure><img src="image_9.png" /><figcaption>**Fig. 6.** Derived BSh.3 scenario under constraints. Blue curves in panels (a) and (c) represent the default configuration corresponding to BSh.3, whereas the red curves depict constrained scenarios that remain within the confidence band (green shaded area). (a) Only temperature was constrained to its minimum level (residual standard deviation = -1). (b) Correlation matrix of subclimate BSh.3. (c) Temperature and precipitation were constrained to their minimum and maximum levels, respectively (residual standard deviation set to -1 for temperature and +1 for precipitation). (d) Correlation matrix of subclimate BSh.3. The red and blue squares highlight the correlation values between temperature and precipitation with the other climate variables, respectively.</figcaption></figure>

the variations observed in the other variables (e.g., relative humidity, shortwave radiation) differ from Example 1 in a manner that reflects the empirical correlations encoded in the correlation matrix of the subclimate BSh.3 (Fig. 6d). For instance, higher precipitation values tend to coincide with increased relative humidity and reduced radiation, in agreement with the typical inter-variable relationships observed.

These examples illustrate how the imposed mean constraints shape the generated scenarios, while the remaining variables are adjusted in a statistically coherent way, preserving the cross-variable relationships captured by the empirical correlation structure.

This flexible statistical framework enables consistent generation of multivariate annual scenarios, suitable for climate impact studies, sensitivity analysis, and synthetic data generation. The code used to generate these scenarios, along with the full set of derivative Köppen-Geiger classes, is available via the Zenodo linked provided in the Data Availability Statement section.

## 6. Transverse classification of the derived Köppen-Geiger classes

The proposed subdivision into subclasses is based on the observation that multiple climatic patterns can exist within a single Köppen-Geiger class. However, at this stage, in the context of a numerical approach — such as assessing the dynamics of a lake system based on a climatic scenario assuming a given Köppen-Geiger class — the method does not yet provide guidance on which subclass should be prioritized. In the context of multidecadal climate projections, it may be considered reasonable to assign, for a given area experiencing a change in climate class, the subclass that is most representative of that class at a macroregional scale. However, such a geographical analogy is inherently limited and lacks underlying physical justification. Moreover, this approach becomes even less relevant when the objective is to simulate a past geological period. In any case, it remains pertinent to analyze the characteristics of each subclass from a different perspective, in order to assess whether there are robust arguments for prioritizing certain

subclasses over others in a given context. To address this gap, we investigated whether a link could be established between the curve patterns and specific climatic contexts, whose descriptions could help refine the characterization of the proposed subclasses.

The approach implemented is again based on a clustering method, but the comparison criterion differs from that used for the initial subclass definition. In the previous phase, the patterns emerged from point clouds in the eigenspace of the principal component analysis (PCA), where the physical parameters were linked to the annualized mean values of climatic variables. In contrast, the goal of this phase is to investigate whether the patterns associated with a subclass of one Köppen-Geiger class also appear in a subclass of another class. To this end, we use the full climatic curves as input features and apply a hierarchical clustering method.

The hierarchical clustering method separates sets of curves using a dendrogram-based structure. The cut-off criteria that determine the number of groups can be defined based on various considerations. We employed several widely-used statistical methods to identify the optimal number of clusters. These methods include:

* **Elbow Method:** It is the same one we used to create the derivative Köppen-Geiger classification. But here, we refined this method by analyzing the second derivative of the curve to objectively detect the inflection point.

* **Silhouette Index:** This index measures the cohesion and separation of clusters by comparing the average intra-cluster distance with the distance to the nearest other cluster. The silhouette coefficient ranges from -1 to 1, and the optimal number of clusters is the one that maximizes the average silhouette score across all points.

* **Gap Statistic:** This method compares the observed within-cluster dispersion to that expected under a reference null distribution (typically uniform). The optimal number of clusters is chosen either as the one with the largest gap or as the smallest k for which the gap is within one standard error of the maximum.
The application of the aforementioned methods to our climatic dataset yielded the following results:

* The elbow method (using second derivative analysis) suggested 8 clusters.

* The silhouette index indicated an optimum at 11 clusters.

* The gap statistic pointed to 8 clusters.

Although the estimates differ slightly, they remain in close agreement, which reflects a consistent underlying structure in the data. However, none of these methods can guarantee a universally optimal clustering, especially in domains such as climatology where interpretability and practical significance are essential (it was the same reason we required 80% and 90% of the total variance in data for the construction of the derivative Köppen-Geiger classifications).

The selection of the number of clusters should be guided by expert knowledge and by the interpretability of the results. This requires assessing whether clusters delineate distinct climatic regimes, whether additional clusters reveal meaningful sub-structures, or whether further partitioning only adds complexity and noise. In our analysis, we retained a division into 8 groups (the dendrogram and the statistical method curves are provided in Supplementary Material, external link in Appendix A). To justify this choice, we first characterize these groups in relation to broadly recognizable macro-environmental contexts (e.g., atmospheric, oceanographic, or topographic). We then discuss the potential implications of applying a criterion that would yield a greater number of groups.

The distribution of subclass groups is shown in Fig. 7 and sum-up Table 2. A preliminary global qualitative observation already reveals that latitude emerges as a key factor influencing group membership. Furthermore, the spatial distribution of subclass groups bears a notable resemblance to that of classical macro-climatic zones (that is reflected by the first letter of each class in the Köppen-Geiger classification) without being strictly similar. This observation still strongly supports the proposition of associating each group with a macro-environmental signature.

Group 1 gathers regions characterized by a predominantly continental climate, potentially with some mountainous influences (the mean elevation is approximately 800 m, but with substantial variability across the area). This group is confined to two major regions: the semi-arid zones of South America, associated with the high Cordillera of the Andes, and the vast semi-arid belt of Central Asia and the Tibetan Plateau. It then extends onto the elevated Ganges plains of northern India, the Terai of Nepal, and across the western Tibetan Plateau into China's Xinjiang. It even reaches southern China, where mountainous highlands (the Yunnan-Guizhou Plateau) influence regional climate conditions. This group gathers several classes belonging to three different macro climatic classes in the Köppen-Geiger meaning: the temperate climate (C) through Cwa and Cfa, the continental climate (D) through Dwa and Dwc. However, these classes cannot be considered as representative of group 1, simply because it involves less than half of their associated subgroups, excepted Dwa class for which two over four subclasses are concerned. The analysis of the dendrogram shows that humidity first, and then cloud cover, are the main discriminant climatic factors to differentiate group 1 to the others, which is compatible with a mountainous context.

Group 2 exhibits characteristics typical of monsoon-influenced tropical climates. It encompasses a broad tropical to subtropical belt extending across Central and South America, West and Central Africa, South and Southeast Asia, and parts of southern China and Oceania. These regions are generally characterized by warm temperatures and substantial precipitation, often governed by monsoonal circulation patterns or persistent tropical humidity. Despite their geographic dispersion, these areas share similar climatic dynamics driven by marked seasonal variations in heat and moisture. This group differs from group 1 by its strong association with some classes which are Aw (three over five subclasses),

Bsh (four over five subclasses), Cwa (three over four subclasses) and all the Cwb subclasses. Therefore this group encompass what the Köppen-Geiger classification identifies as the subtropical highland climate or temperate oceanic climate, almost the humid summer subtropical Cwa, but also the arid climate (which can be cold or hot) located at a climatic transition zone. The analysis of the dendrogram shows that humidity first, and then cloud cover, are the main discriminant climatic factors to differentiate group 2 to the others, which is compatible with a monsoon influence.

Group 1 and group 2 form a set that differ from the other ones by their relative humidity annual variation. Groups 3 to 8 can be further subdivided into two sets based on their annual precipitation variation: groups 3 to 4 and groups 5 to 8. While groups 3 to 4 display a continental signature characterized by dry winters reflecting an anticyclonic influence that acts as a barrier to oceanic moisture, groups 5 to 8 exhibit summer dryness directly related to the dynamics of the subtropical anticyclones.

Group 3 encompasses 25 subclasses or around 25% of the subclasses previously identified. It totally includes the Bsk, Dfa and Dfb classes, and a large part of Cfa (four over five subclasses), Cwa (three over four subclasses) and Bsh (four over five subclasses). It includes a large portion of the United States and southern Canada, much of South America, the Sahara, South Africa, Eastern Europe, Central Asia, and southern Australia. Group 4 also encompasses several classes, totally including Dsd and Dfd classes, and a large part of Dsc, Dfc, Dwc and ET. It broadly encompasses the main portion of the continent at high latitudes in the Northern Hemisphere, as well as the central part of Australia. Main distinction between group 3 and group 4 is the cloud cover annual variability for which the strong and permanent arid signature is expressed as a low and relatively constant cloud coverage over a year. We thus associate group 3 to seasonal variations in contrast with the arid condition maintain for area belonging to group 4.

Groups 5 and 6 are exclusively associated with classes belonging to Category A (tropical climates). When combined with Group 2, they enable a finer differentiation of certain subclasses within Type-A climate classes. Group 5 is equatorial and dominates the central and western Amazon Basin and the Guianas. Its permanent and abundant atmospheric moisture reflects the persistent influence of the Intertropical Convergence Zone; Group 6 can be also categorized as an equatorial climate but with a more pronounced seasonal variability, influenced by seasonal atmospheric circulation (either the monsoon or the Intertropical Convergence Zone). It features high annual precipitation with slight seasonal variation, marked by a drier period linked to wind reversals and reduced convergence. This climate is found in regions like southern Amazonia, the Andean foothills, the Orinoco Basin, and parts of Southeast Asia.

Group 7 corresponds to climate dominated by a marine/oceanic influence. It mainly includes regions of Western Europe, the western coast of North America, and the southern portion of South America. The main signature of this class is present in temperate climates, although it can also be found in arid, continental, or polar climates according to the Köppen-Geiger classification. Group 8 is broadly associated with Mediterranean-like conditions, even if regions may be more attributed to more arid and continental conditions. It fully includes the Dsa class and the majority of the Csa, Csb, and Dsb classes. It encompasses regions surrounding the Mediterranean Basin, as well as parts of California, central Chile, southwestern Australia, and the Western Cape of South Africa, along with more arid areas such as the coastal region of the Arabian Peninsula and parts of Kazakhstan.

Hierarchical clustering thus enables the redistribution of subclasses into 8 groups, each of them being assigned to an identified macro-environment. These groups are as follows: high-altitude continental, monsoon-like tropical, seasonal continental, arid continental, equatorial with permanent Intertropical Convergence Zone (ITCZ) influence, seasonal equatorial, oceanic, and Mediterranean. While some traits of these groups may be redundant in certain classes, it is especially worth
<figure><img src="image_10.png" /><figcaption>Fig. 7. Geographical distribution of the eight groups that emerge from the hierarchical clustering methodology.</figcaption></figure>

noting that some subclasses may now exhibit a more specific signature. Returning to the case study of class BSh, three of its subclasses align naturally with the continental group, while the fourth exhibits a pattern influenced by climatic forcing akin to that observed in regions classified as Mediterranean.

The transversal classification based on these 8 groups allows to 23 classes to be distributed across at least two different groups. For the others classes (BSh, BSk, Dwd, Dfa, Dfb, and Dfd classes) the transversal classification does not provide added value. A legitimate question, therefore, is whether a finer, cross-cutting classification might produce groups of subclasses with easily identifiable distinctions, thereby reducing the number of classes whose subclasses are grouped together. We therefore investigated how the organization of our groups evolves as we reduce the distance threshold, thereby increasing the number of groups in accordance with the constraints imposed by our hierarchical classification.

The ninth group emerging from the hierarchical classification corresponds to a relatively limited surface area in the current configuration and concerns one of the four classes of the Csb class. These zones are characterized by a high average altitude (~2200 meters). Although surrounded by a predominantly equatorial climate, this climate is considered closer to Mediterranean in terms of the annual variations of its climatic parameters. However, it is perfectly acceptable to consider this group as reflecting a high-altitude equatorial climate. This clear and distinct contextual identification of this group compared to the other eight would therefore a priori justify the addition of this group.

The tenth emerging group is associated with a continental context. It separates from Group 4 (arid continental climate). Two subclasses are involved, belonging to the Bsh and Csc classes. Geographically, this group corresponds to central Australia. The criterion justifying the separation of this zone from the rest of the arid continental group is related to wind seasonality. In a context where the user aims to identify a past environmental setting for a paleo-study, we consider this level of detail to be excessive.

The eleventh group is an equatorial group derived from the seasonal equatorial group, and appears to separate areas subject to a marine influence. Specifically, this group includes South Asia, the southern and western parts of Central and South America previously identified as seasonal equatorial, as well as the coastal areas of two West African regions (Sierra Leone, Guinea, Gabon, Cameroon, and Nigeria). This group could be referred as an equatorial marine-influenced.

If we wish to explore the characteristics of groups that emerge beyond the range suggested by statistical criteria, several key steps can be identified. The first occurs when we consider 12 groups.

When choosing an initial number of 12 groups, we observe that the continental Sahara, together with all the continental areas of the Southern Hemisphere, forms a distinctly separate group from the continental clusters encompassing North America, Eastern Europe, and Central and Northern Asia.

Taking into considerations 16 groups allows the emergence of areas such as the Southwestern United States more clearly highlights the Chinese climate and also brings out the distinction between the continental climate of Central Asia and that of Eastern Europe and the central United States.

This illustrates how the hierarchical clustering approach can highlight intra-class climatic variability that would otherwise remain hidden in the standard Köppen-Geiger framework. Several options regarding the number of groups can be considered, and these will inherently depend on the purposes for which this cross-classification is proposed. As an example, we suggest that this cross-classification could be useful for assigning relative weights to climate scenarios in the context of a paleoenvironmental study, where a major challenge is the high level of uncertainty in the available climatic data. In this context, we consider it justified to adopt the smallest number of “meaningful” groups (here, 8).

The choice of a different number of groups remains fully adaptable based on the dendrogram. To facilitate this process, subclass-to-group assignment files corresponding to group numbers ranging from 8 to 20 are available through the Zenodo linked provided in the Data Availability Statement section.

## 7. Application of derivative Köppen-Geiger classes to hydro-environmental modeling of a lacustrine system

The processing of climatic data allows us to generate scenarios that realistically represent both annual means and intra-annual variability at the daily scale for six key climatic parameters.

The objective of this section is to illustrate a potential application of these scenarios in the simulation of a hydro-environmental system. In particular, we highlight two critical considerations:

* First, whether the climatic scenario should be complemented by a high-frequency signal to adequately capture short-term variability;

* Second, the overall degree of realism that can be expected from a numerical simulation based on such a climate scenario.

In this study, we chose to focus on the example of modeling the dynamics of a lake system. Lake water quality models can be categorized
<table><caption>**Table 2**
Correspondence between groups (1 to 8) and Köppen-Geiger climate classes (excluding EF).</caption><thead><tr><th>Climatic group</th><th>1 Cont altitude</th><th>2 Trop</th><th>3 Cont seasonal</th><th>4 Cont arid</th><th>5 Eq humid</th><th>6 Eq seasonal</th><th>7 Oc</th><th>8 Med</th></tr></thead><tbody><tr><td><strong>Tropical</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Af</td><td></td><td></td><td></td><td></td><td>2</td><td>1/3/4/5</td><td></td><td></td></tr><tr><td>Am</td><td></td><td>4/6</td><td></td><td></td><td>1/2</td><td>3/5</td><td></td><td></td></tr><tr><td>Aw</td><td></td><td>1/3/4</td><td></td><td></td><td>2/5</td><td></td><td></td><td></td></tr><tr><td><strong>Arid</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>BWh</td><td></td><td>3</td><td>2/5</td><td>1</td><td></td><td></td><td>4</td><td></td></tr><tr><td>BWk</td><td></td><td></td><td>1/4</td><td></td><td></td><td></td><td></td><td>2/3</td></tr><tr><td>BSh</td><td></td><td>2/3/4/5</td><td></td><td></td><td></td><td></td><td></td><td>1</td></tr><tr><td>BSk</td><td></td><td></td><td>1/2/3/4</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td><strong>Temperate</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Csa</td><td></td><td></td><td></td><td></td><td></td><td></td><td>2</td><td>1/3/4</td></tr><tr><td>Csb</td><td></td><td></td><td></td><td></td><td></td><td></td><td>2</td><td>1/3/4</td></tr><tr><td>Cwa</td><td>1</td><td>2/3/4</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Cwb</td><td></td><td>1/2/3/4</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Cfa</td><td>1</td><td></td><td>2/3/4/5</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Cfb</td><td></td><td>1</td><td>4</td><td></td><td></td><td></td><td>2/3</td><td></td></tr><tr><td>Cfc</td><td></td><td></td><td></td><td>1</td><td></td><td></td><td>3</td><td>2</td></tr><tr><td><strong>Cont</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Dsa</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td>1/2/3</td></tr><tr><td>Dsb</td><td></td><td></td><td></td><td></td><td></td><td></td><td>1</td><td>2/3</td></tr><tr><td>Dsc</td><td></td><td></td><td></td><td></td><td>1/2/4</td><td></td><td></td><td>3</td></tr><tr><td>Dsd</td><td></td><td></td><td></td><td></td><td>1/2</td><td></td><td></td><td></td></tr><tr><td>Dwa</td><td>2/3</td><td></td><td>1/4</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Dwb</td><td></td><td>3</td><td>1/2/4</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Dwc</td><td>2</td><td></td><td>3</td><td></td><td>1/4</td><td></td><td></td><td></td></tr><tr><td>Dwd</td><td></td><td></td><td></td><td></td><td>1/2/3/4</td><td></td><td></td><td></td></tr><tr><td>Dfa</td><td></td><td></td><td>1/2/3</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Dfb</td><td></td><td></td><td>1/2/3/4</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Dfc</td><td></td><td></td><td>4</td><td></td><td>1/2/3</td><td></td><td></td><td></td></tr><tr><td>Dfd</td><td></td><td></td><td></td><td></td><td>1/2/3</td><td></td><td></td><td></td></tr><tr><td><strong>Polar</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>ET</td><td>1</td><td></td><td></td><td>2/4</td><td></td><td></td><td>3</td><td></td></tr></tbody></table>

into (at least) three main types: lumped (0D), 1D, and 2D/3D models. The reader is referred to Piccolroaz et al. (2024) for a comprehensive review of the state of the art in these models. Nevertheless, we provide a brief summary to explain the reasons that led us to choose a lake model.

PCLake (Janssen et al., 2019) is a lumped (0D) model that treats the lake as a single homogeneous unit, focusing on biogeochemical processes like nutrient dynamics and eutrophication, making it suitable for shallow lakes where thermal stratification is minimal. 1D models such as FLAKE (Kirillin et al., 2011) and DYRESM (Imberger and Patterson, 1981) simulate vertical thermal stratification and water column processes, with FLAKE emphasizing thermal exchanges and DYRESM offering detailed seasonal stratification simulations. The General Lake Model (GLM) (Hipsey et al., 2019) also falls into the 1D category but is highly flexible, capable of simulating a broad range of physical and biogeochemical processes and adapted for simulating different lake sizes and conditions. Lastly, 2D/3D models like CE-QUAL-W2[4] and AEM3D (Hodges et al., 2000) provide more complex simulations, incorporating both horizontal and vertical dimensions to model hydrodynamics, nutrient transport, and biogeochemical cycles. These models are ideal for large or complex lakes where spatial variations are significant.

We opted to use a one-dimensional (1D) lake model, which offers a good balance between ease of implementation for a case study and the inclusion of realistic physical processes. More specifically, we selected the General Lake Model (GLM), in order to align our simulation framework with that of the study published by Ladwig et al. (2021).

That study focuses on Lake Mendota, located in the United States, which is typically associated with the Dfa climate class. In our climatic subclassification, this region corresponds to subclass 3. Our approach is as follows: we apply the climatic scenario that reproduces the median

curves of the corresponding subclass. High-frequency components are subsequently incorporated into specific parts of the climatic pattern to more accurately capture the variability observed in real-world climates. In this context, realism is assessed through the variations in gross primary productivity of organic matter generated by the lake system under climatic forcing (revealed by the variable GPP in GLM).

### 7.1. GLM-AED2 and the reference lake mendota study

The combination of GLM and its biochemical library AED2 allows managing an hydrodynamic-ecological model. GLM resolves the vertical thermal structure of lakes by dynamically adjusting water layers in response to changes in stratification, and accounts for key physical processes such as heat exchange, and water mixing by considering energy balances approaches. AED2 builds upon this physical foundation by simulating nutrient cycling, primary production, organic matter decomposition, and oxygen dynamics.

Details of the configuration used by Ladwig et al. to construct the Mendota case study is available in Ladwig et al. (2021). We sum-up here the principal characteristics.

The simulation spans a 37-year period, from 1979 to 2015. The lake water budget is determined through a mass-balance approach that accounts for river and groundwater inflows and outflows, as well as precipitation and evaporation. Inflow dynamics were evaluated using a process-based approach with the PIHM-Lake model. Calibration was achieved by directly comparing simulated streamflow with in-situ measurements provided by the United States Geological Survey (USGS). A similar approach was also deployed to estimate surface nutrient loadings.

The water quality library, AED2, was designed to simulate dissolve oxygen, silica, inorganic carbon, organic matter (refractory, particulate,
<figure><img src="image_11.png" /><figcaption>**Fig. 8.** Comparison between the default climatic scenarios of Dfa and the mean annual curves of climatic parameters recorded near Lake Mendota over a 37-year period. Dfa.3 corresponds to the geographically relevant subclass for Lake Mendota.</figcaption></figure>

and dissolved C, N and P), inorganic matter (refractory, particulate and dissolved C, N, and P), $PO_4$, $NO_3$, $NH_4$, and two phytoplanktons groups (diatoms and cyanobacteria).

Meteorological forcing was provided by the North American Land Data Assimilation System Phase 2 (NLDAS-2), offering hourly resolution to capture climatic variability at fine temporal scales.

Model calibration primarily focused on vertical profiles of temperature and dissolved oxygen, for which observational data are available throughout the entire simulation period.

The authors of this study also made available the parameterization files used for their GLM-AED2 simulation of Lake Mendota. These files serve as the baseline configuration for our own simulations, in which only the climatic forcing parameters are modified.

Our concern is to evaluate the capacity of a climatic scenario emerging from the most relevant subclasses to correctly reproduce the stratifications characteristics of the lake Mendota. From the Köppen-Geiger classification, it falls into the Dfa class. The geographical relevant subclass is the Dfa.3. Our first attempt is thus a comparison between the default climatic scenario of Dfa.3 and the mean annual curves of climatic parameters recorded near Lake Mendota over a 37-year period. As the climate scenario used by Ladwig et al. (2021) considered longwave radiation instead of cloud cover, we reconstructed the longwave signal based on the available climatic variables, specifically cloud cover, air temperature, and relative humidity. Results are shown in Fig. 8. Dfa.3 generally captures the overall shape of the mean curves at daily resolution. Given that high-frequency variability is absent from the synthetic scenario, we dedicate a section to evaluating its potential impact. The objective is to determine whether incorporating high-frequency components into the Dfa.3-based scenario is necessary to accurately reproduce the target dynamics. In this case, the dynamics are entirely driven by primary production of organic matter.

## 7.2. Process for obtaining stationarity

The first step of our workflow consists in establishing a steady-state regime, ensuring that the lake's dynamic response is driven primarily by the annual cycle described by the climate scenario, rather than by the initial state of the lake or the duration of the simulation. To

establish the reference simulation for Lake Mendota, we conducted a multi-year simulation of lake dynamics using the mean annual climate pattern derived from 37 years of meteorological observations (orange curves in the plots shown in Fig. 8). This approach is also extended to inflows, for which a single annualized composite inflow that mimic all the river influx contributions is generated.

Although such a stationary regime simulation is not available in the native version of the GLM code, this approach is not new in 1D lake modeling. For instance, it is implemented in the FLAKE model (Mironov et al., 2010). The long-term simulation, based on the repeated application of the same annual climate cycle, is run until a periodic state — i.e., stationarity — is achieved. In practice, fewer than 20 cycles are generally sufficient. Stationarity is considered to be reached when key lake state variables remain stable from one annual cycle to the next.

To this end, we developed a Python script (available in Zenodo linked provided in the Data Availability Statement section) that wraps the GLM executable and performs multiple annual cycles until a stationary state is reached. Each new cycle is initialized using the final state of the previous cycle. The criteria selected to validate stationarity are (1) the temperature profile and (2) the oxygen profile. The thresholds used to assess the stability of these profiles were empirically defined and tested across various synthetic lake configurations. We readily acknowledge that these criteria are somewhat rudimentary and could be refined. However, for the purposes of the present study, we consider them sufficiently robust and appropriate.

## 7.3. Impact of the daily cycle

As previously noted, the climate scenario derived from a given subclass does not capture daily-scale cyclicity. Nevertheless, several climatic variables, particularly shortwave radiation, exhibit pronounced diurnal variability. To reinforce the realism of the synthetic climatic scenarios, we apply a diurnal modulation to shortwave radiation data. The method involves calibrating a half-sinusoidal modulation aligned with the daylight interval, parameterized by the day of the year and geographic latitude.

The daylight duration *D* (photoperiod) for each day of the year *d* at latitude *φ* is computed based on standard astronomical formulas
defining sunrise and sunset times as follows (Forsythe et al., 1995):

$$t_{\text{sunrise}} = 12 - \frac{D}{2}, \quad t_{\text{sunset}} = 12 + \frac{D}{2}$$

Within this interval, the instantaneous shortwave radiation $SW(t)$ is modulated by a half-sinusoidal function:

$$SW_{\text{mod}}(t) = \begin{cases} SW(t) \times \sin\left(\pi \cdot \frac{t-t_{\text{sunrise}}}{D}\right), & \text{if } t_{\text{sunrise}} \le t \le t_{\text{sunset}} \\ 0, & \text{otherwise} \end{cases}$$

By imposing this physically based diurnal structure, the method enhances the temporal realism of shortwave radiation data, which is a priori crucial for modeling surface energy fluxes and other processes sensitive to sub-daily radiation variability.

Shortwave radiation is not the only climatic parameter exhibiting sub-daily variability. For variables such as precipitation, wind, or cloud cover, fluctuations are inherently stochastic and thus cannot be reliably captured by a simple deterministic signal. By contrast, for parameters such as air temperature or relative humidity, it may be more appropriate to introduce deterministic sub-daily components, as their temporal dynamics often follow more regular and predictable patterns.

In the specific context of this proof-of-concept experiment, we demonstrate that adding high-frequency components to shortwave radiation alone is sufficient to reproduce a realistic climatic dynamics. Therefore, we do not attempt here to develop parameterizations or empirical formulations to mimic sub-daily variations in the other variables.

## 7.4. Results

To assess the sensitivity of lake dynamics to different climatic forcings, six GLM simulations were performed for Lake Mendota using a common model configuration derived from Ladwig et al. (2021), with climate forcing as the only varying component. The reference simulation uses the mean annual climate derived from the original high-frequency meteorological dataset employed in that study (hereafter Original Data, HF), and serves as a benchmark for evaluating alternative scenarios. Five additional climate forcings were considered:

1. Dfa.1, HF: synthetic climate forcing from the first Dfa subclass, including a diurnal cycle in shortwave radiation;

2. Dfa.2, HF: synthetic forcing from the second Dfa subclass, also including diurnal shortwave variability;

3. Dfa.3, LF: synthetic forcing from the third Dfa subclass, without sub-daily variability;

4. Dfa.3, HF: same subclass as above, with an imposed diurnal shortwave cycle;

5. NASA POWER, HF: observational forcing extracted from the NASA POWER dataset at the Lake Mendota location (1991–2020), augmented with a diurnal shortwave cycle.

Fig. 9 summarizes the effects of these forcings on seasonal vertical profiles of temperature and dissolved oxygen, while Fig. 10 illustrates the annual evolution of daily mean gross primary productivity (GPP). Among the tested scenarios, Dfa.3, HF provides the closest agreement with the reference simulation. It accurately reproduces summer thermal stratification, vertical gradients in dissolved oxygen, and the seasonal evolution of GPP. In contrast, the low-frequency forcing Dfa.3, LF underestimates surface temperatures and produces noticeable deviations in summer temperature and oxygen profiles, despite yielding a GPP integral comparable to that of Dfa.3, HF. This highlights the limitations of purely annualized forcing in representing stratification and vertical mixing processes.

The Dfa.2, HF scenario leads to markedly different dynamics: the lake fails to stratify during summer, resulting in nearly uniform vertical

profiles of temperature and dissolved oxygen, and a strong overestimation of summer GPP. This behavior illustrates the high sensitivity of lake dynamics to intra-class climatic variability, even within a single Köppen–Geiger macro-climate. The Dfa.1, HF forcing produces results broadly similar to Dfa.3, HF, but with systematic underestimation of summer surface temperatures and overprediction of annual GPP, yielding the largest relative GPP error.

The simulation driven by NASA POWER, HF reproduces the reference GPP dynamics with the lowest annual error (+2.97%) and captures deep-water temperature and oxygen profiles well. However, it slightly overestimates near-surface temperature and dissolved oxygen during summer, particularly above the thermocline. These results indicate that NASA POWER constitutes a robust and operational alternative to high-resolution regional datasets for lake modeling, albeit with some limitations in surface-layer representation.

Overall, these results demonstrate that lake dynamics are jointly controlled by the choice of climatic subclass and the temporal resolution of radiative forcing. Selecting a geographically relevant subclass alone is insufficient if key high-frequency processes, such as the diurnal cycle of shortwave radiation, are not represented. Conversely, realistic temporal forcing cannot compensate for an inappropriate climatic subclass. Together, these two factors dominate the simulated thermal structure, oxygen dynamics, and primary production of the lake system.

## 8. Discussion and perspectives

This study demonstrates that a set of distinct climatic patterns can be identified, where each pattern is defined by the combination of six annualized climate curves corresponding to key parameters: shortwave solar radiation, air temperature, precipitation, wind speed, relative humidity, and cloud cover. The resulting climate classification, along with the construction of the associated climatic scenarios, has been specifically designed to be compatible with the types of climate data typically required for hydro-environmental modeling.

First of all, it is important to emphasize that the broader objective of improving the Köppen–Geiger classification is not new. Several studies have specifically attempted to introduce a more dynamic perspective into this system, in contrast to its often-criticized static nature. These efforts were primarily aimed at identifying climatic zones most vulnerable to desertification. For instance, Spinoni et al. (2015) proposed to reinforce desertification risk prediction by combining the Köppen–Geiger classification with the FAO aridity index (FAO (Food and Agriculture Organization of the United Nations) and UNESCO and WMO, 1977), defined as the ratio between precipitation and potential evapotranspiration. Based on the value ranges of this index, the FAO defined eight broad climatic macro-classes, allowing for a more nuanced characterization of arid and semi-arid regions. While the present study does not aim to trace climatic transitions over time, the subdivision of classes into internally coherent subclasses provides a framework that could, in principle, be applied to successive time periods. Such an approach could help detect shifts in subclass membership, which may be interpreted as early signals of climatic transitions—complementing, in a different way, previous efforts focused on desertification risk assessment.

In this paper, several steps of the proposed methodological framework for deriving Köppen–Geiger subclasses rely on specific assumptions and simplifications, which are examined and discussed in this section.

A legitimate first critique of this approach lies in the underlying assumption of climatic analogy between the present and the target simulation period. This implies that the climatic characteristics observed today can be meaningfully transferred to a different temporal context. However, in applications focusing on near-term projections — where climate changes are expected to be moderate — this analogy-based strategy may remain sufficiently robust. In such cases, the gradual
<figure><img src="image_12.png" /><figcaption>**Fig. 9.** Results from the GLM simulation for Lake Mendota. (a) Vertical profiles of temperature and dissolved oxygen over time using the annualized original forcing. (b) Legend for subfigures (c) and (d). (c) Mean vertical dissolved oxygen profiles in winter and summer. (d) Mean vertical temperature profiles in winter and summer.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 10.** Gross Primary Productivity (GPP) obtained from the GLM simulation for Lake Mendota using different climatic scenarios.</figcaption></figure>

evolution of climatic conditions could plausibly result in a neighbor-
ing region's Köppen-Geiger subclass expanding into the study area,
making the use of these derivative classes acceptable despite inherent
approximations.

It is particularly important to qualify the analogy with present-day
climate when applying these climatic classes to paleoclimatic contexts,
especially during periods when the global configuration of continents
— and consequently the fundamental drivers of climate — differed
markedly from today. Under such circumstances, it becomes conceiv-
able that certain climatic sub-classes, absent from the present-day
climate system, may have existed. Despite these likely biases, we argue
that the analogy-based approach offers genuine added value, primarily
because no robust alternative currently exists. Moreover, while the
use of the Köppen-Geiger climate classification for specific geological
periods has a long-standing history (e.g. Chaboureau (2012), Chu-
makov et al. (1995), Scotese (2000)), recent advances (Burgener et al.,

2023) have enabled the reconstruction of Köppen-Geiger paleomaps at
remarkably high resolutions, thereby reinforcing the relevance of our
approach in a paleoclimatic context.

Another limitation of the proposed methodology lies in the absence
of explicit diagnostic rules linking the derived subclasses to threshold-
based criteria grounded in the climate variables used in the clustering
process. In the original Köppen-Geiger classification, climatic zones
are defined based on clearly identified and interpretable thresholds—
typically annual or monthly means of temperature and precipitation. In
contrast, our approach relies on unsupervised clustering of a broader
set of variables (e.g., surface shortwave radiation, relative humidity,
wind speed), which enhances the descriptive power of the classification
but complicates the interpretability of the resulting subclasses. While it
would be highly valuable to identify threshold-based criteria (e.g., "if
SW < X W/m² and RH > Y%, then subclass A"), no such systematic
effort was undertaken in this study. As our preliminary explorations
revealed, such a task is far from straightforward: the boundaries be-
tween subclasses often emerge from complex, multivariate relationships
that cannot be easily reduced to simple thresholds. Addressing this
challenge would require a dedicated effort combining dimensionality
reduction, rule extraction methods, and validation against physical or
ecological expectations. We consider this a promising direction for
future work that would contribute to enhancing the interpretability and
operational usability of the derived climate subclasses.

A third key aspect to discuss is the stationarity assumption that is implicitly attached to the definition of climatic patterns. Although this stationarity condition must be examined for each simulation context, the assumption remains consistent, for instance, with lacustrine dynamics under climatic forcing. Several studies have highlighted how rapidly lake systems can adjust to abrupt climate shifts. While a detailed review is beyond the scope of this article, we refer the reader to the work of Woolway et al. (2017), who demonstrated that 20 Central European lakes exhibited significant increases in surface temperatures during the late 1980s, coinciding with an abrupt climatic shift. A similar study (Woolway and Merchant, 2017) conducted on 144 lakes worldwide, using 20 years of satellite data, further confirms that lake response times remain short. While lake depth influences the characteristic time required to reach thermal equilibrium, these adjustments typically occur over the course of a few months. In other words, the assumption of stationarity under a constant climatic forcing remains broadly valid for lake systems. This statement should not be confused with the numerous studies rightly showing that lake systems worldwide undergo year-to-year changes (Livingstone, 2008). These two perspectives can be reconciled by recognizing that lake systems adjust their theoretical stationary regime on an annual basis in response to climatic variations.

Steady-state conditions in a lake system are not solely a mat-
ter of thermal equilibrium but also require a balanced water budget
and a steady material supply, both of which are directly influenced
by watershed dynamics. While sediment (inorganic and organic) and
nutrient fluxes can respond rapidly to seasonal climatic variations,
vegetation cover dynamics may drive much longer transient responses,
often exceeding one year. This phenomenon is examined in studies
investigating tipping point mechanisms in vegetation cover dynamics
(e.g., Wendling et al. (2019)). Moreover, the common assumption of
constant water storage over time has been recently challenged; for
instance, Han et al. (2020) suggested that characteristic timescales for
achieving hydrological steady state can extend over several decades,
particularly in semi-arid regions.

Beyond these considerations, it is important to reflect on how such
studies should be designed. The deterministic nature of the processes
involved must be weighed against the inherent uncertainties associated
with climate scenarios. Consequently, assessing a simulated lacustrine
regime through this framework should rely on the interpretation of
multiple plausible climate scenarios.

Finally, it is essential to clarify the positioning of this approach with
respect to existing frameworks, in particular the ISIMIP project (Lange
and Büchner, 2021). In its most recent phase (ISIMIP3b), this project
provides future climate forcing derived from CMIP6 simulations, statis-
tically bias-corrected and developed under multiple emission scenarios,
with the aim of ensuring consistency and comparability at the global
scale. Our approach, by contrast, relies on the transposition of daily
climatic signatures representative of identified climatic regimes.

By forcing lake dynamics models, and more generally dynamic wa-
tershed models, with climatic regimes observed in analogue regions and
expected to establish locally in the future, the framework developed
here enables the exploration of the dynamic response of the simulated
system to regime transitions. It therefore adopts a mechanistic, process-
oriented perspective that complements the global climate projections
provided by intercomparison frameworks.

9. Conclusion

This study proposes a framework that extends the Köppen–Geiger
climate classification by linking each class to representative, multivari-
ate climatic patterns expressed as annualized daily curves. Rather than
redefining the classification itself, the approach leverages present-day
climate data to capture intra-class variability in a way that is directly
compatible with hydro-environmental modeling requirements.

A key outcome of this work is the demonstration that a single
Köppen–Geiger class may encompass multiple distinct climatic struc-
tures, each associated with different seasonal dynamics and inter-
variable relationships. Accounting for this intra-class variability through
data-driven subclasses proves essential when generating climatic sce-
narios intended to drive numerical models, as these differences can lead
to markedly different system responses.

The scenario-generation framework introduced here allows syn-
thetic climatic forcings to be constructed in a flexible yet statistically
coherent manner. By preserving seasonal structure and cross-variable
dependencies, and by allowing the imposition of user-defined annual
constraints, the method provides a practical alternative to histori-
cal or projected climate datasets in situations where such data are
unavailable, uncertain, or inappropriate.

The application to Lake Mendota illustrates the practical implica-
tions of these choices. Results show that selecting a geographically
relevant climatic subclass, combined with an appropriate represen-
tation of short-term radiative variability, is critical to reproducing
realistic lake dynamics. In particular, the inclusion of a diurnal cycle
in shortwave radiation emerges as a key factor controlling thermal
stratification and primary production, highlighting the sensitivity of en-
vironmental models to both subclass selection and temporal resolution
of climatic forcing.

Beyond this case study, the proposed transversal classification of subclasses into macro-environmental groups offers a pragmatic way to guide scenario selection in contexts where geographical or paleoenvironmental information is incomplete. This additional layer of interpretation is particularly relevant for long-term or paleoenvironmental studies, where uncertainty in climatic forcing is inherently high and multiple plausible scenarios must be considered.

Overall, this work provides a flexible and open framework for
generating climate scenarios tailored to hydro-environmental mod-
eling applications. Future developments will focus on extending its
application to a broader range of modern hydro-environmental sys-
tems, exploring multi-decadal future scenarios, and investigating the
dynamics of paleolakes under reconstructed climatic regimes.

CRediT authorship contribution statement

**Benoît Chauveau:** Writing – original draft, Validation, Supervision, Methodology, Data curation, Conceptualization. **Arnaud Pujol:** Writing – original draft, Validation, Supervision, Methodology, Data curation, Conceptualization. **Ahmed N. Mama:** Validation, Methodology.

Funding

The authors declare that no funds, grants, or other support were
received during the preparation of this manuscript

Declaration of competing interest

The authors declare that they have no known competing finan-
cial interests or personal relationships that could have appeared to
influence the work reported in this paper.

Acknowledgments

The authors are grateful to Javier Vidal Hurtado (EDF) for intro-
ducing them to GLM-AED2 modeling and to Laura Milelli (EDF) for her
careful reading of the manuscript.
# Appendix A. Supplementary data

Supplementary material related to this article can be found online at https://doi.org/10.1016/j.scitotenv.2026.181445.

## Data availability

All data used in this study are publicly available. A part of the data was obtained from National Aeronautics and Space Administration (NASA) Langley Research Center's Prediction Of Worldwide Energy Resources (POWER) project funded through the NASA Earth Science Division. Cloud cover data were obtained from the Moderate Resolution Imaging Spectroradiometer (MODIS) satellite sensor (Platnick et al., 2015). Additional cloud cover data were sourced from the ERA5 re-analysis dataset (Hersbach et al., 2023). The General Lake Model (GLM, version 3.3.0; freeware) was executed using a Python script. The script, all GLM input files, and the polynomial regressions, confidence intervals, and variance data required to reproduce the climatic scenarios for any of the climatic subclasses are archived in the associated Zenodo repository http://dx.doi.org/10.5281/zenodo.16419447.

## References

Abbott, M.B., Bathurst, J.C., Cunge, J.A., O'Connell, P.E., Rasmussen, J., 1986. An introduction to the European hydrological syste- système hydrologique européen, "SHE". J. Hydrol. 87 (1-2), 45-59. http://dx.doi.org/10.1016/0022-1694(86) 90114-9.

Arnold, J.G., Srinivasan, R., Muttiah, R.S., Williams, J.R., 1998. Large area hydrologic modeling and assessment part I: Model development. J. Am. Water Resour. Assoc. 34 (1), 73-89. http://dx.doi.org/10.1111/j.1752-1688.1998.tb05961.x.

Beck, C., Grieser, J., Kottek, M., Rubel, F., Rudolf, B., 2005. Characterizing global climate change by means of Köppen climate classification. In: Klimastatusbericht 2005. pp. 139-149, Available at: http://www.unet.univie.ac.at/~a0551135/pdf/ BeckCo2006.pdf.

Beck, H.E., Zimmermann, N.E., McVicar, T.R., Vergopolan, N., Berg, A., Wood, E.F., 2018. Present and future Köppen-Geiger climate classification maps at 1-km resolution. Sci. Data 5, 180214. http://dx.doi.org/10.1038/sdata.2018.214, URL https://www.nature.com/articles/sdata2018214.

Beck, H.E., Zimmermann, N.E., McVicar, T.R., Vergopolan, N., Berg, A., Wood, E.F., 2023. High-resolution (1 km) Köppen-Geiger climate classification maps for 1901-2099. Sci. Data 10, 1-12. http://dx.doi.org/10.1038/s41597-023-02549-6, URL https://www.nature.com/articles/s41597-023-02549-6.

Belda, M., Holtanová, E., Halenka, T., Kalvová, J., 2014. Climate classification revisited: from Köppen to trewartha. Clim. Res. 59 (1), 1-13. http://dx.doi.org/10.3354/ cr01204.

Boorman, D.B., Hollis, J.M., Lilly, A., 2007. Hydrochemical modelling of integrated surface-subsurface water flow using SHETRAN. J. Hydrol. 234, 101-114. http://dx.doi.org/10.1016/S0022-1694(00)00217-5.

Brutsaert, W., 1975. On a derivable formula for long-wave radiation from clear skies. Water Resour. Res. 11 (5), 742-744. http://dx.doi.org/10.1029/ WR011i005p00742.

Burgener, L., Tabor, N.J., Meyers, S.R., Steier, A., Carroll, A.R., Price, G.D., 2023. Cretaceous climates: Mapping paleo-Köppen climatic zones using a Bayesian statistical analysis of lithologic, paleontologic, and geochemical proxies. Palaeogeogr. Palaeoclimatol. Palaeoecol. 609, 111276. http://dx.doi.org/10. 1016/j.palaeo.2022.111276, URL https://www.sciencedirect.com/science/article/ pii/S0031018222005442.

Chaboureau, A.-C., 2012. Impact du climat et de la tectonique sur la dynamique des systèmes sédimentaires pendant l'ouverture de l'Atlantique Sud (Ph.D. thesis). Université Rennes 1, URL https://theses.hal.science/tel-00830475, Mémoire de Géosciences-Rennes n°147. ISBN : 2-914375-89-1.

Chumakov, N.M., Zharkov, M.A., Herman, A.B., Doludenko, M.P., Kalandadze, N.N., Lebedev, E.L., Ponomarenko, A.G., Rautian, A.S., 1995. Climatic belts of the mid-cretaceous time. Strat. Geol. Correl. 3 (3), 241-260.

Crawford, T.M., Duchon, C.E., 1999. An improved parameterization for estimating effective atmospheric emissivity for use in calculating daytime downwelling long-wave radiation. J. Appl. Meteorol. 38, 474-480. http://dx.doi.org/10.1175/1520- 0450(1999)038<0474:AIPFEE>2.0.CO;2.

DHI, 2012. MIKE SHE User Manual. Technical Report, DHI Water & Environment, Denmark, URL https://www.mikepoweredbydhi.com/.

Ewen, J., Parkin, G., O'Connell, P.E., 2000. SHETRAN: Distributed river basin flow and transport modeling system. J. Hydrol. Eng. 5 (3), 250-258. http://dx.doi.org/10. 1061/(ASCE)1084-0699(2000)5:3(250).

FAO (Food and Agriculture Organization of the United Nations) and UNESCO and WMO, 1977. World Map of Desertification. Food and Agriculture Organization (FAO), Rome.

Forsythe, W.C., Rykiel, E.J., Stahl, R.S., Wu, H.I., Schoolfield, R.M., 1995. A model comparison for daylength as a function of latitude and day of year. Ecol. Model. 80 (1), 87-95. http://dx.doi.org/10.1016/0304-3800(94)00034-3.

Hamman, J.J., Nijssen, B., Bohn, T.J., Gergel, D.R., Mao, Y., 2018. The variable infiltration capacity model version 5 (VIC-5): infrastructure improvements for new applications and reproducibility. Geosci. Model. Dev. 11, 3481-3496. http://dx.doi.org/10.5194/gmd-11-3481-2018.

Han, J., Yang, Y., Roderick, M.L., McVicar, T.R., Yang, D., Zhang, S., Beck, H.E., 2020. Assessing the steady-state assumption in water balance calculation across global catchments. Water Resour. Res. 56, e2020WR027392. http://dx.doi.org/10.1029/ 2020WR027392.

Hersbach, H., Bell, B., Berrisford, P., Biavati, G., Horányi, A., Muñoz Sabater, J., Nicolas, J., Peubey, C., Radu, R., Rozum, I., Schepers, D., Simmons, A., Soci, C., Dee, D., Thépaut, J.-N., 2023. ERA5 hourly data on single levels from 1940 to present. http://dx.doi.org/10.24381/cds.adbb2d47, Copernicus Climate Change Service (C3S) Climate Data Store (CDS). (Accessed on 01 June 2025).

Hipsey, M.R., Bruce, L.C., Boon, C., Busch, B., Carey, C.C., Hamilton, D.P., Hanson, P.C., Read, J.S., de Sousa, E., Weber, M., Winslow, L.A., 2019. A general lake model (GLM 3.0) for linking with high-frequency sensor data from the global lake ecological observatory network (GLEON). Geosci. Model. Dev. 12 (1), 473-523. http://dx.doi.org/10.5194/gmd-12-473-2019, URL https://gmd.copernicus.org/articles/ 12/473/2019/.

Hodges, B.R., Imberger, J., Saggio, A., Winters, K.B., 2000. Modeling basin-scale internal waves in a stratified lake. Limnol. Oceanogr. 45 (7), 1603-1620. http://dx.doi.org/10.4319/lo.2000.45.7.1603.

Holdridge, L.R., 1967. Life Zone Ecology, rev. ed. Tropical Science Center, San José, Costa Rica.

Huth, R., Beck, C., Philipp, A., Demuzere, M., Ustrnul, Z., Cahynová, M., Kyselý, J., Veito, O.E., 2008. Classifications of atmospheric circulation patterns: Recent advances and applications. Ann. New York Acad. Sci. 1146 (1), 105-152. http://dx.doi.org/10.1196/annals.1446.019.

Imberger, J., Patterson, J.C., 1981. A dynamic reservoir simulation model - DYRESM. In: Fischer, H.B. (Ed.), Transport Models for Inland and Coastal Waters. Academic Press, pp. 310-361.

Janssen, A.B.G., Teurlinx, S., Beusen, A.H.W., Huijbregts, M.A.J., Rost, J., Schipper, A.M., Seelen, L.M.S., Mooij, W.M., Janse, J.H., 2019. PCLake+: A process-based ecological model to assess the trophic state of stratified and non-stratified freshwater lakes worldwide. Ecol. Model. 396, 23-32. http://dx.doi.org/10.1016/j. ecolmodel.2019.01.006.

Kirillin, G., Hochschild, J., Mironov, D., Terzhevik, A., Golosov, S., Nützmann, G., 2011. FLake-Global: On line lake model with worldwide coverage. Environ. Model. Softw. 26, 683-684. http://dx.doi.org/10.1016/j.envsoft.2010.12.004.

Köppen, W., 1900. Versuch einer klassifikation der klimte, vorzugsweise nach ihren beziehungen zur pflanzenwelt. Geogr. Z. 6 (12), 657-679.

Köppen, W., 1936. Das geographische system der klimte. In: Köppen, W., Geiger, R. (Eds.), Handbuch der Klimatologie. Gebrüder Rorntraeger, Berlin, pp. 1-44.

Kottek, M., Grieser, J., Beck, C., Rudolf, B., Rubel, F., 2006. World map of the Köppen-Geiger climate classification updated. Meteorol. Z. 15 (3), 259-263. http://dx.doi.org/10.1127/0941-2948/2006/0130.

Ladwig, R., Hanson, P.C., Dugan, H.A., Carey, C.C., Zhang, Y., Shu, L., Duffy, C.J., Cobourn, K.M., 2021. Lake thermal structure drives interannual variability in summer anoxia dynamics in a eutrophic lake over 37 years. Hydrol. Earth Syst. Sci. 25, 1009-1032. http://dx.doi.org/10.5194/hess-25-1009-2021.

Lange, S., Büchner, M., 2021. ISIMIP3b bias-adjusted atmospheric climate input data (v1.1). URL https://coilink.org/20.500.12592/qggxv5. (Accessed 16 December 2025).

Liang, X., Lettenmaier, D.P., Wood, E.F., Burges, S.J., 1994. A simple hydrologically based model of land surface water and energy fluxes for general circulation models. J. Geophys. Res. 99 (D7), 14415-14428. http://dx.doi.org/10.1029/94JD00483.

Livingstone, D.M., 2008. A change of climate provokes a change of paradigm: taking leave of two tacit assumptions about physical lake forcing. Int. Rev. Hydrobiol. 93 (4-5), 404-414. http://dx.doi.org/10.1002/iroh.200811061.

Mironov, D., Heise, E., Kourzeneva, E., Ritter, B., Schneider, N., Terzhevik, A., 2010. Implementation of the lake parameterisation scheme flake into the numerical weather prediction model COSMO. Boreal Environ. Res. 15, 218-230, URL https://www.sciencedirect.com/science/article/abs/pii/S1364815210003245.

NASA POWER Project, 2025. NASA prediction of worldwide energy resources (POWER) project - API access. Data accessed on March 2025 via the POWER API, https://power.larc.nasa.gov/api/.

Neitsch, S.L., Arnold, J.G., Kiniry, J.R., Williams, J.R., 2011. Soil and Water Assessment Tool Theoretical Documentation, Version 2009. Technical Report, Texas Water Resources Institute, Texas A&M University System, URL https://swat.tamu.edu/ media/1130/SWAT2009theory.pdf.
Piccolroaz, S., Zhu, S., Ladwig, R., Carrea, L., Oliver, S., Piotrowski, A.P., Ptak, M., Shinohara, R., Sojka, M., Woolway, R.I., Zhu, D., 2024. Lake water temperature modeling in an era of climate change: data sources, models, and future prospects. Rev. Geophys. 62 (1), e2023RG000816. http://dx.doi.org/10.1029/2023RG000816, URL https://centaur.reading.ac.uk/115202/.

Platnick, S., Ackerman, S.A., King, M.D., Meyer, K., Menzel, W.P., Holz, R.E., Baum, B.A., Yang, P., 2015. MODIS atmosphere L2 cloud product (06L2). http://dx.doi.org/10.5067/MODIS/MOD06_L2.061.

Scotese, C.R., 2000. The paleomap project. URL http://www.scotese.com.

Spinoni, J., Vogt, J., Naumann, G., Carrao, H., Barbosa, P., 2015. Towards identifying areas at climatological risk of desertification using the Köppen-Geiger classification and FAO aridity index. Int. J. Climatol. 35 (9), 2210–2222. http://dx.doi.org/10.1002/joc.4084.

Thornthwaite, C.W., 1948. An approach toward a rational classification of climate. Geogr. Rev. 38 (1), 55–94. http://dx.doi.org/10.2307/210739.

Trewartha, G.T., Horn, L.H., 1980. Introduction to Climate, fifth ed. McGraw-Hill, New York, NY.

Wendling, V., Peugeot, C., Mayor, Á., Hiernaux, P., Mougin, É., Grippa, M., Kergoat, L., Walcker, R., Galle, S., Lebel, T., 2019. Drought-induced regime shift and resilience of a sahelian ecohydrosystem. Environ. Res. Lett. 14 (10), 105005. http://dx.doi.org/10.1088/1748-9326/ab3dde, URL https://hal.science/hal-02344434.

Woolway, R.I., Dokulil, M.T., Marszelewski, W., Schmid, M., Bouffard, D., Merchant, C.J., 2017. Warming of central European lakes and their response to the 1980s climate regime shift. Clim. Change 142 (3–4), 505–520. http://dx.doi.org/10.1007/s10584-017-1966-4.

Woolway, R.I., Merchant, C.J., 2017. Amplified surface temperature response of cold, deep lakes to inter-annual air temperature variability. Sci. Rep. 7, 4130. http://dx.doi.org/10.1038/s41598-017-04058-0.