It is unequivocal that human influence has warmed the atmosphere, leading to more frequent, more intense, longer-lasting, and more widespread hot temperature extremes worldwide (IPCC, 2021). These hot extremes are a major cause of suffering and death as punctuated by some of the deadly cases such as the 2003 European event (∼70,000 mortality) and the 2010 Russian event (∼50,000 mortality, Hoag, 2014). The combination of extremely high temperatures with high humidity, i.e., the so-called humid-heat stress, represents a greater threat to human health, as it lowers the cooling efficiency of sweat, therefore making it difficult or even impossible for the body to prevent overheating (Mora et al., 2017; Buzan and Huber, 2020; IPCC, 2021; Vecellio et al., 2022). Hence, considerable efforts have been devoted to defining, monitoring and understanding extreme humid-heat events, based on diverse metrics configuring temperature and humidity as well as other relevant variables (e.g., winds and radiation) in different ways (Delworth et al., 1999; Diffenbaugh et al., 2007; Stull, 2011; Willett and Sherwood, 2012; Fischer and Knutti, 2013).

Regardless of metrics used, humid heat extremes have exhibited significantly increasing trends for frequency and intensity in most continents (Rogers et al., 2021; Tuholske et al., 2021). When incorporating humidity into the indices in a highly non-linear manner, these increases are markedly larger than their temperature-only counterparts (Delworth et al., 1999; Wang and Zhu, 2020; Rogers et al., 2021). Even at the current level of global warming, several populous regions, including North and South China, Eastern India, and the Middle East (Freychet et al., 2020; Saeed et al., 2021; Raymond et al., 2020; Mora et al., 2017), have experienced humid-heat extremes very close to or even above upper physiological limit. Though the survivability limits are still exceeded only on the rarest of occasions for now, the likelihood of exceedance is expected to grow rapidly in future warmer climates in these hotspots (Dunne et al., 2013; Mora et al., 2017; Kang and Eltahir, 2018; Saeed et al., 2021).

For the bivariate extremes, data quality of both air temperature and relative humidity (RH) matters to the accuracy of estimate for their long-term changes, and hence to assessments of related human health impacts and risks as well as adaptation planning (Sherwood, 2018; Brouillet and Joussaume, 2019). In contrast to high-quality air temperature datasets provided by multiple agencies, observational records of relative humidity are subject to large uncertainties resulting primarily from data inhomogeneity, and thus should be used with caution (IPCC, 2021; Schröder et al., 2019). It has been observed that global mean and boreal midlatitude summer-mean surface RH slightly increased during 1973–2000, followed by a steep decline based on the HadISDH (Willett et al., 2014; Douville and Plazzotta, 2017; Dunn et al., 2017); however, none of the CMIP5 models could capture such evolution. It is not yet clear if this discrepancy is related to the misrepresentation of internal variability in simulations or observational data inhomogeneity. Byrne and O’Gorman (2018) provided a candidate theory suggesting that the observed negative trends for relative humidity over midlatitude continents were linked to warming over the neighbouring oceans. The observational and theoretical uncertainties along with the observation-simulation discrepancy, preclude any robust detection and attribution conclusion with respect to relative humidity changes (IPCC, 2021).

The inhomogeneity issue in relative humidity records stands out starkly over China, the southeastern parts in particular, which may be related to the technical defects of humidity sensors under high temperature and humidity environments there (Yu and Mu, 2008; Zhu et al., 2015). The dielectric parameter of the humicap varies volatilely at high humidity and even more drastically combined with high temperatures, making observed humidity deviating far from real values. Using quality-controlled yet non-homogenized observations, the annual China Climate Bulletin released before 2014 used to report a significant decreasing trend in nationwide surface RH (China Climate Bulletin, 2014); whilst turning to a homogenized RH dataset instead, RH is believed to remain essentially unchanged since 1961 as concluded in the post-2014 Bulletins. Zhu et al. (2015) detected break points in the RH series observed in more than 68% of stations across China, mainly occurring during 2000–2005, coincident with the massive automation of the observing system. Without diving into the history of instrument replacement, conventional homogenization procedures by automatic software are reluctant to accurately detect and properly address potential region-wide inhomogeneities from the source (e.g., homogenized HadISDH, Willett et al., 2014).

Despite improved awareness of RH inhomogeneity (Li et al., 2020c; Wang and Sun, 2021; Zhang et al., 2021), the extent to which such inhomogeneity distorts the estimate for past and future changes in humid-heat extremes has been seldomly quantified (Freychet et al., 2020). The ERA5 (Hersbach et al., 2018) and HadISD (Dunn et al., 2016), are the most widely-used datasets to analyze humid-heat extremes at impact-relevant scales (sub-daily to daily) worldwide (Li et al., 2020b; Raymond et al., 2020; Speizer et al., 2022), with in-situ observations adding values to understanding on the regional heterogeneity (Luo and Lau, 2018; Wang et al., 2019; Li et al., 2020a). If RH is inhomogeneous in these datasets, the derived quantification of long-term changes in extremes would be biased. The recent availability of two homogenized RH datasets covering the domain of mainland China (Zhu et al., 2015; Li et al., 2020a) offers us a unique chance to revisit the existing understanding of humid-heat extremes in the populous region as a showcase, calling for more emphasis on the importance of RH data quality in characterizing the multivariate extreme event worldwide.

We do this by comprehensively comparing homogenized observations with popularly-used ∼10 datasets including raw station-based observations, gridded observations and reanalysis (detailed information see Supplementary Table S1). In addition to the differences in local-to-regional changes in RH and humid-heat extremes, we will also compare the T-RH joint distribution during humid-heat extremes, i.e., the impact-relevant nature of events, in homogenized and non-homogenized datasets.

During 1979–2013, daily mean air temperature exhibits significant warming trends across China in both homogenized and raw site observations (Figures 1A, B), characteristic of highly similar spatial patterns and local-to-regional magnitudes for trends between the two (Figure 1C). This suggests a minor influence of potential temporal inhomogeneity in T on the estimate for long-term changes of the variable. By contrast, the pattern, magnitude, significance and even the sign of trends for RH differ pronouncedly between homogenized and non-homogenized station-based observational datasets. The inhomogeneity issue in relative humidity is particularly stark in southeast China (south of 30N, east of 105E, black rectangle in Figures 1D, E), leading to widespread biases in trend estimates for RH with incorrect signs there (Figure 1F). The region (black rectangle) is typical of humid-hot subtropical climates in favor of the occurrence of oppressive heat (Supplementary Figure S1; Vargas Zeppetello et al., 2022) with the two variables playing equivalently important roles to cause extremes there (Buzan and Huber, 2020).

With respect to the identified RH inhomogeneity hotspot, we further expand the comparison to gridded observations and widely-used reanalysis data at a regional scale. Large interannual to inter-decadal variability could be found in the domain-average RH series in homogenized observations (black curve, OBS-1, Figure 2A), experiencing a relatively wet regime before 2000 followed by a shift to a dry episode. This might be explained via a moist static energy balance perspective in the context of equal fractional changes in specific humidity over land and neighboring oceans (Byrne and O’Gorman, 2018) through remote moisture transport and local evapotranspiration processes (Byrne and O’Gorman, 2016; Douville et al., 2020).

Though the regional variability is generally captured by other datasets, the RH drop around 2000, when coincidentally the observational network across China shifted from manual to automatic practices (2000–2005), is exaggerated in non-homogenized datasets (except for CRA-40). This is mainly because of the systematic difference in manual and automatic observing instruments, with the lower values recorded in the latter period (e.g., Supplementary Figure S2, RH) deemed more accurate and used since then (Zhu et al., 2015). The inhomogeneity gives rise to long-term drying trends of unrealistically high magnitude and significance (Figure 2B). Specifically, the regional drying rate estimated by OBS-1 is around −.47%/decade; while the estimate based on OBS-3 is −1.43%/decade, almost three times sharper than the correct one. In particular, the HadISDH was commonly used as homogenized observations to calculate heat stress indices globally (Raymond et al., 2020; Rogers et al., 2021). Our analysis brings into question the homogeneity of HadISDH over southeastern China, suggesting that detecting region-wide inhomogeneities by automated methods (Dunn et al., 2014) without knowing detailed histories of station information (e.g., relocation and instrument replacement) might be less effective than expected. Another extensively-used dataset—ERA5 also suffers from severe inhomogeneities in RH for the region, and is therefore not suitable to characterize humid-heat extremes at sub-daily to daily scales in the region (Freychet et al., 2020). Encouragingly, incorporating into homogenized station observations (OBS-1), the CRA-40 reanalysis newly released by the China Meteorological Administration performs well in reproducing both variability and trends of RH in southeastern China. Though expanding the study region northward to 35N does not alter the comparison of trend magnitude and significance amongst datasets in any significant manner (figure omitted), it might introduce confounding factors related to differential representation of irrigation in the 30–35N latitudinal band, an intensely irrigated region (Kang and Eltahir, 2018). So, to distinguish influences of data inhomogeneity from other similar confounding factors, southeastern China to the south of 30N is a better candidate study area.

Though air temperature records were also subject to changes in observing sites, instruments and timing, the inhomogeneities from these sources seemed to pose a trivial influence on the temporal variability and long-term trends as well as their significance of domain-averaged series of air temperature (Supplementary Figure S3, OBS-1 and OBS-2 vs. OBS-3). The differences in magnitudes and significances of trends for air temperature amongst datasets are much less conspicuous than those for relative humidity. Notably, despite better performance in representing RH during 1979–2013 in CRA-40, for unknown reasons, its temperature records in 2013 are inconsistent with counterparts from any of other datasets. This year’s error in temperature records has reported back to the data developer for improvement, and the flaw also excludes CRA-40 for participation in follow-up analysis on extremes.

We next quantify extent to which the RH inhomogeneity in observation and reanalysis datasets distorts the trend estimate for humid-heat metrics, similarly starting with a site-by-site comparison (OBS-1 vs. OBS-3). We do this by regressing the differential trends for extremes in homogenized and non-homogenized observations onto the differences of RH trends in them (homogeneous trend minus non-homogeneous one). As indicated by Figures 3A–C, given the same atmospheric warming (local trends for T from OBS-1), an exaggeration of drying trend by 1% decade⁻¹ due to RH inhomogeneity translates into an underestimation of increasing trends for humid-heat extremes’ frequency by more than 1.2 days decade⁻¹, with extreme sWBGT events influenced most. The RH inhomogeneity also significantly abates the intensification rates of extreme humid-heat events, especially those resulting from the highly non-linear combination between T and RH as showcased by the HI index (Figures 3D–F; also Supplementary Figures S4A–C). Locally, an 1% decade⁻¹ drying bias caused by the RH inhomogeneity leads to an underestimation of intensification rates (absolute trends normalized by local index climatology, for the purpose of objective comparison amongst indices) for HI extremes by around .19% decade⁻¹, equivalent to a bias around .07°C decade⁻¹ there (Supplementary Figure S4).

For the domain-average series of frequency, with reference to the trend evaluated via homogenized data, all non-homogenized datasets discernably underestimate the increase, by the magnitudes ranging from .5 to 2 days decade⁻¹, basically proportional to the regional drying biases (both the X-axis and Y-axis represent the homogenized OBS-1 minus other datasets). The limited sample size (∼9 datasets) prohibits us from further evaluating the significance of such relationship. Despite being widely used, the HadISD and ERA5 datasets are amongst the worst that considerably underestimate past increases in the frequency of humid-heat extremes (Figure 4). Their bad performance arises jointly from an overestimation of the RH drying rate and an underestimation of the warming trend (Supplementary Figure S5).

Notably, despite significant negative correlation between air temperature and relative humidity at the regional scale, we did not find the cancellation or say compensation effect between greater warming and stronger RH reduction biases in the same dataset (Supplementary Figure S5) as reported in historical simulations and future projections (Fischer and Knutti, 2013), except for NCEP2 which shows the strongest negative correlation between the two variables. Rather, slightly weaker warming trends and drying bias in RH trends in tandem contribute to an underestimation of regional increases in humid-heat extremes in frequency and intensity, with the “out-of-range” RH bias (outside the grey shadings in Supplementary Figure S5) dominating in most datasets. Even for NCEP2 seemingly most faithfully reproducing extremes’ statistics, such better performance apparent arises from the wrong reason, and therefore should not be further leveraged as a reliable benchmark for model evaluation and impact quantification (Casanueva et al., 2019). Similar “right statistics for wrong reason” might exist at local scales where the T-RH anti-correlation is even stronger (e.g., <−.85).

When referring to extreme humid-heat conditions, the combination in a fashion of high temperature and high humidity comes to one’s mind first and naturally. But as a matter of fact, extreme humid-heat events represent a typical case that not all of its components are necessarily extreme, but their combination leads to an extreme impact (Zscheischler et al., 2018). Different forms of joint distributions of temperature and humidity during extreme humid-heat events are worth further sorting out, because they might cause distinct impacts on human health (Mora et al., 2017; Asseng et al., 2021).

Compared to humid-heat extremes previously identified by non-homogenized observations (OBS-3), counterparts from homogenized observations are similarly warm (mostly within the range of 25–33°C, Supplementary Figure S6) but not necessarily that moist, regardless of indicators selected (Figure 5). Quantitatively, we used to believe that the majority of extremes (around 60%–80% events across the region) occurred in the context of RH wetter than 70%; while subject to the homogenization correction, the understanding is revised into that the conventionally favorable thermodynamic environment host only half of extremes. Such a contrast in the RH constituent to extremes is also ubiquitous amongst other observation and reanalysis datasets (Supplementary Figure S7). In particular, though the bivariate distribution difference is largely determined by RH biases, the temperature difference also plays a non-trivial role in shaping the differential joint distribution structure in some datasets, e.g., ERA5 (Supplementary Figure S7E) and NCEP2 (Supplementary Figure S7G). In addition to the absolute magnitude for humid-heat extremes, the accurate mapping of typical bivariate configuration constituting the extremes, based on homogenized datasets, matters to identification of key thresholds for early warning against health-damaging events (Han et al., 2022; Vecellio et al., 2022). In light of the distinct bivariate distribution structures (Supplementary Figure S7) along with the harder-detectable inhomogeneity in the bivariate combined indices (e.g., Supplementary Figure S2, Tw), a component-wise bias correction technique might be superior to an “one-step” correction scheme with respect to the combined index, for both observations and simulations (Casanueva et al., 2019).

We present the most comprehensive comparison to date on long-term changes in relative humidity across China, based on ∼10 widely-used observational and reanalysis datasets subject to different levels of homogenization. We report an artificially sharp decline in RH around 2000–2005 in most datasets due mainly to the massive transition from manual observation practices to automatic observations with observational instruments replaced at the time. The inhomogeneity issue stands out particularly stark in warm-moist southern China, leading to unrealistically strong and significant drying trends there.

The biased drying trends in these inhomogeneous datasets lead to underestimation of changes in humid-heat extremes, with a drying bias of 1% decade⁻¹ attenuating the magnitude for increases in extremes’ frequency (intensity) by more than 1.2 days decade⁻¹ (.07% decade⁻¹) locally, regardless of metrics considered. For the region as a whole, these inhomogeneous records have underestimated domain-averaged frequency (intensity) of extremes by .47–1.95 days decade⁻¹ (−.04%∼.41% decade⁻¹). Humid-heat extremes identified by homogenized and non-homogenized datasets also differ in the T-RH joint anomalies, with most cases in the latter equivalently hot or slightly cooler yet discernably moister.

Our results call for emphasis on data quality of observed/reanalyzed relative humidity, not only in China but also in other parts of the world, and the need of re-evaluating past changes of humid-heat extremes to improve projections of deadly heat stress worldwide.