Changing rapid weather variability increases influenza epidemic risk in a warming climate

Daily maximum surface air temperatures of Europe and the continental United States for the period January 1, 1997 to February 28, 2018 are derived from the Global Historical Climatology Network (ftp://ftp.ncdc.noaa.gov/pub/data/ghcn/daily/), covering a total number of 7729 days. The data are all the meteorological station-observed daily maximum temperatures. Many individual stations contain non-observed days. In this study, only the daily temperature time series of individual stations that satisfy the following conditions are analyzed: (1) the station-observed temperature time series has no more than 80 missing days for any individual winter half year (from August 1 to February 28); and (2) fewer than 750 days of missing data over the study period (1997–2018). There are 559 stations satisfying these two conditions in the United States and 547 in Europe. The available temporal domain of Italy is from 1999 to 2018. The spatial interpolation using the Matlab function griddata is then applied to fill in missing data for any individual days for all meteorological stations. The homogenized data is then averaged to obtain state-wise (United States) or country-wise (Europe) averaged daily maximum temperatures. The Global Historical Climatology Network contains relatively fewer numbers of stations over China. Therefore, the station-observed daily maximum temperature over China is from the quality-controlled [21] observations assembled by the National Meteorological Information Center of the China Meteorological Administration (http://data.cma.cn/). There are 654 stations over China that satisfy the above two conditions.

Daily near-surface mean absolute humidity is from reanalysis provided by the National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR reanalysis) [22]. This dataset is a spatially gridded one that has a fixed zonal resolution of 1.875 degrees of longitude and a varying Gaussian-shaped meridional resolution with its average close to 1.875 degrees of latitude. The data covers the whole global domain starting from January 1, 1948. In this study, we analyze data covering the same temporal span, from January 1, 1997 to February 28, 2018, to obtain state-wise (United States) or country-wise (Europe and China) averaged absolute humidity fields. The reanalysis can be downloaded from: https://esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.html.

The simulated and projected daily surface maximum temperatures are from the CMIP5 [23] for the periods Jan. 1950–Dec. 2005 and Jan. 2006–Dec. 2100, respectively. The projected future climate data selected are from two future emission scenarios: the RCP4.5 and RCP8.5. Seven widely acclaimed Earth system models from different countries are selected, with their names, institutes, and horizontal resolution listed in supplementary table S1 (available online at stacks.iop.org/ERL/15/044004/mmedia). The simulated and projected data can be downloaded from: https://esgf-node.llnl.gov/projects/cmip5/.

Due to the different settings in the influenza datasets (supplementary table S2), the analyzed data from the United States of America (USA) are the ILI, defined as the percentage of patients with ILI among all patients. For European countries and mainland China, the analyzed datasets contain the weekly percentage of confirmed influenza patients among all tested patient samples, which is referred to as influenza morbidity (IM) in this study. Due to the varying spatiotemporal resolutions of different datasets, interpolation or summation methods are used to obtain the same spatiotemporal resolutions if necessary.

Weekly IM data over Europe and China are obtained from the Influenza Laboratory Surveillance Information of the World Health Organization. Weekly IM data from34 countries are analyzed in this study, with source laboratories and available temporal domains listed in supplementary table S2. The data can be downloaded from: http://apps.who.int/flumart/Default?ReportNo=14.

The ILI data are from the US Center for Disease Control and Prevention. The dataset has weekly temporal resolution and includes three different subsets: (1) the state-wise reported percentage of ILI patients covering the period from the 40th week of 2010 to the present; (2) the standard federal region-wise data covering the period from the 40th week of 1997 to the present; and (3) the whole United States averaged data covering the period from the 40th week of 1997 to the present. The data can be downloaded from: https://gis.cdc.gov/grasp/fluview/fluportaldashboard.html. The reason for using ILI data rather than IM data over the United States is due to the inconsistency in IM data. According to the US Center for Disease Control and Prevention, the methods for collecting IM data changed after 2015–2016, while those for collecting ILI data remained the same for 1997 to 2018.

It has been proposed that sudden large changes in temperature can impair the human immune system and trigger immune evasion [24–29]. The possible mechanism is that the human thermoregulation of immune defense is less adjustable to the sudden large change in temperature [26, 27] and less resistant to various diseases [24, 25]. Recently, studies showed that the sudden large change in temperature tends to cause high respiratory mortality [24–27, 30] and to impact influenza seasonality [31]. It is, in fact, based on this pathological mechanism that we suspect that frequent fluctuating weather, as an additional factor, may play a significant role in influenza epidemics.

To quantify the sudden large change in temperature, here we introduce a climatic quantity called rapid weather variability (RWV) for any given temporal location: the total number of consecutive days with surface air temperature differences [32] larger than 3 K over a three-week period ending at that temporal location (see supplementary table S1). The procedure for quantifying RWV is shown in supplementary figure S1; the number of absolute differences of consecutive-day surface air temperature (figure S1(b)) is first calculated from the raw daily surface air temperature time series (figure S1(a)), and the number of days with differences larger than 3 K over a three-week sliding window is obtained (figure S1(c)). The selection of a temporal window length of three weeks allows us to retain sufficient variability of the RWV and reduce the effect of randomness. It is important to note that the threshold of 3 K is not an arbitrary selection; rather, it is based on the sensitivity of the human immune system to temperature variability [26]. Both positive and negative changes in temperature of more than 3 K would remarkably increase the respiratory mortality [26]. It is noted that the key results presented later in this paper are not sensitive to the selection of window size for defining RWV as long as it is between two weeks to four weeks (see supplementary text and supplementary figures S2–S3).

In this study, we explore the impact of RWV on ILI/IM. Since the former cannot be changed by the latter, we anticipate an either simultaneous or delayed response of ILI/IM to RWV. To characterize this relation, lagged correlations of various delays are calculated for the weekly anomalous weather variable V_i,j and the weekly anomalous value N_i,j of ILI/IM of the same region, where the subscripts i and j represent the ordered year and ordered week of a year, respectively. The lagged correlation between V_i,j and N_i,j is defined as

$$\begin{eqnarray*}&&{r}_{j}\left(\tau \right)=\displaystyle \frac{\displaystyle {\sum }_{i={i}_{s}}^{{i}_{e}}{V}_{i,j}{N}_{i,j+\tau }}{\sqrt{\displaystyle {\sum }_{i={i}_{s}}^{{i}_{e}}{V}_{i,j}{V}_{i,j}}\sqrt{\displaystyle {\sum }_{i={i}_{s}}^{{i}_{e}}{N}_{i,j+\tau }{N}_{i,j+\tau }}}\end{eqnarray*}$$

where i_s is the starting influenza season, and i_e is the ending influenza season. In this study, i_s corresponds to the 1997–1998 influenza season and i_e corresponds to the 2017–2018 influenza season, but with the 2009–2010 influenza season excluded as the ILI/IM in that season was dominated by swine influenza that has abnormal spreading dynamics [33–35]. In addition, the Monte Carlo test was used to examine the statistical significance of lagged temporal correlation. The statistical significances of different countries for the correlation coefficients of ±0.4, ±0.45, ±0.5, and ±0.55 are listed in supplementary table S3.

Figure 1 presents various relations between autumn RWV and ILI/IM. For the four northern mid-latitude countries/regions, i.e. the United States, mainland China, Italy, and France, which have relatively longer ILI/IM data, a common feature that emerged is that the autumn RWV appears to have a long-lasting effect on the influenza epidemic strength in the subsequent months (coherent reddish blocks in the correlation maps displayed in the left column panels of figure 1; statistical significances are shown in supplementary table S3), implying that the RWV has preconditioned the occurrence of the influenza epidemic, although the influenza peaking times were different for all these regions. The case in figure 1(b) in which the January RWV is correlated with influenza mortality in the previous November was largely caused by a few individual influenza seasons, such as the 2010–2011 winter, 2013–2014 winter, and 2015–2016 winter (supplementary figure S8(b)). With a small sample, such a spurious relation can exist, although it is unreal. Another feature common to all four regions is that the peak strength of ILI/IM, P increases with the strength of the RWV in autumn. The normalized change rates, $\left(P-\bar{P}\right)/\bar{P},$ corresponding to the RWV value changed by 1, are 23%, 27%, 14%, and 26% for the United States, mainland China, Italy, and France, respectively (supplementary figure S13). These values indicate that the ILI/IM is highly impacted by RWV, at least from a statistical perspective. It is also noted that the statistical relation between later autumn RWV and the following winter ILI is most robust in populations of ages less than 5 and greater than 65 for the United States (see supplementary text and supplementary figure S4).

Zoom In Zoom Out Reset image size

To further confirm the lasting effect of autumn RWV on the later winter ILI/IM peaking strength, figure 2 presents the spatial patterns of RWV and ILI/IM of the 2017–2018 influenza season for the United States and a large portion of Europe. The influenza epidemic season of 2017–2018 was one of the severest influenza seasons in the United States and Europe [11, 12], causing up to 4064 mortalities a week in the United States alone. The spatial patterns of RWV in mid-October to mid-November of 2017 for the United States and in September of 2017 for Europe share similar spatial structures of their corresponding peak ILI/IM in January and February of 2018, with spatial correlations between RWV in autumn and ILI/IM of 0.32 for the United States and 0.33 for Europe, and both exceeding a 95 percent confidence level against a null hypothesis of random spatial distribution. However, the simultaneous spatial correlations between RWV and peak ILI/IM are much smaller, with values of −0.06 in the United States and −0.23 in Europe, and are not statistically significant.

Zoom In Zoom Out Reset image size

Results from the above analysis suggest a mechanism of ILI/IM temporal evolution: in later autumn, the intensified RWV contributes to the increase in influenza patients. When the mass of patients reaches a critical level in a densely populated region, the direct contacts between influenza patients and healthy persons increase and the rate of persons being infected with the influenza virus reaches a level greater than the rate of influenza recovery, leading to a fast increase in influenza patients and a severe influenza season. The 2017–2018 influenza season of the United States appears to confirm the hidden operation of such a mechanism, although it may not be the most dominant one. Since an influenza epidemic's peak season is mostly in the second half of winter and RWV in autumn contributes significantly to building critical patient levels in densely populated regions, the strength of RWV during autumn may serve as a valuable predictor of the severity of the consequent influenza season, thus facilitating earlier preparation and prevention.

The statistical relationships between winter-time temperature/humidity and influenza epidemics in many highly populated regions of northern mid-latitudes over the past two decades was also reexamined. However, the identified relations between winter temperature/humidity and influenza epidemics are less robust (see supplementary text and supplementary figures S5–S10).

Over the past century, the Earth's climate has been changing at an unprecedented pace [36, 37], especially in the highly populated northern mid-latitudes. This climate change is not limited only to surface temperature; rather, it also includes changes in other climate variables, such as the mid-latitude synoptic variability [23, 38]. Since the change in RWV can cause large changes in the ILI/IM, as revealed above, understanding the spatial patterns of future RWV can help us determine the severity of future ILI/IM threats.

The above identified relation between RWV and ILI/IM provides a chance to estimate the changing risk of the influenza epidemic in a warming climate. Figure 3 presents the triweekly RWV changes for two three-decade spans, 2020–2049 and 2070–2099, based on the diagnoses of seven climate system model (supplementary table S1) outputs from the Representative Concentration Pathway (RCP) 8.5 emission scenario of the Coupled Model Intercomparison Project Phase 5 (CMIP5) [39]. CMIP5 was used to project future climates in the Intergovernmental Panel on Climate Change's fifth assessment report. The displayed RWV changes are the anomalies of the means of triweekly autumn RWV of the above selected two temporal spans with respect to the historical three-decade mean RWV for the temporal span 1970–1999. While the detailed spatial patterns of RWV changes in these seven models are not exactly the same, a common feature emerges: many regions of northern mid-latitudes will have an RWV anomaly increase of more than 0.5 over the next three decades and more than 1 in the last three decades of the 21st century. Most of these models projected that highly populated Europe will have the largest RWV increase, with triweekly RWV increases larger than 1 over the next three decades and more than 2 in the last three decades of the 21st century. This large increase in RWV projects that some regions of Europe will have the risk of IM increase by more than 50% if solely following the statistically robust RWV–IM relation presented in figure 1. However, this number may be an overestimation since RWV is not the only cause of IM. All model-projected RWV changes in the United States show an increase only in its highly populated northeast region with a relatively smaller value. The models also projected that the large majority of China will also have RWV increases by up to 1 by the end of the 21st century, implying that IM will increase by more than 20%. In summary, the increasing RWV is likely to increase the risk of IM in future, although a reliable quantification is out of our reach.

Zoom In Zoom Out Reset image size

The validity of the above projection of the risk of ILI/IM is highly dependent on the accuracy of the projection of RWV. Is this projected RWV physically justified? The historical data shows that the noticeable land warming after the industrial revolution started first in the land neighboring the Arctic and the subtropical northern hemisphere [40]. These two bands intensified and generated an anomalous temperature gradient that led to the increased anomalous heat transport from polar/subpolar regions and subtropical regions to northern mid-latitude regions in the autumn season. By the end of the 21st century, the maximum warming was at about 50 °N (figure 4), leading to the increased meridional temperature gradient north to that latitude [40]. The increased meridional temperature gradient provided a more favorable environment for the synoptic system to develop [41–43], explaining the intensification of RWV. This maximum warming latitude is still shifting southward as the globe continues to warm, see figure 3 of Ji et al [40]. The climate projection based on the above-mentioned model outputs show that this zone of increased meridional temperature gradient is further south, with its southern edge located at about 45 °N. Thus, the increased RWV in climate model projection is consistent with previous understandings and can be anticipated. It is noted that, for the Asian region, the zone of the intensified meridional temperature gradient locates south to the zonally averaged one, while for the United States and Europe, this zone locates north to the zonally averaged one.

Zoom In Zoom Out Reset image size

The similar spatial pattern changes in RWV and zonally averaged temperatures are also seen in the projections under a more moderate emission scenario, RCP4.5. Supplementary figure S11 displays the triweekly autumn RWV changes in seven climate system models for a more moderate emission scenario, RCP4.5. The spatial patterns of RWV increase under RCP4.5 and resemble those of RCP8.5 (figure 3), but with reduced amplitude. Under the RCP4.5 emission scenario, many regions of northern mid-latitudes will have an RWV anomaly increase more than 0.3 over the next three decades, and more than 0.4 in the last three decades of the 21st century. Under the RCP4.5 scenario, the IM over Europe and China may increase by more than 30% and 15%, respectively, in the last three decades of the 21st century. In addition, model projections under the RCP4.5 emission scenario are characterized by maximum warming over northern mid-latitudes (supplementary figure S12).