Association between Cold Spells and Mortality Risk and Burden: A Nationwide Study in China.

BACKGROUND: Few multicity studies have evaluated the association between cold spells and mortality risk and burden.
OBJECTIVES: We aimed to estimate the association between cold spells and cause-specific mortality and to evaluate the mortality burden in China.
METHODS: We conducted a time-series analysis with a nationally representative Disease Surveillance Points System database during the cool seasons spanning from 2013 to 2015 in 272 Chinese cities. We used 12 cold-spell definitions and overdispersed generalized additive models with distributed lag models to estimate the city-specific cumulative association of cold spells over lags of 0-28 d. We controlled for the nonlinear and lagged effects of cold temperature over 0-28 d to evaluate the added effect estimates of cold spell. We also quantified the nationwide mortality burden and pooled the estimated association at national and different climatic levels with meta-regression models.
RESULTS: For the cold-spell definition of daily mean temperatures of ≤5th percentile of city-specific daily mean temperature and duration of ≥4 consecutive d, the relative risks (i.e., risk ratios) associated with cold spells were 1.39 [95% confidence interval (CI): 1.15, 1.69] for non-accidental mortality, 1.66 (95% CI: 1.20, 2.31) for coronary heart disease mortality, 1.49 (95% CI: 1.12, 1.97) for stroke mortality, and 1.26 (95% CI: 0.85, 1.87) for chronic obstructive pulmonary disease mortality. Cold spells showed a maximal lagged association of 28 d with the risks peaked at 10-15 d. A statistically significant attributable fraction (AF) of non-accidental mortality [2.10% (95% CI: 0.94%, 3.04%)] was estimated. The risks were higher in the temperate continental and the temperate monsoon zones than in the subtropical monsoon zone. The elderly population was especially vulnerable to cold spells. DISCUSSION: Our study provides evidence for the significant relative risks of non-accidental, cardiovascular, and respiratory mortality associated with cold spells. The findings on vulnerable populations and differential risks in different climatic zones may help establish region-specific forecasting systems against the hazardous impact of cold spells. https://doi.org/10.1289/EHP9284.

Introduction

Global climate change has been associated with the increased occurrence of extreme weather events (McMichael et al. 2006) and has impacted human health throughout the life cycle from infancy and adolescence to adulthood and old age (reviewed by Watts et al. 2019). Epidemiological studies have provided evidence for an association between low temperature and increased mortality (Gasparrini et al. 2015) or morbidity (Ye et al. 2012). The latest Global Burden of Disease report added low temperature as a new risk factor for premature deaths, which ranked as one of the top-10 risk factors among the oldest age group (GBD 2019 Risk Factors Collaborators. 2020). Most of the temperature-related mortality, such as cardiovascular disease (CVD), was attributable to low temperature (Chen et al. 2018b; Gasparrini et al. 2015; Yang et al. 2015). Low temperature may trigger CVD by increasing platelet counts, red cells, blood viscosity, plasma cholesterol, fibrinogen, and blood pressure (Huynen et al. 2001; Xie et al. 2013). Cold spells are extreme cold weather events that last for a few days in cool seasons and have been associated with increased mortality rates throughout the world (Ryti et al. 2016). Previous studies indicated that the impact of cold spells differs from just low temperature (Lee et al. 2018; Sartini et al. 2016). Thus, the association between cold spells and cause-specific mortality risk deserves further investigation.

Frequent cold spells with high intensity (i.e., the Siberian cold current) occur during cool seasons in China (Song and Wu 2019). During January and February 2008, the majority of China suffered a severe continuous cold spell lasting for nearly a month (10 January–6 February). It caused substantial excess mortality (Xie et al. 2013) and tremendous economic losses exceeding USD (Zhou et al. 2011). Epidemiological evidence suggested that cold spells may significantly increase the risks of CVD, coronary heart disease (CHD), and stroke (Gao et al. 2019; Guo et al. 2013b; Sartini et al. 2016; Yang et al. 2015). A systematic review indicated that cold spells were associated with an approximately 11% increase in CVD mortality rates worldwide (Ryti et al. 2016). Previous epidemiological studies have reported the association between cold spells and mortality of respiratory diseases and chronic obstructive pulmonary disease (COPD) (Chen et al. 2019; Wang et al. 2016). However, with nonstandardized cold-spell definitions and limited cities or regions included in these studies, the pattern of association and potential vulnerable population of cold spells are still not determined. Moreover, few multicity studies have evaluated the mortality burden of cold spells (Gao et al. 2019; Lee et al. 2018). Complicating matters further, without considering the diversity of climatic and socioeconomic characteristics, previous studies might have exaggerated or underestimated the association between cold spells and cause-specific mortality risk (Davídkovová et al. 2014; Lee et al. 2018).

We used a nationally representative database including 272 main cities in different climatic zones of China during 2013–2015 to evaluate the patterns of association between cold spells and cause-specific mortality. We estimated the added effect estimates of cold spells by controlling for the main effect of low temperature. Moreover, we calculated the attributable fractions (AFs) attributable to cold spells to estimate the mortality burden. We then identified vulnerable populations with stratified analyses of the estimated associations of cold spells in different climatic zones, sex, age groups, and other potential association modifiers.

Methods

Data Collection

The data source of our study was based on a nationally representative Disease Surveillance Points System (DSPS) database spanning from 2013 to 2015 in 272 Chinese cities, as described in previous publications (Chen et al. 2017, 2018a; Yin et al. 2017). Cities were selected if they experienced an average of more than three non-accidental deaths occurring per day according to the DSPS death registry in China (Chen et al. 2018b), which is administered by the Chinese Center for Disease Control and Prevention. The DSPS covered 605 surveillance points of 31 provinces in China (Liu et al. 2016). Based on the classification of climatic types in China, the 272 cities can be divided into five climatic zones: 98 cities in the temperate monsoon zone, 27 cities in the temperate continental zone, 141 cities in the subtropical monsoon zone, 3 cities in the tropical monsoon zone, and 3 cities in the alpine zone. The locations of these cities and climatic zone divisions are displayed in Figure S1.

Daily death records were extracted from the DSPS. We assessed a range of causes of death based on the sole primary diagnosis coded by the International Classification of Diseases (ICD), 10th Revision (ICD-10; WHO 2016), including non-accidental mortality (ICD-10 codes A00–R99), overall CVD (ICD-10 codes I00–I99), CHD (ICD-10 codes I20–I25), stroke (ICD-10 codes I60–I69), overall respiratory diseases (ICD-10 codes J00–J98), and COPD (ICD-10 codes J41–J44).

The daily mean temperature and mean relative humidity in each city were derived from the China Meteorological Data Sharing Service System (Homepage: http://data.cma.cn/). We obtained air pollution data from National Urban Air Quality Real-time Publishing Platform of China (Homepage: http://106.37.208.233:20035) for a sensitivity analysis. We collected the gross domestic product (GDP) per capita and urbanization rates for each city from statistic yearbooks at national or provincial levels for a meta-regression analysis. The data were reported in our previous studies (Chen et al. 2018b; Yin et al. 2018). Our study protocol was approved by the institutional review board at the School of Public Health, Fudan University, with a waiver of informed consent (no. 2014-07-0523) because all data were aggregated at the city level and no subjects were contacted.

Cold-Spell Definition

Because of the lack of a universally accepted definition of cold spells, previous studies have proposed several cold-spell definitions with different cold thresholds (absolute or relative) and durations (different number of consecutive days) (Cheng et al. 2019; Kim et al. 2018; Liang et al. 2018). Given that China is a country with multiple climatic zones and that the adaptability of residents differs, we used city-specific relative thresholds (10th, 7.5th, 5th, and 2.5th percentiles) based on the distributions of the daily mean temperature during the study period, which was consistent with the definition’s criteria of extreme weather (i.e., heat waves) in previous studies (Guo et al. 2017; Yin et al. 2018). At the same time, different durations of cold spells were also considered: (2 or more consecutive d), (3 or more consecutive d), and (4 or more consecutive d). In total, we used 12 cold-spell definitions comprising the above four thresholds and three durations (Table 1). The days during the study period were divided into cold-spell or non–cold-spell days according to the definitions to calculate the relative risk (RR; i.e., risk ratio) of cause-specific mortality in each city.

Table 1

Descriptive statistics of cold-spell definitions and number of cold spells days in 272 Chinese cities from 2013 to 2015.

Name	Cold-spell definitions	Number of cold spells per year					Total number of cold-spell days
		Mean	SD	Min	Median	Max
P10_≥2d	≤P10 with ≥2-d duration	34	2	23	34	37	27,413
P10_≥3d	≤P10 with ≥3-d duration	29	2	21	29	35	23,693
P10_≥4d	≤P10 with ≥4-d duration	25	3	13	25	33	20,501
P7.5_≥2d	≤P7.5 with ≥2-d duration	25	1	18	25	28	20,402
P7.5_≥3d	≤P7.5 with ≥3-d duration	21	2	16	21	26	17,016
P7.5_≥4d	≤P7.5 with ≥4-d duration	18	3	8	18	26	14,496
P5_≥2d	≤P5 with ≥2-d duration	16	1	12	16	19	13,086
P5_≥3d	≤P5 with ≥3-d duration	13	2	8	13	19	10,588
P5_≥4d	≤P5 with ≥4-d duration	11	3	3	11	18	8,746
P2.5_≥2d	≤P2.5 with ≥2-d duration	8	1	5	8	9	6,333
P2.5_≥3d	≤P2.5 with ≥3-d duration	6	1	3	6	9	4,937
P2.5_≥4d	≤P2.5 with ≥4-d duration	5	2	0	5	9	3,677

Note: A total of 123,216 days of 272 cities from 2013–2015 (November to March) were analyzed in this analysis. Max, maximum; Min, minimum; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution; SD, standard deviation.

Statistical Analysis

We estimated the association between cold spells and cause-specific mortality with a two-stage analytical framework. Our statistical analyses were restricted to the cool seasons (i.e., from November to March) and were conducted for each cold spell definition. In the first stage, city-specific associations were estimated. The regional or national average associations were obtained in the second stage.

First-stage analysis.

In the first-stage analysis, we compared the cold-spell days with the non–cold-spell days and examined the RRs of cause-specific mortality in each city with standard time-series regression models (Bhaskaran et al. 2013). The overdispersed generalized additive models linked with a quasi-Poisson link function were used to estimate the association between cold spells and mortality counts given that the overdispersion of the daily mortality counts distribution. We applied distributed lag models (DLMs) to capture the lagged associations of cold spells (Gasparrini 2014). Previous studies indicated that the estimated associations of cold could last up to 3 or 4 wk (Gasparrini et al. 2015; Guo et al. 2014). Thus, we a priori set the lag period to 28 d. The exposure–response association of cold spells was modeled by DLM with a linear function of cold-spell days (as a binary variable). The lag–response curve was modeled with a natural cubic spline with 4 degrees of freedom (df) and 2 internal knots (plus an intercept) placed at equally spaced values in the log scale. To estimate the independent association of cold spells from temperature, our model included a distributed lag nonlinear model (DLNM) for daily mean temperature to control for the main effect of temperature by establishing a cross-basis function of temperature with the same lag period as cold spells (28 d). The DLNM for temperature included a quadratic B spline with 2 df and 2 internal knots placed at equally spaced values in the log scale and a lag–response curve with a natural cubic spline with 3 internal knots (plus an intercept) placed at equally spaced values in the log scale (Chen et al. 2018b). Other covariates include a natural cubic B spline of calendar day with 6 df per cool season (5 months) to control for seasonality and long-term trends in mortality, a natural cubic B spline of the present-day relative humidity with 3 df, and a categorical variable of the day of week (Yang et al. 2015). We used 6 df in adjusting for calendar day because it resulted in the best model fit according to the generalized cross-validation (GCV) for non-accidental mortality when comparing multiple models with 4–7 df per cool season (Table S1). We then reduced the city-specific cold spell–mortality associations to the overall association that was cumulated throughout the lag period by reducing the number of parameters to be pooled in the meta-regression models (Gasparrini and Armstrong 2013).

Second-stage analysis.

In the second-stage analysis, we used univariate meta-regression models to pool the cumulative RRs for each definition of cold spells and each cause of death at the national or the regional (climatic zones) level. We used multivariate meta-analysis models to pool the lag–response associations, which allow the analysis of complex nonlinear and delayed associations in multicity studies (Chen et al. 2018b; Gasparrini and Armstrong 2013). The potential between-city heterogeneity was evaluated with I-square () statistics and -values from the Cochran -test in meta-regression models. The value of represents the percentages of variability in RRs attributable to the heterogeneity factors, including climatic zones, latitude, longitude, average temperature, average humidity, GDP per capita, and urbanization rates of 2014.

In addition, we conducted stratification analyses by potential individual-level association modifiers, including sex (male and female), age (5–64 and ), and climatic zones. Stratification analysis was restricted to non-accidental mortality because of to the limited numbers of daily deaths from cause-specific mortality. We tested the statistical significance of differences between the estimated associations of the strata of a potential association modifier by calculating the 95% confidence interval (CI) by the formula used in the previous study (Chen et al. 2012), where and are the estimates for the two categories, and and are their respective standard errors.

Estimation of AFs.

To estimate the AFs of mortality attributable to cold spells, we calculated the overall cumulative RRs in each city by comparing the cold-spell days with the non–cold-spell days. We used pooled RRs derived from univariate meta-regression models rather than city-specific RRs in calculating national or regional burden according to recent large-scale studies (Gasparrini et al. 2015; Lee et al. 2018). Thus, we used the overall lag-cumulative RRs corresponding to cold-spell definition to calculate the AFs and 95% CI for cause-specific mortality. The formula is displayed as follows: where indicates the nationwide AFs under different cold-spell definitions (i.e., daily AFs). is the pooled cumulative RRs obtained from the second-stage analysis. is the total number of deaths in the cold-spell days during cool seasons. The total AFs of cold spells () were calculated by dividing the total number of deaths during cool seasons () by the total number of deaths attributable to cold spells in the cold-spell days ().

To compare the AFs of mortality attributable to cold spells and cold temperature during the cool season, we estimated the AFs of cold [all temperatures below the minimum-mortality temperature (16.6°C, the 97th percentile)], moderate cold (temperatures from 0.2°C to 16.6°C), and extreme cold [temperatures between the minimum () and the 10th percentile (0.2°C)] during the cool seasons, according to our previous study based on the same data set (Chen et al. 2018b). Briefly, the overall cumulative RR was calculated by comparing each day’s temperature to the minimum-mortality temperature (16.6°C; the 97th percentile). The total counts of deaths attributable to cold temperatures were counted by summing the contributions from all the days in the series and gained the total AFs by dividing the total number of deaths by the total number of attributable deaths (Chen et al. 2018b).

Sensitivity analysis.

For the sensitivity analysis, we controlled for the lagged associations of daily 24-h average particulate particles with an aerodynamic diameter of () and maximum 8-h average ozone () concentrations by adding a DLM with a duration of lag 0, 1, and 2 d for and in the above two-stage analysis models. We conducted this sensitivity analysis among 69 cities with complete meteorological and air pollution data during the study period.

All statistical analyses were conducted in R software (version 4.0.2; R Development Core Team) with the dlnm package for the first-stage analysis and the mvmeta package for the second analysis. A two-tailed was considered statistically significant for main analysis and meta-regression analysis and a 95% CI was used to test the statistical significance of the difference between groups. Using the time-series data for the 10 regions of England and Wales publicly available in a previous study (Gasparrini and Armstrong 2013), we have provided an example of R code for the identification of cold spells, first-stage analysis, and second-stage analysis (see “Code S1 Example of code in the statistical analysis” in the Supplemental Material).

Results

Descriptive Data

The 272 cities included in this study had a daily median temperature of 7.2°C and an interquartile range of 9.3°C during the cool season in 2013–2015 (Table S2). As shown in Table 1, during the study period (123,216 d in all 272 cities), we identified between 3,677 and 27,413 cold-spell days according to different cold-spell definitions. Generally, the annual mean number of cold-spell days for each city displayed a decreasing tendency with more stringent definitions of cold spells. No cold-spell day was filtrated under the most stringent definition ( with a duration) in some cities. Thus, we did not estimate the association of this cold-spell definition. During the study period, the mean concentrations of 24-h and 8-h were () and , respectively (Table 2).

Table 2

Descriptive statistics on annual average daily deaths and air pollutant concentrations in 272 cities of Chinese during the cool season of 2013 to 2015.

Variables	Mean	SD	Min	P25	P50	P75	Max	Total number of deaths
Deaths
Non-accidental	18	18	2	7	13	22	177	2,185,744
Male	10	10	1	4	8	13	105	1,253,143
Female	8	8	1	3	6	9	73	932,601
5–64 years of age	4	4	1	2	3	5	44	532,056
>65 years of age	13	14	1	5	10	16	132	1,653,688
Cold-spell day	19	19	2	8	14	23	193	510,219
Non–cold-spell day	17	17	2	7	13	22	173	1,675,525
Temperate monsoon zone	20	16	4	8	15	27	107	884,857
Temperate continental zone	9	5	3	5	8	11	24	98,162
Subtropical monsoon zone	19	20	3	9	14	23	177	1,191,650
Tropical monsoon zone	4	1	4	4	4	4	4	5,413
Alpine zone	4	2	2	3	4	5	6	5,662
CVD	9	8	1	4	7	11	72	1,101,220
CHD	3	4	0	1	2	4	33	399,616
Stroke	4	4	0	2	3	5	36	527,450
Respiratory disease	2	3	0	1	2	3	39	301,144
COPD	2	3	0	1	1	2	34	225,777
Air pollutants (μg/m^3)
PM_2.5 (24-h)	72	27	23	52	69	87	173	—
Ozone (8-h)	55	12	25	46	55	62	93	—

Note: —, not applicable; CHD, coronary heart disease; COPD, chronic obstructive pulmonary disease; CVD, cardiovascular diseases; Max, maximum; Min, minimum; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution; , particulate matter with an aerodynamic diameter of ; SD, standard deviation.

We also summarized the annual mean daily deaths of different causes in 272 cities during the cool season from 2013 to 2015 (Table 2). Briefly, we recorded 2,185,744 non-accidental deaths and 1,101,220 deaths from CVD, 399,616 from CHD, 527,450 from stroke, 301,144 from respiratory diseases, and 225,777 from COPD. The annual average daily non-accidental death on cold spell days was and on non–cold-spell days during the cool seasons. The unadjusted crude non-accidental mortality in males was higher than that in females. Higher mortality was recorded in elderly populations (Table 2).

Regression Results

The lag pattern for the associations of cold spells varied by different definitions of cold spells (Figure 1; Figures S2–S6). As shown in Figure 1, cold spells were associated with increased non-accidental mortality risk and displayed a nonlinear lagged association. In general, the RR estimates were positive on lag 1–4 d and peaked on lag 10–15 d after cold-spell exposure. In most cases, there was no cold spell–mortality association on lag 28 d. The lag patterns for the estimated associations of cold spells on different causes of deaths were generally consistent. When a lenient temperature threshold (10th percentiles of temperature) was used, the lag pattern tendencies were different from other cold-spell definitions. In the same temperature threshold categories, the lag patterns were generally consistent by duration of cold spell. In most cases, more stringent thresholds of temperature generally resulted in more marked lag patterns and higher RR estimates.

Figure 1.

The pooled lag structure for the relative risks (RRs) of non-accidental mortality associated with different cold-spell definitions in 272 Chinese cities using overdispersed generalized additive models (GAMs) and distributed lag models (DLMs) to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature with distributed lag nonlinear models (DLNMs). Estimates (95% CIs) were generated using overdispersed GAMs and DLMs to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from DLNMs), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). Solid lines represent the mean RRs of mortality (cold-spell days vs. non–cold-spell days); shaded areas represent the 95% CIs. Note: CI, confidence interval; df, degrees of freedom; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

Figure 2 and Table S3 show the pooled cumulative associations of cold spells over lag 0 to 28 d on cause-specific mortality. In general, the association between cold spells on cause-specific mortality was positive and statistically significant. The more stringent thresholds of cold spells were associated with higher cumulative RR estimates, whereas prolonging the duration was not associated with higher cumulative RR estimates. We present the results of cold spell with the definition of as the main analysis because that definition resulted in the best model fit according to the GCV value for non-accidental mortality (Table S4). Specifically, for the cold-spell definition of , the cumulative RRs were 1.39 (95% CI: 1.15, 1.69) for non-accidental mortality, 1.50 (95% CI: 1.19, 1.90) for CVD mortality, 1.66 (95% CI: 1.20, 2.31) for CHD mortality, 1.49 (95% CI: 1.12, 1.97) for stroke mortality, 1.35 (95% CI: 0.95, 1.92) for respiratory diseases mortality, and 1.26 (95% CI: 0.85, 1.87) for COPD mortality (Figure 2; Table S3).

Figure 2.

The pooled cumulative relative risks (RRs) of cause-specific mortality associated with different cold-spell definitions in 272 Chinese cities using overdispersed generalized additive models (GAMs) and distributed lag models (DLMs) to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature with distributed lag nonlinear models (DLNMs). Estimates (95% CIs) were generated using overdispersed GAMs and DLMs to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from DLNMs), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). The corresponding numeric data are presented in Table S3. Points represent the estimated RRs of mortality (cold-spell days vs. non–cold-spell days); lines represent the 95% CIs. Note: CHD, coronary heart disease; CI, confidence interval; COPD, chronic obstructive pulmonary disease; CVD, cardiovascular disease; df, degrees of freedom; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

Figure 3 and Table S5 show the cumulative RR estimates of different cold-spell definitions and non-accidental mortality by sex and age. In most cases, the cumulative RR estimates in both males and females were positive and statistically significant. The cumulative RR estimates in females were similar or higher than those in males, whereas no statistical significance was observed for the difference between males and females (Figure 3A). Among those of age, the cumulative RR estimates were positive and statistically significant; the cumulative RR estimates were similar in the population of age and in younger populations in lenient temperature thresholds (P10, P7.5) of cold-spell definitions (Figure 3B). However, the cumulative RR estimates in the elderly population were higher and statistically significant between-strata in some stringent definitions of cold spells. In the cold-spell definition of , the RR estimate in the population of age [1.77 (95% CI: 1.43, 2.20)] was higher than the younger population [1.41 (95% CI: 1.06, 1.88)], and the between-strata difference was statistically significant [ (95% CI: , )].

Figure 3.

The pooled cumulative relative risks (RRs) of non-accidental mortality associated with cold spells among (A) sex and (B) age in 272 Chinese cities using overdispersed generalized additive models (GAMs) and distributed lag models (DLMs) to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature with distributed lag nonlinear models (DLNMs). Estimates (95% CI) were generated using overdispersed GAMs and DLMs to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from DLNMs), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). The corresponding numeric data are presented in Tables S8 and S9. Points represent the estimated RRs of mortality (cold-spell days vs. non–cold-spell days); lines represent the 95% CIs. Note: CI, confidence interval; df, degrees of freedom; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

Figure 4 and Table S6 show the cumulative RR estimates of cold spells with various definitions in different climatic zones. Because of the limited number of cities in the tropical monsoon zone (3 cities) and the alpine zone (3 cities), and several cold-spell definitions in the temperate continental climatic zone, the cumulative RR estimates were of high statistical uncertainty. Despite that, the cumulative RR estimates in the temperate monsoon and temperate continental zones were generally higher than those in the subtropical monsoon zone except for the most stringent thresholds of cold spells. Specifically, for the cold-spell definition of , the estimated RRs of non-accidental mortality were 2.14 (95% CI: 1.51, 3.02) in the temperate monsoon zone and 1.33 (95% CI: 1.02, 1.74) in the subtropical monsoon zone, respectively, and the difference in RR estimates was statistically significant [0.81 (95% CI: 0.37, 1.25)]. The lag patterns of cold spells in different climatic zones are shown in Figures S7–S9.

Figure 4.

The pooled cumulative relative risks (RRs) of non-accidental mortality associated with cold spells in in 272 Chinese cities and three main climatic zones for different cold-spell definitions using overdispersed generalized additive models (GAMs) and distributed lag models (DLMs) to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature with distributed lag nonlinear models (DLNMs). Estimates were generated using overdispersed GAMs and DLMs to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from DLNMs), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). The corresponding numeric data are presented in Table S10.The estimates in the tropical monsoon zone and the alpine zone are not presented because of high statistical uncertainty caused by the limited number of cities and mortality in these climatic zones. Points represent the estimated RRs of mortality (cold-spell days vs. non–cold-spell days); lines represent the 95% CIs. Note: CI, confidence interval; df, degrees of freedom; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

In the sensitivity analysis, the cold spell with the definition of was used in models after controlling for the associations of and concentrations in 69 cities. As shown in Table 3, the pooled cumulative RRs of cause-specific mortality associated with different cold-spell definitions were not appreciably changed by this adjustment.

Table 3

The pooled cumulative relative risks (95% confidence intervals) of cause-specific mortality associated with the cold-spell definition of in 69 Chinese cities estimated using overdispersed generalized additive models and distributed lag models.

Categories	Without adjustment	With adjustmenta	With adjustmentb	With adjustmentc
Non-accidental	1.43 (1.11, 1.83)	1.41 (1.03, 1.94)	1.38 (1.01, 1.88)	1.42 (1.05, 1.91)
CVD	1.59 (1.08, 2.34)	1.57 (1.02, 2.41)	1.55 (1.01, 2.39)	1.56 (1.05, 2.30)
CHD	1.67 (0.98, 2.87)	1.67 (0.90, 3.12)	1.48 (0.78, 2.79)	1.30 (0.65, 2.60)
Stroke	1.46 (0.92, 2.32)	1.50 (0.91, 2.47)	1.52 (0.91, 2.56)	1.43 (0.83, 2.47)
Respiratory disease	1.59 (0.94, 2.72)	1.44 (0.82, 2.52)	1.31 (0.71, 2.40)	1.25 (0.63, 2.46)
COPD	1.37 (0.75, 2.53)	1.20 (0.65, 2.22)	1.15 (0.59, 2.24)	1.19 (0.51, 2.81)

Note: Relative risk estimates were generated using overdispersed generalized additive models and distributed lag models to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from a distributed lag nonlinear model), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). CHD, coronary heart disease; COPD, chronic obstructive pulmonary disease; CVD, cardiovascular disease; df, degrees of freedom; P5, 5th percentiles of the daily mean temperature distribution; , particulate matter with an aerodynamic diameter of .

Adjustment was conducted with the lag 0 d of 24-h and 8-h ozone concentrations on the same day in 69 cities with complete meteorological and air pollution data in 2013–2015.

Adjustment was conducted with the lag 1 d of 24-h and 8-h ozone concentrations on the same day in 69 cities with complete meteorological and air pollution data in 2013–2015.

Adjustment was conducted with the lag 2 d of 24-h and 8-h ozone concentrations on the same day in 69 cities with complete meteorological and air pollution data in 2013–2015.

Table 4 shows the results of AFs of cause-specific mortality associated with different cold-spell definitions nationwide. In general, the overall AFs of CVD were higher than respiratory diseases except for the stringent definitions of cold spells (, ). For the cold-spell definition of , the AFs were 2.10% (95% CI: 0.94%, 3.04%) for non-accidental mortality, 2.56% (95% CI: 1.20%, 3.63%) for CVD mortality, 2.93% (95% CI: 1.22%, 4.17%) for CHD mortality, 2.46% (95% CI: 0.82%, 3.70%) for stroke mortality, 2.11% (95% CI: , 3.89%) for respiratory disease mortality, and 1.68% (95% CI: , 3.79%) for COPD mortality. In comparison, cold temperature contributed to the largest AF of cause-specific mortality. Specifically, the AF of cold temperature was 20.29% (95% CI: 17.59%, 20.95%) for non-accidental mortality, 20.11% (95% CI: 16.84%, 20.77%) for CVD mortality, 21.04% (95% CI: 16.39%, 21.88%) for CHD mortality, 15.30% (95% CI: 13.08%, 16.03%) for stroke mortality, 18.52% (95% CI: 14.73%, 19.14%) for respiratory disease mortality, and 19.31% (95% CI: 15.98%, 20.18%) for COPD mortality during cool seasons (Table S7).

Table 4

National average fractions of cause-specific mortality attributable to different cold-spell definitions in 272 Chinese cities from 2013 to 2015.

Categories	Name	Cold-spell definitions	Attributable fractions {% [AF (95% CI)]}a
Non-accidental	P10_≥2d	<P10 with ≥2-d duration	3.87 (0.46, 6.78)
	P10_≥3d	<P10 with ≥3-d duration	3.45 (1.11, 5.50)
	P10_≥4d	<P10 with ≥4-d duration	3.18 (1.14, 4.97)
	P7.5_≥2d	<P7.5 with ≥2-d duration	4.27 (1.98, 6.23)
	P7.5_≥3d	<P7.5 with ≥3-d duration	4.69 (2.93, 6.20)
	P7.5_≥4d	<P7.5 with ≥4-d duration	3.14 (1.62, 4.46)
	P5_≥2d	<P5 with ≥2-d duration	4.07 (2.36, 5.44)
	P5_≥3d	<P5 with ≥3-d duration	3.66 (2.46, 4.64)
	P5_≥4d	<P5 with ≥4-d duration	2.10 (0.94, 3.04)
	P2.5_≥2d	<P2.5 with ≥2-d duration	1.86 (0.64, 2.78)
	P2.5_≥3d	<P2.5 with ≥3-d duration	0.89 (−0.28, 1.76)
CVD	P10_≥2d	<P10 with ≥2-d duration	5.94 (2.38, 8.91)
	P10_≥3d	<P10 with ≥3-d duration	5.94 (3.34, 8.13)
	P10_≥4d	<P10 with ≥4-d duration	4.89 (2.59, 6.83)
	P7.5_≥2d	<P7.5 with ≥2-d duration	5.73 (3.08, 7.89)
	P7.5_≥3d	<P7.5 with ≥3-d duration	5.51 (3.42, 7.21)
	P7.5_≥4d	<P7.5 with ≥4-d duration	3.87 (2.07, 5.37)
	P5_≥2d	<P5 with ≥2-d duration	4.78 (2.78, 6.32)
	P5_≥3d	<P5 with ≥3-d duration	4.27 (2.84, 5.38)
	P5_≥4d	<P5 with ≥4-d duration	2.56 (1.20, 3.63)
	P2.5_≥2d	<P2.5 with ≥2-d duration	1.85 (0.27, 2.97)
	P2.5_≥3d	<P2.5 with ≥3-d duration	0.61 (−1.07, 1.78)
CHD	P10_≥2d	<P10 with ≥2-d duration	4.97 (−0.60, 9.25)
	P10_≥3d	<P10 with ≥3-d duration	5.66 (2.22, 8.44)
	P10_≥4d	<P10 with ≥4-d duration	5.67 (3.01, 7.84)
	P7.5_≥2d	<P7.5 with ≥2-d duration	5.37 (1.47, 8.32)
	P7.5_≥3d	<P7.5 with ≥3-d duration	5.49 (2.70, 7.63)
	P7.5_≥4d	<P7.5 with ≥4-d duration	4.73 (2.63, 6.39)
	P5_≥2d	<P5 with ≥2-d duration	5.16 (2.48, 7.02)
	P5_≥3d	<P5 with ≥3-d duration	3.71 (1.68, 5.18)
	P5_≥4d	<P5 with ≥4-d duration	2.93 (1.22, 4.17)
	P2.5_≥2d	<P2.5 with ≥2-d duration	2.06 (−0.12, 3.39)
	P2.5_≥3d	<P2.5 with ≥3-d duration	0.72 (−1.67, 2.16)
Stroke	P10_≥2d	<P10 with ≥2-d duration	3.52 (−1.26, 7.39)
	P10_≥3d	<P10 with ≥3-d duration	4.38 (1.10, 7.08)
	P10_≥4d	<P10with ≥4-d duration	3.90 (1.27, 6.10)
	P7.5_≥2d	<P7.5 with ≥2-d duration	6.15 (3.13, 8.54)
	P7.5_≥3d	<P7.5 with ≥3-d duration	6.07 (4.00, 7.74)
	P7.5_≥4d	<P7.5 with ≥4-d duration	3.91 (1.97, 5.49)
	P5_≥2d	<P5 with ≥2-d duration	5.20 (3.03, 6.79)
	P5_≥3d	<P5 with ≥3-d duration	4.49 (2.95, 5.64)
	P5_≥4d	<P5 with ≥4-d duration	2.46 (0.82, 3.70)
	P2.5_≥2d	<P2.5 with ≥2-d duration	1.85 (−0.23, 3.18)
	P2.5_≥3d	<P2.5 with ≥3-d duration	0.82 (−1.23, 2.12)
Respiratory disease	P10_≥2d	<P10 with ≥2-d duration	0.83 (−8.90, 7.75)
	P10_≥3d	<P10 with ≥3-d duration	—
	P10_≥4d	<P10 with ≥4-day duration	—
	P7.5_≥2d	<P7.5 with ≥2-d duration	2.95 (−2.63, 7.06)
	P7.5_≥3d	<P7.5 with ≥3-d duration	1.72 (−2.83, 5.15)
	P7.5_≥4d	<P7.5 with ≥4-d duration	0.65 (−3.28, 3.66)
	P5_≥2d	<P5 with ≥2-d duration	3.66 (−0.08, 6.23)
	P5_≥3d	<P5 with ≥3-d duration	2.46 (−0.61, 4.62)
	P5_≥4d	<P5 with ≥4-d duration	2.11 (−0.42, 3.89)
	P2.5_≥2d	<P2.5 with ≥2-d duration	2.71 (0.06, 4.14)
	P2.5_≥3d	<P2.5 with ≥3-d duration	2.12 (0.49, 3.12)
COPD	P10_≥2d	<P10 with ≥2-d duration	0.84 (−9.46, 8.02)
	P10_≥3d	<P10 with ≥3-d duration	—
	P10_≥4d	<P10 with ≥4-d duration	—
	P7.5_≥2d	<P7.5 with ≥2-d duration	3.18 (−3.07, 7.62)
	P7.5_≥3d	<P7.5 with ≥3-d duration	1.26 (−4.37, 5.30)
	P7.5_≥4d	<P7.5 with ≥4-d duration	0.54 (−4.08, 3.95)
	P5_≥2d	<P5 with ≥2-d duration	3.25 (−1.20, 6.18)
	P5_≥3d	<P5 with ≥3-d duration	1.94 (−1.96, 4.52)
	P5_≥4d	<P5 with ≥4-d duration	1.68 (−1.45, 3.79)
	P2.5_≥2d	<P2.5 with ≥2-d duration	2.94 (0.26, 4.33)
	P2.5_≥3d	<P2.5 with ≥3-d duration	2.19 (0.09, 3.35)

Note: Estimates were generated using overdispersed generalized additive models and distributed lag models to estimate the associations of cold spell over lags of 0–28 d after controlling for the main effect of daily temperature (cross-basis function for temperature lagged for 0–28 d from a distributed lag nonlinear model), adjusted for calendar day (natural cubic spline with 6 df), day of the week, and humidity (lag 0, natural smooth function, 3 df). —, not applicable; AF, attributable fraction; CHD, coronary heart disease; CI, confidence interval; COPD, chronic obstructive pulmonary disease; CVD, cardiovascular diseases; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

The attributable fractions of cold spells in cool seasons.

Unable to calculate the attributable fractions due to relative risk estimates below one.

Table 5 shows the results of heterogeneity tests in the cumulative RRs of non-accidental mortality associated with four representative cold-spell definitions using meta-regression models. There was moderate heterogeneity for the associations between cold spells and non-accidental mortality risks with values ranging from 42.29% to 49.93%. We found some statistically significant geographical and climatic modifiers (i.e., latitude, average temperature, and humidity) that contributed to appreciable proportions of the non-accidental heterogeneity () in the meta-regression models (Table 5). The results from meta-regression analyses showed stronger associations with non-accidental and cardiovascular mortality in cities with higher latitude, lower daily average temperatures, and lower average humidity.

Table 5

The heterogeneity and potential heterogeneity factors in the cumulative relative risks of non-accidental mortality associated with four representative cold-spell definitions in meta-regression models in 272 Chinese cities from 2013 to 2015.

Cold-spell definitions	Predictor	p-Value for Q-test	I^2 (%)	Estimate for predictors	p-Value for predictors
P10_≥3d	None	<0.001	48.40	—	—
	Climatic zone	<0.001	45.39	—	0.68
	Latitude	<0.001	47.76	0.021	0.044
	Longitude	<0.001	48.47	−0.007	0.47
	Average temperature	<0.001	47.63	−0.019	0.048
	Average humidity	<0.001	47.57	−0.015	0.023
	GDP per capita	<0.001	48.43	<0.001	0.23
	Urbanization rates	<0.001	48.59	0.004	0.42
	All predictors	<0.001	45.24	—	0.55
P7.5_≥4d	None	<0.001	49.74	—	—
	Climatic zone	<0.001	47.12	—	0.056
	Latitude	<0.001	48.68	0.031	0.008
	Longitude	<0.001	49.40	0.018	0.088
	Average temperature	<0.001	49.07	−0.020	0.053
	Average humidity	<0.001	47.76	−0.024	0.001
	GDP per capita	<0.001	49.87	<0.001	0.45
	Urbanization rates	<0.001	49.93	0.003	0.60
	All predictors	<0.001	47.42	—	0.22
P5_≥4d	None	<0.001	48.15	—	—
	Climatic zone	<0.001	47.99	—	0.19
	Latitude	<0.001	48.31	−0.006	0.72
	Longitude	<0.001	48.31	−0.006	0.65
	Average temperature	<0.001	48.32	0.005	0.75
	Average humidity	<0.001	47.93	−0.016	0.082
	GDP per capita	<0.001	48.31	<0.001	0.62
	Urbanization rates	<0.001	48.27	0.007	0.35
	All predictors	<0.001	46.88	—	0.36
P2.5_≥3d	None	<0.001	44.22	—	—
	Climatic zone	<0.001	42.29	—	0.33
	Latitude	<0.001	44.16	−0.041	0.098
	Longitude	<0.001	44.43	0.001	0.96
	Average temperature	<0.001	44.13	0.039	0.083
	Average humidity	<0.001	44.20	−0.007	0.64
	GDP per capita	<0.001	44.21	<0.001	0.28
	Urbanization rates	<0.001	44.38	−0.012	0.30
	All predictors	<0.001	43.85	—	0.48

Note: The results of heterogeneity and potential heterogeneity factors in the associations of non-accidental mortality with four representative cold-spell definitions were generated using meta-regression models. represents the percentage of variability in RRs attributable to the heterogeneity factors; -values were generated by meta-regression models, and two-tailed was considered statistically significant. —, not applicable; GDP, gross domestic product; P2.5, 2.5th percentiles of the daily mean temperature distribution; P5, 5th percentiles of the daily mean temperature distribution; P7.5, 7.5th percentiles of the daily mean temperature distribution; P10, 10th percentiles of the daily mean temperature distribution.

Discussion

Our study investigated the association between cold spells and cause-specific mortality, and estimated the mortality burden using a nationwide meteorological and cause-specific mortality database in China. Our results showed significant associations between cold spells of various definitions and cause-specific mortality risk with a lagged association peaked at 10–15 d and attenuated at 28 d. The stratified analyses showed that the elderly population was more vulnerable to cold spells, whereas no statistically significant difference was observed by sex. Our observed differentiated associations of cold spells in different climatic zones highlighted the necessity to establish the region-specific forecasting systems against the impact of cold spells.

In our study, the impact of cold spells on non-accidental and cause-specific mortality showed a notable lag pattern up to 28 d, indicating a relatively durable impact of cold spells. Consistent with our results, previous studies also showed that the estimated association of cold spells may last up to 3 or 4 wk (Chen et al. 2019; Gao et al. 2019; Wang et al. 2016). Most definitions of cold spells resulted in positive and statistically significant RRs for non-accidental mortality. The associations between cold spells and increased mortality have been well documented by previous studies in Asia (Chen et al. 2019; Lee et al. 2018; Wang et al. 2016) and in other regions around the world (Huynen et al. 2001; Ryti et al. 2016). Similar to the associations of heat wave and mortality risk in our previous study (Yin et al. 2018), we found that, for the same temperature threshold of cold spells, prolonging duration did not correspond to larger RR estimates, which could be partly due to the long lag period of cold spells. Thus, identification of the lag pattern may have important public health significance in extending the preventive measures of cold spell to longer duration.

Consistent with previous epidemiological studies (Gao et al. 2019; Sartini et al. 2016; Xie et al. 2013), we found that cold spells were associated with increased mortality risks of CVD, CHD, and stroke. Similar associations between cold spells and CVD were observed in Guangdong Province (Xie et al. 2013) and Shanghai (Ma et al. 2013) of China. Several epidemiological studies have also documented the associations between increased CHD mortality risk and cold spells (Davídkovová et al. 2014; Revich and Shaposhnikov 2008; Wolf et al. 2009). It was also reported that cold spell may be a risk factor for stroke (Gao et al. 2019; Ma et al. 2013), which coincided with our results. For a cold spell with definition of , the RR for CVD in our study was comparable with the RR [1.69 (95% CI: 1.48, 1.89)] estimated in a previous epidemiological study of 31 cities in China (Chen et al. 2019). The RR of CVD mortality due to cold spells in our study was higher than the RR [1.11 (95% CI: 1.03, 1.19)] estimated in a previous systematic review and meta-analysis (Ryti et al. 2016). The disparity could be explained by the different methods of controlling for daily temperature and different climatic zones among these studies.

We also estimated the added and independent association of cold spells, controlling for the main effect of low temperature, which is consistent with a previous study (Sartini et al. 2016). The association of cold weather and mortality is largely attributable to the direct effects of cold (Huynen et al. 2001), and the association of cold on the cardiovascular system is often due to increased cardiovascular risks in relation to changes in the autonomic nervous system, blood pressure, blood coagulation system, inflammatory response, and oxidative stress (Cai et al. 2016; Ryti et al. 2016). As an added and independent association from cold temperature (Sartini et al. 2016), the possible mechanism for the estimated association of cold spells could be the higher morbidity of respiratory infections, such as during influenza epidemics (Chen et al. 2019; Huynen et al. 2001). Moreover, contractility of veins, plasma fibrinogen, blood viscosity, and blood pressure may further trigger thrombosis and induce sudden mortality (Davídkovová et al. 2014). In addition, cold spells may also contribute to excess mortality because health services are less available during cold spells (Wang et al. 2016).

We estimated statistically significant RRs of respiratory disease and COPD associated with cold spells in the stringent definitions of cold spells (i.e., P5, P2.5), which is consistent with a previous epidemiological study (Xie et al. 2013). The RR for respiratory disease estimated in our study was comparable with the RR for respiratory diseases [1.21 (95% CI: 0.97, 1.51)] in a previous systematic review and meta-analysis from several Asian and European cities (Ryti et al. 2016). Conversely, previous studies in Shanghai (Liang et al. 2018) and the Netherlands (Huynen et al. 2001) did not observe a statistically significant association between cold spells and respiratory mortality. The disparity could be explained by the different climatic zones and cold-spell definitions among these studies. The estimated association of cold spells on the respiratory system might be due to increased respiratory infections, bronchoconstriction, and susceptibility of the immune system during cold days (Huynen et al. 2001; Liu et al. 2015; Xie et al. 2013). The independent respiratory associations of cold spells may be also explained by the reinforced cold effects and reduced or limited availability of medical services during cold spells.

Few studies have evaluated the mortality burden of cold spells. In our study, we estimated small but statistically significant AFs of mortality attributable to cold spells after controlling for the lagged and nonlinear effects of daily temperature. The results of our study were comparable with the AFs estimated in Hefei City [ischemic stroke mortality, (95% CI: 0.47%, 4.18%)] (Gao et al. 2019). For the same cold-spell definition of , the AF of non-accidental mortality [406% (95% CI: 2.37%, 5.44%)] estimated in our study was higher than the AF in Korea and Japan [1.44% (95% CI not provided)] (Lee et al. 2018). The disparity of AFs attributable to cold spells in different regions may be due to the different methods of controlling for daily temperature. In addition, the climatic diversity and cold-spell definitions could also contribute to the difference among studies. Our results suggest that the mortality burden of cold spells in China cannot be neglected and is worthy of attention. Besides, it has been indicated that the risk of cold spells in mid-latitude continents might not decrease in the future (Cohen et al. 2018). Considering the aging population of China (Zeng and Hesketh 2016), the potential mortality burden of cold spells should be taken into account.

The stratified analysis showed that the RRs of cold spells in females were similar to or higher than those in males in certain cold-spell definitions, but this difference was not statistically significant. Several studies also observed no significant difference in the cold spell–mortality association across sex (Chen et al. 2019; Ma et al. 2013). Conversely, a study conducted in 66 communities of China reported lower cumulative RR estimates for females than males (Wang et al. 2016). The disparity could be explained by the differences in sample size and cold-spell definitions, as well as by whether there was an adjustment of daily temperature. Consistent with previous studies (Huynen et al. 2001; Xie et al. 2013), our study indicates that the elderly population may be more susceptible to cold spells.

In this nationwide analysis, our meta-regression analysis revealed moderate heterogeneity in city-specific risk estimates. We found some significant geographical and climatic modifiers (i.e., latitude, average temperature, and humidity) contributing to appreciable proportions of the heterogeneity. Different demographical structure, socioeconomical levels, and long-term climatic adaptability could also contribute to the heterogeneity (Zhang et al. 2014). It is noteworthy that annual average temperature was a significant heterogeneity factor in the risk estimates of cold spell with a lenient threshold (). This could be interpretable with the relatively weak estimated association of cold spells in the subtropical monsoon zone, which is characterized by warm temperature. Studies also indicated that the estimated association of cold spells could be modified by different climatic types (Guo et al. 2013b; Ryti et al. 2016). Correspondingly, the RR estimates of cold spells in high-latitude and low-humidity areas (i.e., the temperate monsoon zone and the temperate continental zone) were generally higher than those in lower-latitude and high-humidity areas (i.e., the subtropical monsoon zone) in our study. Considering the diverse features in different regions of China, it is necessary to establish region-specific health forecasting systems for cold spells.

This study has several notable strengths. First, we applied the largest database of long-term, multiple cities, and cross-regionals collected from 2013 to 2015 of 272 cities in China. Second, we used the DLM methodology and various cold-spell definitions to estimate the cause-specific mortality risk and mortality burden of cold spells. Third, this investigation also provided robust evidence of the risks in different climatic zones, vulnerable populations, and socioeconomic characteristics. All the factors mentioned above contributed to the reliable external representativeness and robustness of our results and conclusion.

Potential limitations of this study should be noted. First, as done in most previous studies, we evaluated the exposure to cold spells relying on fixed-site measurements of temperature, which would lead to inevitable exposure measurement errors. However, a previous study had indicated that the exposure misclassification was likely to be randomly distributed and to typically cause an underestimation of the risk (Guo et al. 2013a). Second, this analysis is inherently an ecological study and individual-level risk factors cannot be controlled. Third, the statistical power of this study could be attenuated, especially for stratified and region-specific analyses, because of the relatively short study period (2013–2015) at the city level, but this may not substantially affect our results at the national level.

Conclusion

Our study provides evidence for the association between increased risk of non-accidental, cardiovascular, and respiratory mortality and cold spells that were independent of low temperature and accounted for a small but statistically significant mortality burden. The potential risks could last for up to 28 d, with a peak at . The RRs varied appreciably by cold-spell definition and a range of demographical, geographical, and climatic factors. Our findings highlight the importance of establishing region-specific forecasting systems for protecting vulnerable populations from the hazardous impact of cold spells.

Click here for additional data file.

Click here for additional data file.