Reversed asymmetric warming of sub-diurnal temperature over land during recent decades

Domain: earth_science
Published: 2023-11-08
Editorial outcome: major_revision

Abstract
In the latter half of the twentieth century, a significant climate phenomenon “diurnal asymmetric warming” emerged, wherein global land surface temperatures increased more rapidly during the night than during the day. However, recent episodes of global brightening and regional droughts and heatwaves have brought notable alterations to this asymmetric warming trend. Here, we re-evaluate sub-diurnal temperature patterns, revealing a substantial increase in the warming rates of daily maximum temperatures (Tmax), while daily minimum temperatures have remained relatively stable. This shift has resulted in a reversal of the diurnal warming trend, expanding the diurnal temperature range over recent decades. The intensified Tmax warming is attributed to a widespread reduction in cloud cover, which has led to increased solar irradiance at the surface. Our findings underscore the urgent need for enhanced scrutiny of recent temperature trends and their implications for the wider earth system.

Introduction
The surface air temperature (SAT) is a commonly used measure of surface climate change due to its ability to represent terrestrial energy exchange with reasonable accuracy1,2. In addition to daily average temperatures, the diurnal temperature range (DTR), defined as the difference between the daily maximum temperature (Tmax) and daily minimum temperature (Tmin), provides useful information about the climate3,4. Changes in DTR have received considerable attention because it is closely linked with crop yields5–7, plant growth8–11, animal wellbeing12,13 and human health14–16. Existing studies found that the land surface warming since the 1950s has been associated with larger increases in Tmin than in Tmax, i.e., decreases in DTR, which is commonly known as night warming or asymmetric warming17–20.

Changes in DTR are complex; the changes are subject to many factors, including cloud cover21–23, solar radiation24,25, aerosols26, precipitation27,28, planetary boundary layer height29, land use change30,31 and deforestation32,33. For example, an increase in total cloud cover reduces DTR due to a decrease in daytime surface insolation and an increase in night-time downward longwave radiation34. Increases in precipitation and soil moisture can reduce Tmax and therefore DTR through increased evaporative cooling35. However, considering known changes in these processes, the higher warming rate of Tmin than that of Tmax is seemingly inconsistent with two recent phenomena. One contradiction is the increase in surface solar radiation that has occurred since the late 1980s, which is referred to as brightening after dimming36–38. Solar radiation affects Tmax more than Tmin26, thus, brightening should contribute to further warming of Tmax. The second contradiction is the increased occurrence of drought events and heatwaves, especially in spring and summer39–41. This phenomenon indicates that the cooling effect of soil moisture may have been weakened, which should lead to a faster increase in Tmax. Motivated by these different changes, we re-evaluated the warming rates of Tmax and Tmin over the period of 1961–2020 and investigated the causes behind the associated changes in DTR.

Results
Observed reversing asymmetric warming
With changes in recent years included, we evaluated the surface warming rates of Tmax and Tmin using two station observation-based datasets from Berkeley Earth Surface Temperatures (BEST)42 and the Climatic Research Unit Time-Series version 4.07 (CRU TS)43. Both datasets were gridded and gap-filled over the land masses, and BEST was used to detect global surface temperature changes in the Intergovernmental Panel on Climate Change (IPCC) sixth assessment report44. As shown in Fig. 1a, b, the trends in the global area-weighted average of Tmax calculated using a 30-year moving window increased faster than that of Tmin during 1961–2020. The warming rate in global average Tmax reached the warming rate in Tmin in recent decades, with an earlier surpassing moment in BEST than in CRU TS. In the last 30-year window (1991–2020), both datasets exhibit a slight (CRU TS dataset) or substantially more pronounced (BEST dataset) rise in the global average of Tmax compared to that of Tmin. Spatially, a stronger warming rate of Tmax was found in almost 25% of the land area in the earlier time window (1961–1990); the area expanded rapidly, and in the recent time window (1991–2020) the area of stronger warming rate of Tmax was found in at least half of the area (approximately 52% in CRU TS and 70% in BEST, Fig. 1e, f, with consistent results between the two datasets in 61% of the total land area). Observations from stations based on Global Surface Summary of the Day (GSOD) support the finding of a greater increase in Tmax than in Tmin from 1991–2020 (Fig. 2). Analyzes of all observed monitoring sites indicate that 63% of the sites exhibited an overall upward trend in DTR, while more than one-third (35%) of all sites displayed a statistically significant increase in DTR during 1991–2020. Spatially, a widespread decreasing trend in DTR was detected (81% and 76% land area fraction in the BEST and CRU TS temperature datasets, respectively) for the period 1961–1990, with the exceptions of southern Africa, Southern Europe, and some regions in Northern America (Fig. 1c, d). However, during the recent decades of 1991–2020, both temperature datasets show a consistent increase in DTR over more than half of the land area, in particular over the western United States, southern Europe, West Africa, inner East Asia, and Australia.

Our analyzes reveal that the reversal of asymmetric warming was more pronounced in the BEST dataset compared with CRU TS. We assessed the accuracy of these two sets of gridded DTR data by comparing them against GSOD station DTR data (see Methods). The findings revealed a significantly stronger correlation between BEST’s DTR data and the GSOD station dataset (one-tailed t-test, p < 0.001) compared to the correlation between CRU TS’s DTR data and the GSOD station dataset (Supplementary Fig. 1). Thus, in the subsequent analysis we used the BEST dataset to study annual and seasonal trends in Tmax, Tmin and DTR for the periods 1961–1990 and 1991–2020.

A significant decline in the global average of DTR (−0.08°C decade−1, p < 0.05) was found during 1961–1990. For the period 1991–2020, the globally averaged DTR increased at a rate of 0.06°C decade−1 (p < 0.05) (Supplementary Fig. 2c). This reversal in DTR occurred mainly due to a marked increase in Tmax for the latter period (0.35°C year−1 decade−1, p < 0.05 during 1991–2020 vs. 0.13°C year−1 decade−1, p < 0.05 during 1961–1990, Supplementary Fig. 2a, b). Seasonally, the largest decline in DTR was detected in winter during 1961–1990, which is consistent with previous findings35. Significant increases (p < 0.05) in DTR were detected in spring, summer and winter, which were mainly caused by a significant increase in Tmax (Supplementary Fig. 2a). All these evidences point to a reversing asymmetric warming over land in recent decades.

Potential mechanisms behind reversing asymmetric warming
Earlier studies suggested that the recent global warming is mainly forced by greenhouse gases45, while changes in DTR are influenced largely by cloud cover23,46. In addition, aerosols and soil moisture may also have affected the variability of DTR47,48. These potential drivers of DTR exhibit a high degree of correlation. For instance, the phenomenon of aerosol-cloud interactions, such as the influence of aerosols on cloud albedo49 and lifetime50, has gained growing attention51–53. Moreover, clear associations can be observed between cloud cover and soil moisture, as increasing cloud cover is associated with enhanced precipitation54, which subsequently leads to soil wetting. The intricate interactions among these factors highlight the presence of high multi-collinearity when performing regression analysis with them as independent variables, posing challenges in discerning the dominant drivers of DTR change. To address the impact of multi-collinearity and accurately identify relationships, we conducted ridge regression analyzes55 at individual grid points. This analysis (Methods) employed as independent variables monthly total cloud cover from the fifth-generation ECMWF reanalysis (ERA5)56 dataset, aerosol optical depth from Modern-Era Retrospective analysis for Research and Applications, version 2 (MERRA-2)57, and soil moisture from the Global Land Evaporation Amsterdam Model (GLEAM)58 dataset, while DTR was calculated from the BEST dataset as the dependent variable spanning 1981 to 2020. Ridge regression is a linear regression method with regularization that can effectively address the issue of multi-collinearity. By introducing a penalty term in the cost function that discourages overly large parameter values, it improves upon the ordinary least squares regression model in scenarios where there exist strong correlations among independent variables. An evaluation of the ridge regression (Methods) shows that the ridge regression model can capture the majority of the explained variance in DTR (Supplementary Fig. 3a), except in certain regions in particularly Africa, South America, and Northern Hemisphere high-latitudes.

The analysis using ridge regression suggests a worldwide negative response of DTR to changes in total cloud cover and soil moisture (Supplementary Fig. 3b, d), with generally stronger negative impacts of cloud cover than soil moisture. This result aligns with previous findings21,23. Negative responses of DTR to aerosol concentrations are detected in Western Africa, the Arabian Peninsula, India and southern China (Supplementary Fig. 3c). To make a quantitative comparison, Fig. 3 shows a cyan-magenta-yellow (CMY) composite map (see Methods) of the relative contributions of total cloud cover, aerosol and soil moisture to DTR. Generally, changes in DTR during 1981–2020 were dominated by cloud cover variations over 83.5% of the global land area; this phenomenon is widely detected in most of North America, southern South America, Europe, southern Africa, Central Asia, and East Asia. Moreover, the dominant effect of cloud cover on DTR over land is further confirmed when using an independent total cloud cover dataset obtained from the Moderate Resolution Imaging Spectroradiometer (MODIS) satellite during 2003–2020 (84.4% of the total land area; Supplementary Fig. 4).

The ridge regression analysis assumes linear relationships and may overlook the non-linear associations among cloud, aerosol, soil moisture, and DTR. To address this limitation, we repeated the regression analysis using the Random Forest algorithm59 (Methods). Random Forest is a machine learning technique based on decision trees and can capture non-linear relationships. The results from the Random Forest regression consistently corroborate our prior findings, indicating that total cloud cover is indeed the dominant driving force behind DTR fluctuations in terrestrial regions, encompassing a substantial proportion (79.1%) of the land area (Supplementary Fig. 5). DTR is significantly negatively correlated with total cloud cover in most regions of the terrestrial land surface, suggesting that DTR is greatly influenced by cloud cover.

To further elucidate the influence of different environment variables on DTR, we analyzed the changes of DTR, total cloud cover, aerosol optical depth and soil moisture over the global land area during recent decades (Supplementary Fig. 6). During the 1960s to the early 1970s, there was a strong and significant upward trend in total cloud cover, which later stabilized. However, a remarkable downward trend emerged from the mid-1980s onwards. This decline in cloud cover corresponded with an increase in DTR during the same period. Between 1981 and 2020, a significant (p < 0.05) negative correlation was found between the annual average global surface total cloud cover and DTR. The mean bootstrapped partial correlation coefficient (−0.47; −0.54 to −0.41, 95% confidence interval (CI)) between these two variables was higher than that of the mean correlation between DTR and soil moisture (−0.31; −0.39 to −0.23, 95% CI), as well as the correlation between DTR and aerosol optical depth (−0.18; −0.28 to −0.07, 95% CI).

Cloud cover can influence DTR through two primary mechanisms. On the one hand, diminished cloud cover leads to greater daytime solar radiation, resulting in elevated daytime maximum temperatures and, consequently, an expanded DTR. Conversely, reduced cloud cover leads to decreased night-time outgoing radiation at the surface, causing lower minimum temperatures and further increasing DTR23,46. To discern the predominate effect, we calculated the partial correlations between both cloud cover and Tmax and cloud cover and Tmin (Supplementary Fig. 7). The results indicate that, overall, the prevalent association is a negative correlation between cloud cover and maximum temperatures. This suggests that decreased cloud cover has led to an increased DTR in recent decades, primarily due to the associated increase in incoming solar radiation, while the impact of night-time cooling has played a secondary role.

To reveal the spatial patterns of cloud and radiation changes, we compared the trends in total cloud cover and solar radiation during the earlier decades (1961–1990) and the recent decades (1991–2020) (Fig. 4). Specifically, a global (over 69% of land) increase in cloud cover corresponded to a decline in solar radiation at the surface over almost three-quarters (74%) of the land area for the period 1961–1990; the dimming has disappeared since the 1990s, and increases in solar radiation were detected in over 63% of global land, which is partly due to an extensive (over 66% of land) decrease in cloud cover. Over the period 1961–1990, central North America, southern South America, the Mediterranean, and Australia saw weak decreases in solar radiation, with the local DTR narrowing or experiencing little change. However, during 1991–2020, a substantial increasing trend in solar incident radiation was detected over central North America, southern South America the Mediterranean, and Australia, corresponding to a local general increase of DTR.

To further elucidate the relationship between incident solar radiation and DTR, we conducted a partial correlation analysis between DTR and incident shortwave radiation extracted from ERA5 (1961–2020), CERES (2001–2020), and MERRA-2 (1981–2020), while adjusting for the influence of soil moisture (Supplementary Fig. 8). The results show a significant positive correlation between DTR and surface solar radiation over most land areas of the world in all four datasets, confirming a credible positive influence of solar radiation on DTR, consistent with previous studies24,25. The above results indicate that the recent global decrease in cloud cover60–62 has increased incoming solar radiation at the surface, resulting in a greater increase in Tmax than in Tmin, ultimately expanding DTR.

Methods
Data
The Tmax and Tmin data utilized in this study were obtained from the Climatic Research Unit Time-Series version 4.07 (CRU TS)43 and Berkeley Earth Surface Temperatures (BEST) datasets42. The CRU TS temperature dataset is derived from a blend of data sources including weather station records, ship logs, and more recent satellite observations. This dataset undergoes meticulous calibration to account for biases and variations in measurement methods, making it an extensively used resource in climate research for documenting long-term temperature trends and variations. It provides spatiotemporal resolutions of 0.5° × 0.5° on a monthly basis, covering the period from 1901 to 2022.

The BEST dataset encompasses a larger sample of approximately 39,000 records. It employs advanced statistical techniques to quantify and adjust measurement biases, ensuring a highly accurate representation of global temperature trends. This dataset offers spatiotemporal resolutions of 1° × 1° on a monthly basis, spanning from 1850 to the present.

In situ Tmax and Tmin records were acquired from the Global Surface Summary of the Day (GSOD) database, downloaded in September 2021 from https://data.nodc.noaa.gov/cgi-bin/iso?id=gov.noaa.ncdc:C00516. The GSOD dataset originates from the Integrated Surface Hourly (ISH) database developed by the US National Climatic Data Center. This dataset contains meteorological variables, including temperature, precipitation, and wind speed, from over 9000 weather stations. The original observations undergo rigorous quality control procedures to ensure their accuracy84. Specifically, we selected daily maximum and minimum temperature data with relatively complete records from all stations covering the period 1991 to 2020. We employed stringent selection criteria to exclude incomplete data series, only including stations with missing values not exceeding one year of the analysis period (1991 to 2020) and with complete records for all 12 months. Furthermore, we required each monthly value to be derived from at least 15 days of data. Consequently, our trend analysis encompassed a total of 2557 stations.

Monthly data on total cloud cover were obtained from the fifth-generation ECMWF reanalysis (ERA5) dataset56 on a 0.25° × 0.25° regular latitude–longitude grid (the native resolution of ERA5 is about 31 km). In addition, monthly total cloud cover was obtained from the MODIS MCD06 COSP data at a spatial resolution of 1° × 1° after 2003. ERA5 total cloud cover was used for the main analyzes.

The monthly incident shortwave radiation data in all-sky conditions were obtained from the ERA5 dataset on a 0.25° grid, the Cloud and the Earth’s Radiant Energy System energy balanced and filled edition 4.1 (CERES–EBAF)85 dataset at a spatial resolution of 1° × 1° after March 2000, and the Modern-Era Retrospective analysis for Research and Applications, version 2 (MERRA-2)57 dataset with a spatial resolution of 0.625° × 0.5° after 1980.

The monthly average aerosol optical depth and surface albedo at a spatial resolution of 0.625° × 0.5° was obtained from the MERRA-2 dataset beginning in 1980. The monthly surface soil moisture at a spatial resolution of 0.25° was obtained from the Global Land Evaporation Amsterdam Model (GLEAM) version 3.5 dataset58, which is a global dataset spanning 40 years from 1981 to 2020 and based on satellite and reanalysis data. The cloud cover, radiation, aerosol optical depth, and soil moisture data sets were aggregated to a common 0.5° grid.

Assessment of temperature-gridded data
We assessed the accuracy of the temperature gridded data using GSOD station temperature dataset. The evaluation period spanned from 1978 to 2020, during which the GSOD data records were relatively complete. Following a similar methodology used in analyzing station-based DTR trends, we selected a total of 2,058 stations with missing observations comprising less than 5% of the evaluation period for temperature-gridded data assessment. In the evaluation process, we extracted the time series of DTR from the gridded BEST and CRU TS datasets corresponding to each station’s location. Daily data were averaged to obtain yearly data and missing values at these stations were filled using linear interpolation. We then computed Pearson correlation coefficients between the annual station observation series and the observed DTR sequence at each station, serving as an indicator of accuracy. The significance of the correlation coefficient differences between the DTR station data and the two sets of DTR gridded data was tested by the one-tailed Student’s t-test.

Seasonal analysis
To analyze DTR variations by season, we defined the seasons as follows: March, April, May for spring (autumn), June, July, August for summer (winter), September, October, November for autumn (spring), and December, January, February for winter (summer) in the Northern Hemisphere (Southern Hemisphere). Correspondingly, the annual average value of a variable in one year is defined as the average of 12 months from December in the preceding year to November of that year.

Ridge regression
Ridge regression is a method for estimating coefficients in multiple linear regression models in scenarios characterized by high correlations among the independent variables. In ridge regression, a regularization term is introduced to the standard least squares objective function, aiding in stabilizing the estimated coefficients. This regularization term is controlled by a tuning parameter denoted as λ. The ridge regression objective function can be expressed as follows:

\hat{\beta} = \sum_{i=1}^{n} (y_i - \beta_0 - \sum \beta_i x_i)^2 + \lambda \sum \beta_i^2 (1)

where \hat{\beta} represents the estimated regression coefficients, y is the dependent variable, \beta_0 is the intercept term, and \beta_i signifies the regression coefficient for the independent variable x_i. The tuning parameter \lambda governs the extent of shrinkage applied to the coefficients. When \lambda is set to zero, the regularization term exerts no effect, and ridge regression reduces to ordinary least squares regression. However, as \lambda increases, the penalty term gains influence. The larger \lambda becomes, the more pronounced is the shrinkage applied to the coefficients. Consequently, ridge regression effectively mitigates the impact of multi-collinearity.

Prior to performing the ridge regression analysis in this study, all time series were detrended by subtracting the linear trend and transformed into z-scores by subtracting the monthly climatology means dividing by the monthly climatological standard deviations from 1981 to 2020. For validation, the z-scores were randomly divided into an 80% calibration dataset and a 20% validation dataset. The dataset in the training set was used to train the ridge regression model, while the dataset in the validation set was employed to assess the performance of the ridge regression model. Afterwards, all the datasets were merged to determine the ridge regression coefficients.

The ridge regression was performed on the variables in each grid cell, and the tuning parameter \lambda was determined for each grid cell based on the Variance Inflation Factor (VIF) of the independent variables within each individual grid cell. VIF serves as a measure of multi-collinearity among the independent variables in the regression model. The following formula was used to calculate the VIF:

VIF_i = \frac{1}{1 - R_i^2} (2)

where R_i^2 denotes the coefficient of determination between the ith independent variable and all other independent variables. A higher VIF value indicates stronger multi-collinearity within the regression model’s independent variables. Here, a VIF value less than 3 suggests an acceptable level of multi-collinearity86.

Throughout the regression analysis process, the initial value of \lambda was set to 0 and incremented by a step size of 0.01. As \lambda increased, the degree of multi-collinearity decreased, subsequently resulting in a decline in the VIF value. The incrementation of \lambda ceased when the VIF value dropped below 3, with this \lambda value being determined as the tuning parameter at this grid point. The significance of the ridge regression analysis was assessed utilizing an F-test at a significance level of 0.05. The accuracy of the ridge regression model was evaluated by employing the coefficient of determination on the validation dataset.

CMY composite
The color of the cyan-magenta-yellow (CMY) composite was determined by the relative contributions from the magnitudes of the ridge regression coefficients (R.c), which can be expressed as follows:

C = \frac{|R.c_{x1}|}{|R.c_{x1}| + |R.c_{x2}| + |R.c_{x3}|} (3)

M = \frac{|R.c_{x2}|}{|R.c_{x1}| + |R.c_{x2}| + |R.c_{x3}|} (4)

Y = \frac{|R.c_{x3}|}{|R.c_{x1}| + |R.c_{x2}| + |R.c_{x3}|} (5)

Here, R.c_{x1}, R.c_{x2}, R.c_{x3} represent the R.c of DTR to variables x_1, x_2, and x_3, respectively. We denoted the relative contributions of variables x_1, x_2, and x_3 to DTR as C, M, and Y, respectively. These contributions then served as the brightness values for the cyan, magenta, and yellow channels, respectively, culminating in the generation of a CMY image.

Random Forest regression analysis
We also utilized the Random Forest algorithm for regression analysis59. The input to the Random Forest was standardized in the same way as for the ridge regression analysis, and the data was also randomly partitioned into an 80% calibration dataset and a 20% validation dataset for the validation. During the modeling process at individual grid cells, the Random Forest algorithm leveraged the provided sequences of independent variables and dependent variables in the calibration dataset to train 100 decision trees. Each decision tree independently predicted DTR values based on the given predictor variables. Moreover, we employed out-of-bag (OOB) prediction error estimation, an inherent capability of the Random Forest algorithm. Additionally, we assessed the importance of predictors by enabling the OOB Predictor Importance feature, providing insights into the relative contributions of total cloud cover, aerosol optical depth, and soil moisture in predicting DTR. Similar to the ridge regression coefficients, these Predictor Importance features were subsequently utilized in generating a CMY image.

Partial correlation analysis and bootstrap method for assessing relationships
To determine the magnitude of the relationship between annual solar radiation and DTR, Tmax or Tmin at individual grid points, while controlling for the effect of soil moisture, we conducted a partial correlation analysis. Partial correlation analysis is a statistical technique used to assess the relationship between two variables while controlling for the influence of one or more additional variables. Here, both the variables involved in the correlation analysis and the variable being controlled for were specific to each individual grid cell. The significance of the partial correlations was evaluated at a threshold of p < 0.05.

The uncertainties of the partial correlations between global annual DTR and total cloud cover, soil moisture, or aerosol optical depth were assessed using the bootstrap method87. Specifically, we generated 1000 bootstrap samples through random samples with replacements from the original data to create samples of equal size to the original dataset. For each bootstrap sample, we computed the partial correlation coefficient between DTR and an environmental factor while adjusting for the potential confounding effect of other variables. Subsequently, we calculated the mean correlation coefficient across all bootstrap samples and determined the 95% confidence interval using the 2.5th to 97.5th percentile of the bootstrap distribution.

Data availability
All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. The source data underlying Figs. 1–4 are provided as Source Data files and have been deposited in the Figshare repository available at https://doi.org/10.6084/m9.figshare.24310699.v1. The Berkeley Earth Surface Temperature dataset is from https://berkeleyearth.org/data/. The Global Surface Summary of the Day dataset is from https://data.nodc.noaa.gov/cgi-bin/iso?id=gov.noaa.ncdc:C00516. The ERA5 temperature dataset is from https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-single-levels-monthly-means?tab=form. The MODIS cloud dataset is from https://ladsweb.modaps.eosdis.nasa.gov/. The CERES incident shortwave radiation dataset is from https://asdc.larc.nasa.gov/project/CERES/CERES_EBAF_Edition4.1. The MERRA-2 downwards shortwave radiation, aerosol optical depth and surface albedo dataset is from https://disc.gsfc.nasa.gov/datasets?project=MERRA-2. The GLEAM soil moisture data is from https://www.gleam.eu/#datasets.

Code availability
The code for the analysis and mapping can be obtained from the Figshare repository (https://doi.org/10.6084/m9.figshare.24310699.v1)88.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Review of “Reversed asymmetric warming of sub-diurnal temperature over land during recent decades\nby Zhong et al.\nThe paper documents a shift in the trends from the established diurnal asymmetric warming linked to\nearlier found processes that makes the climate change more efficient of increasing the nighftime\ntemperatures than during dayfime. The currently presented analysis shows a shift in that behavior, i.e.\ntoward processes that increase the maximum temperature more than earlier, thus a reversed\nasymmetric warming and even an overall increase in the diurnal temperature range (DTR). Such a result\nis significant enough to merit publicafion in a high-profile journal.\nThe authors aftribute the change in maximum temperature to a widespread decline in cloud cover,\ndecreased aerosol opfical depth in Western Europe and decreased soil moisture in southern South\nAmerica. Here the story gets a bit complicated as the observafions of these quanfifies are not as well\nestablished as the maximum and minimum temperatures, which are not as well established as the mean\ntemperature. The reliability of the datasets and other possible processes to aftribute the change to is not\ndiscussed in the manuscript, although some similar analyses of alternafive datasets are presented.\nEarlier published datasets have large white areas and other issues and crucial references discussing this\nare missing. I propose that the authors reduce the number of datasets used, discusses the uncertainfies\nin a robust way and sort out what reliably be concluded from the data. These aspects need to be\nimproved before considering publicafion.\nComments\nPage 1, Line 54: While I agree that the maximum and minimum temperatures provide useful\ninformafion, care is needed as it is much harder to produce a gridded data set of these fields. In a state-\nof-the-art dataset that I used five years ago, several land areas were not even reported due to these\nuncertainfies, such as most of S. America and Africa. What is the quality of these datasets now? How\ndoes the observafional coverage vary over the studied periods? Reading in AR6 there is a discussion on\nthe quality and reference to a study by Thorne et al. (2016a,b) that I do not see referenced here. It is also\nevident from Fig 1 that any regional conclusion are sensifive to the quality of the original data which\nvaries substanfially (would be good for the reader if the data source and period was clear in the figure\nand not only in the capfion).\nPage 2, Line 60: Here you list a number of factors, however, you mix physical parameters with manmade\nchanges which is confusing. I miss surface albedo for example, or do you consider that to be in the land\nuse change? The land use change covers more things though and I think that albedo might be befter\nobserved than soil moisture which I assume is supposed to cover one aspect of the drought you discuss\nbelow but changes in albedo are more widespread. It is also clear that the factors vary regionally and\nseasonally.\nPage 6, Line 138: The informafion on the analysis here is very limited and quite crypfical as the reader\ndoes not know which data that is used for the total cloud cover and what is a ridge regression and why\ndo you use that?\nPage 6, Line 140: Where is the total cloud cover from? It that the best variable to use?\nPage 7, Line 145: This has been seen in earlier studies, please refer to them.\nPage 7, Line 154: What type of cloud data goes into the CRU TS data, is it synop data? Can you be sure\nthat there is no fime dependent observafional bias in the datasets? The only one covering the earlier\nperiod is CRU TS. Is it independent or does it make use of the satellite data for the later period? Would\nbe good to have difference plots in Fig S3.\nPage 8, Figure 3. In some of your figures, North Africa is grey which I assume is of lack of data but it is\nnot commented and not clear which dataset that is lacking informafion and why.\nPage 8, Line 173: While I understand it is interesfing to try to see the role of different types of clouds, but\nusing cloud cover is a blunt method and using reanalysis, where the observafional material varies over\nfime, is tricky. Do we really trust the very strong trend in low-level clouds? Did you try to do the whole\nanalysis only in the reanalysis world? That would be interesfing to see if it shows the change in DTR\ntrends and the relafionship with clouds, soil moisture and aerosols.\nPage 9, Figure 4: What is the meaning of the dots?\nPage 10, Line 206: What does “controlling for PDSI or soil moisture” mean?\nPage 10, Line 208: We already know that surface solar radiafion influences the DTR!\nPage 10, Line 224: Do you mean that the surface of the air cools?\nPate 10, Line 227: The sentence starfing with “However”, states that there should be a limited role for\nthe drought and aerosols, that is not reflected in the abstract which claims that these two factors has a\nrole. The only conclusion is that the solar radiafion has changed due to changes in clouds and that\nresults in a changed trend for DTR.\nThere are a number of issues in the Methods secfion, please rewrite and make sure that everything is\nlogical and refer to the right dataset and explain all symbols and do not introduce symbols that are not\nused.\nI will not comment any details in the supplement at this stage.\nReferences\nThorne, P.W. et al., 2016a: Reassessing changes in diurnal temperature range: Intercomparison and\nevaluafion of exisfing global data set esfimates (2016b). Journal of Geophysical Research: Atmospheres,\n121(10), 5138–5158, doi:10.1002/2015jd024584.\nThorne, P.W. et al., 2016b: Reassessing changes in diurnal temperature range: A new data set and\ncharacterizafion of data biases (2016a). Journal of Geophysical Research: Atmospheres, 121(10), 5115–\n5137,", "reviewer_2": "The manuscript considers how the diurnal temperature range (DTR) has increased during the last 30\nyears, in contrast to the previous 30 years where the DTR decreased. The authors explain this reversal in\nthe DTR trend by comparing the contribufion due to aerosols, cloud cover and drought, and suggest that\nthe trend is primarily driven by changes in cloud cover and therefore incoming solar radiafion at the\nsurface, which increases Tmax more than Tmin.\nThe paper addresses an important trend and the results appear robust, even given various different\ndatasets. However, there are some aspects of the stafisfical methods that are not explained in sufficient\ndetail and I have two main concerns. Firstly, I expect more jusfificafion for the choice of regression and\nthe features used. Secondly, the ridge regression model used frequently throughout the paper is not\nvalidated and I find it difficult to trust the results without evidence that ridge regression is a suitable and\naccurate model. Assuming this to be true, the remainder of the stafisfical analysis on the regression\ncoefficients is valid and robust. The authors carefully address each potenfial driver of DTR changes and\nensure robustness by repeafing analyses on different datasets. The paper is clearly wriften and well\npresented. If the authors can show that the regression model is 1) jusfified by previous studies and 2)\naccurate in terms of predicfion, I expect this paper to be an interesfing contribufion to the community.\nMain Comments\nL137: To explore\n138 possible factors behind the changes in DTR in recent decades, we performed a ridge\n139 regression analysis49, reducing the impact of mulficollinearity on the esfimated\n140 relafionships during 1981–2020. …\nFurther explanafion of this method should be included here. Are the independent variables defined on\nevery single grid cell? Why does one expect mulfi-collinearity? The dependent variable is DTR, is this a\nsingle value such as the global mean DTR? Or is this the DTR at every grid point?\nL135-144\nBefore showing the values of the ridge regression coefficients, I think it is very important that the ridge\nregression model is validated. This can be done by leaving aside some samples from the training dataset\nand performing tesfing on these samples. Metrics such as R^2, RMSE, etc. should be used to ensure the\npredicfion from the ridge regression matches the test data. Metrics or plots highlighfing the accuracy of\nthis method should be included in the Supplementary. Then the full dataset can be used for the ridge\nregression in the remainder of the paper.\nLimitafions to the choice of model?\nAre there any possible limitafions to this study? I think there should be more jusfificafion for the model\nchoices made in the introducfion secfion, as well as comments on any limitafions towards the end of the\npaper.\nThe paper focuses on aerosols, clouds and drought, which I understand are likely to be the dominant the\ndominant drivers of temperature change, but I would like to see more jusfificafion. Have previous studies\nconfirmed that these are the main three drivers? I would also expect there to be addifional contributors\ntowards surface temperature changes that are not commented on. Can we be enfirely confident that\nthe regression does not ignore other important variables? What addifional variables could be explored or\nprobed in future studies?\nAlso, the ridge regression model used is a linear model. The authors do not comment on this\nassumpfion. Can this be jusfified based on previous studies?\nFig. 3.\nThis is a very interesfing composite plot that nicely summarizes the coefficients. I have one query which\nwould affect the output of the plot. Is the composite esfimated based on the relafive contribufion from\nall three coefficients? i.e.\nR = |coefficient_1|/(coefficient_1 + coefficient_2 + coefficient_3)\nG =| coefficient_2|/(coefficient_1 + coefficient_2 + coefficient_3)\nB = |coefficient_3|/(coefficient_1 + coefficient_2 + coefficient_3)\nIf so, this would mean that grid cells which have low coefficient values everywhere would sfill appear\nwith equal “amount” of color. I highlight this because I nofice that areas that have low magnitude values\nfor all three coefficients in Fig S2. e.g. Southern Western Africa, parts of Canada, Siberia and Brazil, are\nsfill “bright” in Fig 3. They appear a liftle noisy. This could be misleading. Could the authors clarify the\nmethod for construcfing the composite?\nI think the confusion may come from the menfion of “absolute” in the capfion which states “The colour\nof the composite was determined by the absolute value of the ridge 172 regression coefficients”, given\nthere are two possible meanings of absolute (“magnitude” or “absolute” vs. relafive). For clarity this\ncould be reworded, e.g. “The colour of the composite was determined by the relafive contribufion from\nthe magnitude of the ridge regression coefficients”.\nMethods\nIn some cases, the descripfion of the stafisfical methods are limited. I would expect some more detail on\nthe following:\n1. Ridge regression: how is the penalizafion term added (L2 norm)? How is the penalizafion term tuned?\nThe equafion would be beneficial to readers here. This would also help clarify what the independent and\ndependent variables are and how many there are (i.e. global mean or individual grid cells?)\n2. Parfial correlafion: this is not explained much for those who are not familiar with it and to clarify any\nchoices made. My main quesfions are: how do you “control” for PDSI? Do you take it to be the mean\nvalue everywhere?\nAlso references could be added to the methods for those who would like further details.\nMinor comments, grammar, typos\nL86: “The change in global average Tmax reached the warming rate in\nShould this be “the warming rate in Tmax” rather than “the change in global average Tmax” or are these\nequivalent?\nL87 In the last 30-year window (1991–2020), both two datasets\n88 show a stronger increase in the\n89 global average of Tmax than that of Tmin.\nThis is not clear to me for CRU TS in Fig. 1a, it appears they have a the same increase\nL133 – “All these evidences indicate” -> All this evidence indicates\nL173 – which dataset is used for the low, medium and high-level cloud cover?\nL276 – “resolufion. on a 0.25° grid” remove period.\nFig. S1 – are these global temperature trends?\nFig. S13, 14, etc. – “The black stripes mark ….” dots rather than stripes.", "reviewer_3": "The authors present analysis of the changes in the diurnal temperature range from two stafion-based\ngridded datasets of surface air temperature (BEST and CRU). They find that there is a clear change in the\nlate 1980s from a period where the DTR was decreasing due to the more rapid increase in diurnal\nminimum temperatures than maximum temperatures, to a period where either they are both increasing\nat the same rate or the diurnal maximum temperatures have been increasing more than diurnal\nminimum temperatures. Changes to the diurnal temperature range are extremely important to a range\nof impacts, and so this kind of analysis is very important for many groups. While the analysis seems\nreasonable and supports the conclusions I have a major concern about the inconsistency in the DTR\nchanges between the CRU and BEST datasets which would need to be addressed prior to publicafion.\nMajor concern:\nThere is a large difference in the results from CRU and BEST (Figure 1) so why not also apply the rest of\nthe analysis to BEST as well? I don’t agree with the reasoning to focus on CRU as the change is more\npronounced in BEST since the spafial pafterns and hence many of the conclusions that follow in the rest\nof the paper are going to depend on which dataset is used. This large systemafic uncertainty in changes\nin the DTR needs to be addressed to check for consistency in the findings prior to publicafion.\nShort comments:\nWhy are you showing the non-significant trends in Figure 2? Can’t be sure these are increases/decreases\nrather than noise. I think it’s more appropriate to put these together as non-significant trends to give a\nfuller picture of the surface area where we are seeing significant changes."}, "author_response": "Response to reviewer: 1\nThe authors attribute the change in maximum temperature to a widespread decline\nin cloud cover, decreased aerosol optical depth in Western Europe and decreased soil\nmoisture in southern South America. Here the story gets a bit complicated as the\nobservations of these quantities are not as well established as the maximum and\nminimum temperatures, which are not as well established as the mean temperature.\nThe reliability of the datasets and other possible processes to attribute the change to\nis not discussed in the manuscript, although some similar analyses of alternative\ndatasets are presented. Earlier published datasets have large white areas and other\nissues and crucial references discussing this are missing. I propose that the authors\nreduce the number of datasets used, discusses the uncertainties in a robust way and\nsort out what reliably be concluded from the data. These aspects need to be\nimproved before considering publication.\nRes: Thanks for your constructive suggestions. In the revised manuscript, we have\naddressed the above-mentioned points as follows:\n(1) We have provided a more comprehensive discussion regarding the uncertainties\nassociated with the data and the overall findings of the study.\n(2) We have reduced the number of datasets used to enhance clarity and reliability.\n(3) The issue of spatial gaps in some datasets is acknowledged and explained in the\nrevised manuscript.\nMoreover, we have made major revisions to the Methods section to provide more\nmethodological details and clarify several points raised by you and the other\nreviewers.\nMore details can be found in our point-by-point responses below.\nprovide useful information, care is needed as it is much harder to produce a gridded\ndata set of these fields. In a state-of-the-art dataset that I used five years ago, several\nland areas were not even reported due to these uncertainties, such as most of S.\nAmerica and Africa. What is the quality of these datasets now? How does the\nobservational coverage vary over the studied periods? Reading in AR6 there is a\ndiscussion on the quality and reference to a study by Thorne et al. (2016a,b) that I\ndo not see referenced here. It is also evident from Fig 1 that any regional conclusion\nare sensitive to the quality of the original data which varies substantially (would be\ngood for the reader if the data source and period was clear in the figure and not only\nin the caption).\nRes: You raised an important point regarding the difficulties in producing gridded\ndatasets for maximum and minimum temperatures. We agree that there are\nuncertainties and challenges in obtaining complete coverage, especially in regions like\nSouth America and Africa where observational coverage remains limited. In our study,\nwe used state-of-the-art datasets available at the time, but it is still worth noting that\nthe quality and coverage of these datasets vary over the studied periods.\nIn the revised manuscript, we have discussed the limited observational coverage in\nSouth America and Africa, specifically referencing the study by Thorne et al. (2016a,b)\nfor further context. Please see lines 271-286.\nMoreover, we have conducted additional analyses, which revealed that even after\nexcluding regions with higher uncertainties (most of South America and Africa), the\naverage warming rate of maximum temperatures over the past three decades remains\nhigher than that of minimum temperatures. Please see Fig. R1.\nin Tmax (red) and Tmin (blue). The trends were calculated using a 30-year moving window over\n1961–2020. The region highlighted in red in the insets delineate the study area here.\nFurthermore, we have revised Fig. 1 following your suggestion to clearly indicate the\ndata source and period in the figure itself, please see Fig. 1 in the revised manuscript.\nPage 2, Line 60: Here you list a number of factors, however, you mix physical\nparameters with manmade changes which is confusing. I miss surface albedo for\nexample, or do you consider that to be in the land use change? The land use\nchange covers more things though and I think that albedo might be better\nobserved than soil moisture which I assume is supposed to cover one aspect of the\ndrought you discuss below but changes in albedo are more widespread. It is also clear that\nthe factors vary regionally and seasonally.\nRes: Changes in the physical parameters that were included can be manmade and\nnatural. This study does not go into the causes of these changes. Regarding your\nquestion about the influence of surface albedo and soil moisture on DTR, we\nacknowledge that the original analysis overlooked the potential influence of surface\nalbedo changes on DTR. Therefore, we have conducted additional partial correlation\nanalyses in the revised manuscript, as shown in Fig. R2 below. Fig. R2a shows the\nspatial distribution of partial correlation coefficients between soil moisture and DTR,\n3\ncontrolling for the effects of cloud cover, aerosols, and surface albedo. Further, Fig.\nR2b shows the partial correlation coefficients between surface albedo and DTR,\ncontrolling for the effects of cloud cover, aerosols, and soil moisture. Overall, we\nfound a significant negative correlation between soil moisture and DTR on a global\nscale, while the correlation between surface albedo and DTR was relatively weak and\nvaried in different regions. A possible reason for this relationship between surface\nalbedo and DTR is that a decrease in surface albedo can lead to increased absorption\nof solar radiation and sensible heat flux at the surface, resulting in higher maximum\ntemperatures and DTR. Meanwhile, a decrease in surface albedo is often associated\nwith an increase in vegetation greenness. The cooling effect of vegetation\nevapotranspiration can reduce maximum temperatures and DTR, resulting in the\nobserved regional variations in the correlation between surface albedo and DTR. We\nhave included a detailed discussion on this topic in the main text, as mentioned in\nlines 259-265 and 301-308.\nFig. R2. Spatial distribution of partial correlation coefficient between yearly diurnal temperature\nrange (DTR) and soil moisture (a) or surface albedo (b) during 1981-2020. The soil moisture data\nwas obtained from the Global Land Evaporation Amsterdam Model dataset and the surface albedo\ndata was obtained from MERRA-2. The black dots mark the areas where correlations are\nsignificant at the p<0.05 level.\ncryptical as the reader does not know which data that is used for the total cloud\ncover and what is a ridge regression and why do you use that?\nRes: Sorry that the source of information for clouds were not provided. The total\ncloud cover data used in our study is derived from the ERA5 dataset, which is now\nclarified in lines 152-157 in the revised manuscript.\nAs for the use of ridge regression, it is a statistical method employed to reduce the\nimpact of multi-collinearity among independent variables that influence the variation\nin DTR. Given the presence of multi-collinearity among the explanatory variables\nwe considered, we deemed it necessary to use ridge regression to mitigate the impact\nof this issue. We have incorporated an explanation of ridge regression in the revised\nmanuscript, see lines 140-161 and the expanded “Ridge regression” subsection in\nMethods.\nPage 6, Line 140: Where is the total cloud cover from? It that the best variable to\nuse?\nRes: The total cloud cover data is sourced from the ERA5 dataset. Among the\ncommonly used reanalysis datasets for cloud cover, only two datasets cover the\nperiod from 1961 to 2020: ERA5 and CRU TS. We chose the ERA5 dataset because\nthe generation process of the cloud cover dataset in CRU TS involves the use of\nDTR. Specifically, DTR station anomalies were used to estimate ‘synthetic’ CLD\nstation anomalies since 2003, by a linear transformation with a scaling factor and\nmean offset calculated from CRU TS gridded CLD and DTR values for each latitude\nband (Harris et al. 2014). Therefore, to ensure the integrity of the main analysis on the\nrelationship between DTR and cloud cover, we did not use the CRU TS dataset.\nConsequently, in the main text of the manuscript, we rely on the cloud cover data\n5\nfrom the ERA5 dataset for our primary analysis, which is now mentioned in lines\n152-154.\nIn addition to total cloud cover, the diurnal temperature range is also directly related\nto surface incoming radiation, as we will discuss in later sections. However, since the\nmagnitude of surface incoming radiation is predominantly determined by total cloud\ncover, we have decided not to include both variables as independent variables in this\nanalysis.\nReference cited:\nHarris, I., Jones, P.D., Osborn, T.J. and Lister, D.H. (2014), Updated high-resolution grids of\nmonthly climatic observations – the CRU TS3.10 Dataset. Int. J. Climatol., 34: 623-642.\nhttps://doi.org/10.1002/joc.3711\nPage 7, Line 145: This has been seen in earlier studies, please refer to them.\nRes: Point taken. We have cited relevant references at the specified location in the\nrevised manuscript.\nPage 7, Line 154: What type of cloud data goes into the CRU TS data, is it synop\ndata? Can you be sure that there is no time dependent observational bias in the\ndatasets? The only one covering the earlier period is CRU TS. Is it independent or\ndoes it make use of the satellite data for the later period? Would be good to have\ndifference plots in Fig S3.\nRes: The cloud cover data in CRU TS is a fused product derived from daily sunshine\nduration and diurnal temperature range observations, using empirical formulas.\nHowever, we cannot guarantee the absence of time-dependent observational bias\nduring the data collection process. Therefore, we avoided using the cloud cover data\nfrom the CRU TS dataset.\nprimarily relies on the ERA5 reanalysis dataset, which covers the period from 1961\nto 2020. The ERA5 dataset assimilates multiple satellite data sources and\nobservational records. Therefore, we chose to analyze the ERA5 dataset due to its\ncomprehensive coverage and assimilation of various data sources.\nPage 8, Figure 3. In some of your figures, North Africa is grey which I assume is\nof lack of data but it is not commented and not clear which dataset that is lacking\ninformation and why.\nRes: Sorry for not making this clear. The gray color in the North Africa region is due\nto the absence of PDSI data specifically in the Sahara region. In the revised\nmanuscript, we have included an explicit explanation in the figure caption to clarify\nthe reason for the missing information.\nPage 8, Line 173: While I understand it is interesting to try to see the role of\ndifferent types of clouds, but using cloud cover is a blunt method and using\nreanalysis, where the observational material varies over time, is tricky. Do we\nreally trust the very strong trend in low-level clouds? Did you try to do the whole\nanalysis only in the reanalysis world? That would be interesting to see if it shows\nthe change in DTR trends and the relationship with clouds, soil moisture and\naerosols.\nRes: As you rightly pointed out, the trends in high, medium, and low-level cloud\ncover have a considerable amount of uncertainty, and we cannot fully trust the strong\ndownward trend in low-level clouds in the reanalysis dataset. Therefore, our revised\nanalysis is focused on total cloud cover. The issue of inhomogeneity (and thus\n7\nquality) of reanalysis is relevant. But on a positive note, the latest IPCC report (Chen\net al., 2021) has assessed this issue and the improvement of reanalysis over the last\nyears and its usefulness for climate change studies have been highlighted. Thus, we\nhave some confidence in the ERA5 data. Based on your advice, we have revised this\nsection to analyze the global long-term trends of DTR, cloud cover, soil moisture,\nand aerosols, as well as their relationships using reanalysis data only. Please see Fig.\nS6 and lines 198-212 in the revised manuscript.\nReference cited:\nChen, D., M. Rojas, B.H. Samset, K. Cobb, A. Diongue Niang, P. Edwards, S. Emori, S.H. Faria,\nE. Hawkins, P. Hope, P. Huybrechts, M. Meinshausen, S.K. Mustafa, G.-K. Plattner, and A.-M.\nTréguier, 2021: Framing, Context, and Methods. In Climate Change 2021: The Physical Science\nBasis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental\nPanel on Climate Change\" [Masson-Delmotte, V., P. Zhai, A. Pirani, S.L. Connors, C. Péan, S.\nBerger, N. Caud, Y. Chen, L. Goldfarb, M.I. Gomis, M. Huang, K. Leitzell, E. Lonnoy, J.B.R.\nMatthews, T.K. Maycock, T. Waterfield, O. Yelekçi, R. Yu, and B. Zhou (eds.)]. Cambridge\nUniversity Press, Cambridge, United Kingdom and New York, NY, USA, pp. 147–286,\ndoi:10.1017/9781009157896.003.\nPage 9, Figure 4: What is the meaning of the dots?\nRes: Sorry for having missed the explanation. The black dots represent the regions\nwhere the trends reached a level of significance of p<0.05. We have clarified this in\nthe caption of Fig. 4 in the revised manuscript.\nPage 10, Line 206: What does “controlling for PDSI or soil moisture” mean?\nRes: It means that in the partial correlation analysis, we considered and adjusted for\nthe influence of either PDSI or soil moisture on the relationship between DTR and\nsolar radiation. We have revised the imprecise wording in the revised manuscript.\nDTR!\nRes: Point taken. We have modified the inaccurate wording and provided a more\nprecise statement regarding the influence of surface solar radiation on DTR.\nAdditionally, we have included relevant citations in the revised manuscript.\nPage 10, Line 224: Do you mean that the surface of the air cools?\nRes: Yes, we intended to convey the message that there is a cooling effect of soil\nmoisture evaporation on air temperature. We have provided further clarification on\nthis point in the revised manuscript.\nPate 10, Line 227: The sentence starting with “However”, states that there should\nbe a limited role for the drought and aerosols, that is not reflected in the abstract\nwhich claims that these two factors has a role. The only conclusion is that the\nsolar radiation has changed due to changes in clouds and that results in a changed\ntrend for DTR.\nRes: In the revised manuscript, we have modified the expression in the abstract to\nensure consistency with the statement in this section. Please refer to the revised\nabstract for the updated wording.\nThere are a number of issues in the Methods section, please rewrite and make sure\nthat everything is logical and refer to the right dataset and explain all symbols and\ndo not introduce symbols that are not used..\nMethods section to ensure logical flow, accurate referencing of datasets, and clear\nexplanations of symbols. We have also removed any symbols that are not used in the\nstudy.\nThe paper addresses an important trend and the results appear robust, even given\nvarious different datasets. However, there are some aspects of the statistical methods\nthat are not explained in sufficient detail and I have two main concerns. Firstly, I\nexpect more justification for the choice of regression and the features used. Secondly,\nthe ridge regression model used frequently throughout the paper is not validated\nand I find it difficult to trust the results without evidence that ridge regression is a\nsuitable and accurate model. Assuming this to be true, the remainder of the\nstatistical analysis on the regression coefficients is valid and robust. The authors\ncarefully address each potential driver of DTR changes and ensure robustness by\nrepeating analyses on different datasets. The paper is clearly written and well\npresented. If the authors can show that the regression model is 1) justified by\nprevious studies and 2) accurate in terms of prediction, I expect this paper to be an\ninteresting contribution to the community.\nRes: Thank you for your constructive and positive feedback. Following your\nsuggestions, we have significantly expanded the justification and explanation for\nselecting ridge regression as our main analysis method. We have also conducted\nadditional validation to demonstrate the accuracy of the ridge regression model in\nterms of prediction. Finally, we have conducted additional analyses using a Random\nForest regression model to confirm that the linear assumption in the ridge regression\nis reasonable. These improvements have been incorporated into the revised manuscript,\nproviding a more comprehensive understanding of the statistical methods used.\nIn addition to these new analyses, the revised manuscript also includes additional\ncitations of previous studies to better support our choices of methods and findings.\n138 possible factors behind the changes in DTR in recent decades, we performed a\nridge\n139 regression analysis49, reducing the impact of multicollinearity on the estimated\n140 relationships during 1981–2020. …\nFurther explanation of this method should be included here. Are the independent\nvariables defined on every single grid cell? Why does one expect multi-collinearity?\nThe dependent variable is DTR, is this a single value such as the global mean DTR?\nOr is this the DTR at every grid point?.\nRes: The ridge regression analysis was indeed performed on a grid-based level, where\nboth the independent variables and the dependent variable (DTR) were considered at\nindividual grid points. This has been clarified in the revised manuscript.\nAs for the issue of multi-collinearity, we have discussed it in the specified paragraph\nof the revised manuscript at lines 142-150. Specifically, there is a strong correlation\namong the potential drivers of DTR, which include total cloud cover, aerosols, and\nsoil moisture. For example, the phenomenon of aerosol-cloud interactions has received\nsignificant support from numerous published reports, with aerosols influencing cloud\nalbedo and cloud lifetime. Additionally, clear associations exist between cloud cover\nand soil moisture. Increased cloud cover has the potential to enhance precipitation,\nleading to soil wetting. The complex interactions among these factors indicate a high\ndegree of multicollinearity when conducting regression analysis, making it\nchallenging to identify the dominant drivers of DTR changes.\nL135-144\nimportant that the ridge regression model is validated. This can be done by leaving\naside some samples from the training dataset and performing testing on these\nsamples. Metrics such as R^2, RMSE, etc. should be used to ensure the prediction\nfrom the ridge regression matches the test data. Metrics or plots highlighting the\naccuracy of this method should be included in the Supplementary. Then the full\ndataset can be used for the ridge regression in the remainder of the paper.\nRes: We agree with you. Following your advice, we have partitioned all the data into\nan 80% training set and a 20% validation set to validate the ridge regression model.\nMeanwhile, the full dataset was employed to identify the ridge regression coefficients.\nRegions with dense distributions of temperature measurement stations, including\nNorth America, Europe, East Asia, and Australia, exhibit high R-squared values in\nthe validation datasets. Please see Fig. S3a in the revised Supplementary Information.\nLimitations to the choice of model?\nAre there any possible limitations to this study? I think there should be more\njustification for the model choices made in the introduction section, as well as\ncomments on any limitations towards the end of the paper?\nRes: In the revised manuscript, we have included more justification for the choice of\nridge regression analysis in the introduction section of the manuscript, specifically in\nthe context of \"Potential mechanisms behind reversing asymmetric warming\". We\nexplain why this modeling approach was employed and its relevance to our research\nobjectives. Please see lines140-161.\nthe end of the paper, discussing the uncertainties associated with its application. Please\nsee lines 287-297.\nThe paper focuses on aerosols, clouds and drought, which I understand are likely\nto be the dominant the dominant drivers of temperature change, but I would like to\nsee more justification. Have previous studies confirmed that these are the main three\ndrivers? I would also expect there to be additional contributors towards surface\ntemperature changes that are not commented on. Can we be entirely confident that\nthe regression does not ignore other important variables? What additional variables\ncould be explored or probed in future studies?\nRes: While there is extensive literature suggesting that cloud cover, aerosol optical\ndepth, and soil moisture can significantly influence DTR, they are not the exclusive\ndrivers. In our revised manuscript, we have provided additional discussion on the\ncomplex influence of precipitation, atmospheric water vapor, albedo and land use/land\ncover changes (LULCC) on DTR, recognizing that it may be other important\ncontributors to surface temperature changes that were not explicitly mentioned. We\nalso suggested that the effect of LULCC on DTR could further be explored or probed\nto gain a more comprehensive understanding of DTR variability. Please see lines 294-\n308.\nAlso, the ridge regression model used is a linear model. The authors do not comment\non this assumption. Can this be justified based on previous studies?\nRes: As far as we know, there is no previous study that explicitly justified the linear\nassumption for this application. Thus, we have added an analysis using a Random\nregression model is capable of capturing complex, non-linear relationships. Our\nanalysis using the Random Forest model yields results consistent with the ridge\nregression findings, which indicates that neglecting nonlinear relationships does not\nhave a substantial impact on the research findings. This is discussed in the revised\nmanuscript in lines 181-188.\nFig. 3.\nThis is a very interesting composite plot that nicely summarizes the coefficients. I\nhave one query which would affect the output of the plot. Is the composite estimated\nbased on the relative contribution from all three coefficients? i.e.\nR = |coefficient_1|/(coefficient_1 + coefficient_2 + coefficient_3)\nG =| coefficient_2|/(coefficient_1 + coefficient_2 + coefficient_3)\nB = |coefficient_3|/(coefficient_1 + coefficient_2 + coefficient_3)\nIf so, this would mean that grid cells which have low coefficient values everywhere\nwould still appear with equal “amount” of color. I highlight this because I notice\nthat areas that have low magnitude values for all three coefficients in Fig S2. e.g.\nSouthern Western Africa, parts of Canada, Siberia and Brazil, are still “bright” in\nFig 3. They appear a little noisy. This could be misleading. Could the authors clarify\nthe method for constructing the composite?\nI think the confusion may come from the mention of “absolute” in the caption\nwhich states “The colour of the composite was determined by the absolute value of the\nridge 172 regression coefficients”, given there are two possible meanings of absolute\n(“magnitude” or “absolute” vs. relative). For clarity this could be reworded, e.g.\n15\n“The colour of the composite was determined by the relative contribution from the\nmagnitude of the ridge regression coefficients”?\nRes: You have correctly understood our analysis process. In the revised manuscript,\nwe have revised the corresponding description to avoid any potential\nmisunderstanding. Please see lines 192-197. Additionally, we have provided further\nexplanations in the Methods section to clarify the calculation process, please see lines\n426-433. These revisions aim to improve the clarity and understanding of our\nmethodology.\nMethods\nIn some cases, the description of the statistical methods are limited. I would expect\nsome more detail on the following:\n1. Ridge regression: how is the penalization term added (L2 norm)? How is the\npenalization term tuned? The equation would be beneficial to readers here. This\nwould also help clarify what the independent and dependent variables are and how\nmany there are (i.e. global mean or individual grid cells?)\nRes: Thank you for your valuable feedback on the statistical methods described in our\nmanuscript. For each grid point in the regression analysis, the penalization term λ was\ninitially set to 0 and incremented by a small step size of 0.01. As λ increased, multi-\ncollinearity decreased, leading to a corresponding decline in the Variance Inflation\nFactor (VIF) value. The incrementation of λ ceased when the VIF value dropped below\n3, and this specific λ value was identified as the penalization term for the grid point.\nIn the revised manuscript, we have provided a more comprehensive explanation of the\npenalization term, including the equation itself, to help readers better understand the\nmethodology.\nindividual grid cell, with cloud cover, aerosol optical depth, and soil moisture as the\nindependent variables, and DTR as the dependent variable. These details are now\nelaborated upon in the revised Methods section. Please see lines 386-425.\n2. Partial correlation: this is not explained much for those who are not familiar with\nit and to clarify any choices made. My main questions are: how do you “control”\nfor PDSI? Do you take it to be the mean value everywhere?\nRes: In our revised manuscript, we have provided a more comprehensive explanation\nof the partial correlation analysis, clarifying that it aims to assess the relationship\nbetween two variables while accounting for the influence of other variables. This\ntechnique helps isolate the associations between the variables of interest.\nRegarding the control variable PDSI, we would like to clarify that we did not consider\nit as the mean value everywhere. Instead, we collected the specific values of solar\nradiation, DTR, and PDSI for each grid cell within our study area. These values were\nthen used in the partial correlation analysis to evaluate the direct relationship between\nsolar radiation and DTR while controlling for the potential confounding effects of\nPDSI.\nAlso references could be added to the methods for those who would like further\ndetails.\nRes: We have followed your suggestion by adding appropriate references to the\nMethods section.\nTS”\nShould this be “the warming rate in Tmax” rather than “the change in global\naverage Tmax” or are these equivalent?.\nRes: Yes, it should indeed be“the warming rate in Tmax”. We have corrected it in the\nrevised manuscript.\nL87 In the last 30-year window (1991–2020), both two datasets\n88 show a stronger increase in the\n89 global average of Tmax than that of Tmin.\nThis is not clear to me for CRU TS in Fig. 1a, it appears they have a the same\nincrease.\nRes: We have revised the corresponding statement to make it clearer. Please see lines\n88-90.\nL133 – “All these evidences indicate” -> All this evidence indicates.\nRes: Revised as suggested.\nL173 – which dataset is used for the low, medium and high-level cloud cover?\nanalysis. However, based on the suggestion from Reviewer 1, in the revised version,\nwe have replaced this section with an analysis of the annual variations of DTR and\nother environmental variables due to the large uncertainty in low, medium and high-\nlevel cloud cover in the reanalysis, especially in the earlier period.\nL276 – “resolution. on a 0.25° grid” remove period.\nRes: Fixed.\nFig. S1 – are these global temperature trends?\nRes: Yes. We now mention this in Fig. S2 in the revised Supplementary Information.\nFig. S13, 14, etc. – “The black stripes mark ….” dots rather than stripes.\nRes: Done.\nMajor concern:\nThere is a large difference in the results from CRU and BEST (Figure 1) so why not\nalso apply the rest of the analysis to BEST as well? I don’t agree with the reasoning\nto focus on CRU as the change is more pronounced in BEST since the spatial\npatterns and hence many of the conclusions that follow in the rest of the paper are\ngoing to depend on which dataset is used. This large systematic uncertainty in\nchanges in the DTR needs to be addressed to check for consistency in the findings\nprior to publication.\nRes: We agree with you. The selection of DTR data needs to be handled with caution\nas it can significantly impact the subsequent analysis outcomes. In the revised\nmanuscript, we have taken measures to address this concern. We validated both sets\nof DTR data against DTR observations and found that BEST data exhibited higher\noverall accuracy with the DTR observations. Consequently, we chose to use BEST\ndata for the analysis of factors influencing DTR changes. Please see Fig. S1 and lines\n121-127 in the revised manuscript. The conclusions drawn from this analysis align\nclosely with those derived from the previous analysis using CRU data. Additionally,\nwe have provided a comprehensive discussion on the uncertainties associated with\nboth datasets. Please sees line 267-286.\nShort comments:\nWhy are you showing the non-significant trends in Figure 2? Can’t be sure these\nare increases/decreases rather than noise. I think it’s more appropriate to put these\ntogether as non-significant trends to give a fuller picture of the surface area where\nwe are seeing significant changes?\nrevised manuscript.\nReviewer #1 (Remarks to the Author):\nThe authors have done a good job on revising the manuscript based on my and the other reviewers’\ncomments. I am mostly happy with the current version although I have some smaller suggesfions for\nclarity.\nLine 94: Would be interesfing to know the fracfion of the area that overlap between CRU TS and BEST. Is\nthere a resolufion or data treatment difference that makes out part of this difference? Seems like the\nCRU is much smoother field.\nLines 97 and 98: There is not enough confidence in the data to give the percentages with one decimal,\nplease write 63 and 35 % instead.\nLine 210: What about the influence of the clouds during nighftime? Less clouds during the night leads to\nmore outgoing radiafion and possibly lower minimum temperatures.\nLine 240: It is kind of obvious that you found a stronger correlafion between solar radiafion during\ndayfime, during nighftime there is no sunlight! However, there could be an influence on the nighftime\ntemperatures that remain from the day before, but that is much harder to show and does not make\nmuch sense to go into that discussion here.\nLine 287: Even though you have menfioned regression analysis I think you need to write regression\nmodels as it is easily misunderstood to be other types of models (e.g. climate models) is you only write\nmodels.\nLines 349 and 353: Why do you use two different resolufions of the ERA5? Or is it just two ways to write\nthe same resolufion? Please be consistent and remember that there is no exact translafion of a lat-lon\ngrid to km, it is not equal in north-south nor constant with lafitude.\nFigure 2: It is very difficult to disfinguish between the small points with no trend and the blue dots.\nPlease use a different color scheme. What does the * indicate, please provide informafion in the capfion.\nFigure 3: Please change the color scale to something which is possible for color blind to disfinguish\nbetween, green and red is the worst combinafion. There are a number of websites that provide color\npallets with enough contrast for everyone.\nFigure 4: I think it is good if you menfion the data source for each figure. At least in this one I am missing\nthat informafion.\nReviewer #2 (Remarks to the Author):\nThe authors have done a nice job addressing my previous comments. I greatly appreciate the addifional\ndetail on ridge regression in both the main text and the methods. The newly added random forest\nanalysis also strengthens the paper. The authors have also carefully outlined the limitafions of the study\nin the new paragraph starfing L289, which is addresses many of my comments regarding limitafions and\nchoice of features. Aside from the following two very minor comments, I believe the manuscript is\npublicafion ready.\nL406: I am not familiar with the meaning of “z-score anomalies”\nL158 “Ridge regression is a linear regularizafion method …”\nTo be more precise, I would refer to ridge regression as a regression method that uses regularizafion,\nrather than a regularizafion method in itself. I suggest something along the lines of “Ridge regression is a\nlinear regression method with regularizafion …”\nReviewer #3 (Remarks to the Author):\nI am safisfied that the authors have addressed the concerns I raised.\nResponse to reviewer: 1\nThe authors have done a good job on revising the manuscript based on my and the\nother reviewers’ comments. I am mostly happy with the current version although I\nhave some smaller suggestions for clarity.\nRes: Thank you, we very much appreciate your comments on how to further improve\nthe manuscript and have revised the text accordingly, as detailed below.\nLine 94: Would be interesting to know the fraction of the area that overlap between\nCRU TS and BEST. Is there a resolution or data treatment difference that makes\nout part of this difference? Seems like the CRU is much smoother field.\nRes: Good suggestion! We have included additional information regarding the\noverlap of areas for both datasets. In response to your concern about the potential\nimpact of resolution, we have conducted additional analyses. We resampled the CRU\nTS temperature data to a 1-degree latitude-longitude grid (matching the resolution of\nthe BEST data) and recomputed the DTR, as shown in Fig. R1 below. The results are\nconsistent with those obtained at the original resolution (0.5-degree latitude-longitude\ngrid), indicating that resolution differences may not be a significant factor\ncontributing to the variations between the two datasets. We believe that the reason\nwhy CRU is smoother is how the interpolation/extrapolation was done (Xu et al. 2018).\nXu, W., Li, Q., Jones, P. et al. A new integrated and homogenized global monthly land surface air\ntemperature dataset for the period since 1900. Clim Dyn 50, 2513–2536 (2018).\nhttps://doi.org/10.1007/s00382-017-3755-1\n(b). Here, the CRU TS temperature data was resampled to a 1-degree latitude-longitude grid.\nLines 97 and 98: There is not enough confidence in the data to give the percentages\nwith one decimal, please write 63 and 35 % instead.\nRes: We have made the suggested change accordingly.\nLine 210: What about the influence of the clouds during nighttime? Less clouds\nduring the night leads to more outgoing radiation and possibly lower minimum\ntemperatures.\nRes: Good question. We have included an analysis of the influence of cloud cover on\nboth maximum and minimum temperatures, as shown in the new Fig. S7, which is also\nincluded below. Our findings indicate that, in some regions, there is indeed a positive\ncorrelation between cloud cover and minimum temperatures, likely due to the\ngreenhouse effect of clouds that you mentioned. However, overall, a negative\ncorrelation between cloud cover and maximum temperatures appears to be more\nprevalent. This suggests that the effect of cloud cover on daytime solar incoming\nradiation generally plays a more important role in controlling the DTR. This is now\nmentioned on Lines 209–219.\ncover and maximum temperature (a) or minimum temperature (b) during 1981-2020, while\ncontrolling for the influence of soil moisture. The black dots mark the areas where correlations are\nsignificant at the p < 0.05 level. The total cloud cover was obtained from ERA5 dataset. (c) to (d),\nsame as (a) to (b), but the total cloud cover was obtained from MODIS (2003-2020) dataset.\nLine 240: It is kind of obvious that you found a stronger correlation between solar\nradiation during daytime, during nighttime there is no sunlight! However, there\ncould be an influence on the nighttime temperatures that remain from the day\nbefore, but that is much harder to show and does not make much sense to go into\nthat discussion here.\nRes: Another effect is that the minimum temperature often occurs shortly after sunrise\nand thus could be affected by the amount of solar radiation around sunrise, but we\nagree that these effects are harder to show using our data sources and have therefore\nremoved this part from the manuscript.\nLine 287: Even though you have mentioned regression analysis I think you need\nto write regression models as it is easily misunderstood to be other types of models\n(e.g. climate models) is you only write models.\nLines 349 and 353: Why do you use two different resolutions of the ERA5? Or is it\njust two ways to write the same resolution? Please be consistent and remember that\nthere is no exact translation of a lat-lon grid to km, it is not equal in north-south\nnor constant with latitude.\nRes: The ERA5 data is produced and archive on a reduced Gaussian grid with a\nresolution of N320 and decreasing number of longitude points toward the poles.\n(Some upper air data is archived as spectral coefficients, but we will not go into\nthese details here.) The native resolution of ERA5 is thus around 31 km. However,\nERA5 data is commonly provided on a 0.25° × 0.25° regular lat/lon grid (through\nbilinear interpolation for continuous data), which is also what we used here.\nNevertheless, we think it is important to also point out the native resolution of the\ndata. We agree that the previous text could cause confusion and have therefore\nrevised it as follows:\n“Monthly data on total cloud cover were obtained from the fifth-generation ECMWF\nreanalysis (ERA5) dataset56 on a 0.25°×0.25° regular latitude–longitude grid (the\nnative resolution of ERA5 is about 31 km).\n…\nThe monthly incident shortwave radiation data in all-sky conditions were obtained\nfrom the ERA5 dataset on a 0.25° grid …”\nFigure 2: It is very difficult to distinguish between the small points with no trend\nand the blue dots. Please use a different color scheme. What does the * indicate,\nplease provide information in the caption.\neasier to distinguish between the small points with no trend and the blue dots. We\nhave also updated the caption for improved clarity.\nFigure 3: Please change the color scale to something which is possible for color\nblind to distinguish between, green and red is the worst combination. There are a\nnumber of websites that provide color pallets with enough contrast for everyone.\nRes: Good suggestion! In the revised manuscript, we have replaced the red-green-\nblue (RGB) color scale with a cyan-magenta-yellow (CMY) composite plot, which\nenhances the figure's readability and ensures it is distinguishable for individuals with\nthe most common color vision impairments. Please see the updated plots in Figure 3,\nSupplementary Figure 4, and Supplementary Figure 5.\nFigure 4: I think it is good if you mention the data source for each figure. At least\nin this one I am missing that information.\nRes: In the revised manuscript, we have included the data source in the caption for\nall figures.\nL406: I am not familiar with the meaning of “z-score anomalies”.\nRes: Z-score are calculated by taking the difference between a data point and the mean\nof the dataset and then dividing it by the standard deviation. This standardization helps\nin comparing and interpreting data in a consistent manner, as it provides a measure of\nhow unusual or extreme a data point is within the context of the dataset. We have added\nan explanation about z-score anomalies in the revised manuscript.\nL158 “Ridge regression is a linear regularization method …”\nTo be more precise, I would refer to ridge regression as a regression method that\nuses regularization, rather than a regularization method in itself. I suggest\nsomething along the lines of “Ridge regression is a linear regression method with\nregularization …”\nRes: We agree with you and have made the suggested change accordingly. Thank you\nfor your suggestions!\n6"}]

Decision letter