Main

Due to the outbreak of COVID-19, China activated the First Level Public Health Emergency Response (here called a quarantine), which required local governments to carry out strict restrictions on travel1,2. The entire country was under this quarantine, which lasted for one month and was arguably unprecedented regarding its spatial coverage, duration, strictness and effectiveness for preventing the spread of COVID-19 (ref. 2). There was an observed improvement in ambient air quality during the quarantine, probably due to the limited industrial and transportation activities coupled with favourable meteorological conditions3,4,5,6,7,8. Some expected that the air quality improvement may have reduced the exposure of the population to air pollutants, such as NO2 (refs. 6,7,9). If coupled with reductions in the ambient levels of fine particulate matter with a diameter <2.5 μm (PM2.5), for which the best information is available on health impacts, the quarantine may have yielded an inadvertent health benefit during the COVID-19 pandemic. How actual population exposure changed, however, depends not only on the ambient air quality but also on the air quality indoors, and the mobility and daily activity patterns of individuals, such as the time spent in different locations10,11,12,13.

The quarantine triggered by the outbreak started on 25 January 2020, which coincided with the start of the 2020 Spring Festival. Just before the start, reportedly 125 million migrant workers had moved from urban to rural areas to reunite with their families14. Normally, they would have returned at most one month after the start of the Festival15. The nationwide returning-to-work tide, however, was frozen by travel restrictions under the quarantine16,17. Thus, an extra 9% of the Chinese population were kept in rural areas longer because of the COVID-19 outbreak, areas where household air pollution (HAP) is more severe due to the prevalent use of solid fuels (coal and biomass) for cooking13,18. Also, during that season, there was still substantial space heating in households over much of the country, which is even more likely to be done with solid fuels than is cooking19.

The question we ask is how the overall PM2.5 exposure of the Chinese population changed during the COVID-19 quarantine, taking into account the changes in indoor and outdoor concentrations, time spent indoors and outdoors, and large-scale migration patterns. Such an assessment is of interest because of the health impacts of short-term exposure to PM2.5 (refs. 20,21,22) and the reported associations between PM2.5 and the spread and severity of the COVID-19 infection23,24. Here, we use real-time migration data, during-pandemic activity survey data, national census data, advanced air quality modelling techniques, and an indoor exposure model to track the dynamic changes in the population exposure to PM2.5 across China before and during the nationwide quarantine (Methods).

Results

Overall change in PM2.5 exposure before and after the COVID-19 quarantine

We focus our comparison on two periods: P1, one month preceding the Spring Festival spanning from 25 December 2019 to 24 January 2020; and P2, one month following the start of the Spring Festival, that is, the quarantine period, spanning from 25 January 2020 to 25 February 2020 (Fig. 1). Using surveys on time-activity patterns of the Chinese population both in normal days and during the quarantine (Supplementary Table 1), population time use is parsed between indoors and outdoors (Methods). Data fusion using an ensemble deep learning method to integrate the ground-level measurements of the national monitoring network with the outputs of a chemical transport model25 (Methods) shows a decrease of 16.6 μg m−3 (15.3–17.9 μg m−3, uncertainty is expressed as 95% confidence interval (CI) throughout) in the population-weighted average of ambient (outdoor) PM2.5 concentrations between P1 (64 (57–70) μg m−3) and P2 (47 (41–54) μg m−3). In contrast, the population-weighted exposure (PWE) that considers both indoor and ambient concentrations shows an increase of 5.7 μg m−3 (1.2–11.0 μg m−3) in P2 (101 (87–117) μg m−3) compared to P1 (95 (83–109) μg m−3) (Fig. 1), suggesting important roles of other factors, including population migration and the time spent indoors, in the PWE change under the quarantine. Decomposition analysis, by changing the factors severally, attributes the changes of –10.5 (–12.8 to –8.5), 10.8 (8.2–14.2) and 5.4 (3.2–8.1) μg m−3 in PWE to the changes in ambient PM2.5, population migration and time spent indoors, respectively (Fig. 2). Note that changes in ambient PM2.5 affect indoor concentration through infiltration, which is included in our assessment. Population migration alone offsets the effect of the ambient air quality improvement on PWE (10.8 μg m−3 due to migration versus −10.5 μg m−3 due to ambient air) (Fig. 2).

Fig. 1: Daily trends of PWE in the real case and under different counterfactual scenarios during the study period.
figure 1

The dark- and light-shaded areas represent the interquartile ranges and the 95% CIs of the time series, respectively. Compared to the real case, the ‘2019 migration’ scenario assumes that there was no COVID-19 outbreak such that the migration followed the pattern of the 2019 Spring Festival (instead of 2020) and the time spent indoors was not affected by the quarantine. The ‘no migration’ scenario assumes no COVID-19 outbreak and no Spring Festival migration. Weather conditions and ambient PM2.5 levels in the two counterfactual scenarios remain the same as in the real case. The difference between the real case and the 2019 migration scenario reflects the impacts of the quarantine-induced freezing of the migration and the change in time spent indoors on PWE. The difference between the ‘2019 migration’ and ‘no migration’ scenarios reflects the impact of the Spring Festival migration on PWE in normal years. Our uncertainty analysis shows that despite the overlapping of the PWE uncertainty ranges between the real and counterfactual cases, the effects of migration and time-activity patterns on PWE are significant (Supplementary Methods 3 and 4).

Full size image
Fig. 2: Decomposition analysis of the PWE change between P1 and P2.
figure 2

The overall change in the PWE of the Chinese population is decomposed into the changes in PWE due to the changes in ambient air quality, population relocation, household energy consumption and time spent indoors. Note that the migration had two phases with opposite directions: the first one (during P1) preceded the Spring Festival when people returned to their hometowns; the second (during P2) followed the first as people travelled back to work. The quarantine froze the second phase of the migration, leading to a net difference in the migration impact on PWEs between P1 and P2, as marked in the figure. The impact of the quarantine-induced freezing of the migration in response to COVID-19 in P2 is evaluated by comparing with the PWE under 2019 migration pattern and is also marked in the figure. PWEs are in μg m−3.

Full size image

The effects of migration on PWE

The dynamic cross-province migration dataset we established is based on the national census data26 and official reports14 and is temporally allocated using the Baidu real-time mobility data27 (Methods). The direction of the migration is characterized on a province-to-province basis and further divided into four categories: (1) urban-to-rural, (2) urban-to-urban, (3) rural-to-rural, and (4) rural-to-urban. The migration started about 25 days before the Spring Festival and had two phases with opposite directions—one occurred in P1, the other in P2. Before the Spring Festival (P1), there were estimated 236 million people returning to their hometowns, accounting for one-sixth of the total population. Urban-to-rural migration contributed 53% of the total, of which most were reportedly rural migrant workers14. Urban-to-urban migration contributed 34% and other two types of migration were relatively minor (10% for rural-to-rural and 3% for rural-to-urban). After the Spring Festival (P2) when people would normally move back to cities, however, the nation was under quarantine in response to the outbreak of COVID-19 and the migration froze. The effect of the quarantine on the migration in P2 is clearly illustrated by the day-by-day comparison in the migration intensity in 2020 with the previous year (Fig. 3). The migration in P2 was close to completion within 25 days after the 2019 Spring Festival, by which time this year the migration was only 18% complete (Fig. 3).

Fig. 3: The population migration around the Spring Festivals of 2019 and 2020.
figure 3

The shaded areas illustrate the temporal trend of the number of people migrating each day. The solid lines show the temporal trends of the fraction of rural population (the population residing in rural areas) in the total population. The black dashed line marks the Spring Festival. The x axis represents the calendar date in 2020.

Full size image

The migration led to a shift in the fraction of population residing in rural areas. The fraction reached its maximum of 47.6% during the Spring Festival as in normal years but decreased at a pace one-seventh that of normal years afterwards due to the quarantine (0.05% per day in 2020 versus 0.34% per day in 2019; Fig. 3).

Two main consequences of such a change in the migration for population exposure were: (1) a larger fraction of people exposed for a longer time to HAP in rural households which is usually more severe than in urban households13,18; and (2) increased rural energy consumption to meet the demand of the immigrants, both of which further worsened HAP. On the basis of a recently conducted national survey on rural household energy consumption19 and the indoor exposure model10,12 (Methods), we estimate that by increasing the fraction of rural population, the migration enhanced the nationwide PWE by 3.1 (2.3–4.0) μg m−3 and 7.7 (5.9–10.0) μg m−3 in P1 and P2, respectively, compared to a baseline scenario assuming no migration (Fig. 2), while by increasing the household energy consumption, the migration further increased the PWE by 3.6 (2.7–4.6) μg m−3 and 9.7 (7.4–12.5) μg m−3, respectively, in P1 and P2 (Fig. 2). This amounts to total increases of 6.6 (5.0–8.6) μg m−3 and 17.4 (13.2–22.6) μg m−3 in PWE in P1 and P2, respectively. To isolate the impact of the COVID-19-induced freezing of the migration on PWE, we substitute migration in 2019 for that experienced in 2020 while keeping all other factors equal (outdoor air quality, time spent indoors, baseline energy mix, and so on). The results show a comparable increase in PWE in P1 (6.5 μg m−3 in 2019 versus 6.6 μg m−3 in 2020) but a much smaller increase in P2 (7.2 versus 17.4 μg m−3), suggesting an enhancement of 10.2 (7.7–13.3) μg m−3 (17.4 minus 7.2 μg m−3) in PWE due to the freezing of the migration under the national quarantine.

The contribution of HAP on PWE and the inequality of the PWE change

Focusing on the quarantine period (P2), we consider the changes in HAP and other sectors (transportation, industry and power generation) and assess the overall impacts of the quarantine on indoor and ambient air quality and on PWE. Our assessment shows an estimated decrease of 15.6 (10.0–20.6) μg m−3 in the population-weighted average ambient PM2.5 due to the quarantine, which is similar in magnitude to the PM2.5 reduction before and after the COVID-19 outbreak (16.7 μg m−3) (Fig. 4). We note, however, that unlike our fused PM2.5 field which is the best guess of the real-world PM2.5 concentrations, our impact assessment on ambient PM2.5 using chemical transport model is limited by the uncertainty in the estimation of quarantine-induced emission reduction (Methods) and the capability of the model to reproduce the actual PM2.5 change in response to the emission reduction28, both of which warrant further investigation. The indoor PM2.5 concentration is estimated to increase by 3.1 (−2.4–8.8) μg m−3 due to the quarantine (Fig. 4b) which is a result of the competition between the exacerbation of HAP (12.2 μg m−3) and the improvement in ambient PM2.5 that infiltrates indoors (−9.1 μg m−3). Incorporating the changes in indoor and ambient PM2.5 with population migration and human activities, we estimate that the COVID-19 quarantine led to a net increase of 5.9 (0.2–11.6) μg m−3 in PWE (Fig. 4b).

Fig. 4: The impacts of the responses to COVID-19 on PWEs and the use of solid fuels as a driving factor.
figure 4

a, The PWEs in the real case and in the ‘no COVID-19’ scenario. b, The changes in PWEs due to the responses to COVID-19. c, The shares of solid fuel use in household energy mix. Panels a–c show effects in China, in indoor vs outdoor environments, in urban vs rural areas, in heating vs non-heating regions and in provinces with different per-capita income levels. The shares of solid fuel use in household energy mix in indoor vs outdoor environments are the same as in national total and thus are not shown in c. Error bars in a and b represent 95% CIs.

Full size image

We calculate the contribution of HAP on PWE, which includes the direct contributions to indoor and outdoor PM2.5 and the indirect contribution of outdoor HAP to indoor PM2.5 through infiltration (Supplementary Fig. 1). HAP dominated the PWE in P2 regardless of whether there was a quarantine, whereas the COVID-19-induced quarantine increased the HAP contribution to PWE from 74% (no quarantine or no COVID-19) to 82% (in the real case) (Supplementary Fig. 2).

The contribution of HAP to PWE during this period was higher than that before the COVID-19 quarantine (68%), or in a counterfactual scenario where there was no migration (70%), or for the annual average (62%) (Supplementary Fig. 2). The leading cause of HAP is the use of solid fuels (for example, coals and biomass) for cooking and heating, which is much more prevalent in rural areas (67.5% as the share of solid fuels in the household energy mix) than in urban areas (4.7%) (Fig. 4c). Further investigation shows a clear tendency toward a stronger positive effect of the quarantine on PWE as the share of solid fuel use increased (Fig. 4b,c and Supplementary Fig. 3) and that the PWE in rural areas was estimated to increase by 19.2 (13.5–25.6) μg m−3 due to the quarantine, while the urban PWE decreased by 14.0 (8.6–20.7) μg m−3.

The change in exposure associated with solid fuel use and the contrary changes in rural and urban PWEs are due primarily to the interaction between HAP and the human activities: the longer time spent indoors during the pandemic increased the time length for people being exposed to HAP and thus increased the PWE among rural residents; the freezing of the migration increased the rural household energy consumption and subsequently increased the severity of HAP. On the other hand, in urban areas where indoor air quality is often better than outdoor29, the increase in the time spent indoors reduced PWE.

The association between PWE and HAP led to the spatial heterogeneity (Fig. 5) and population inequality in the quarantine-induced changes in PWE (Fig. 4c) which ranged from −19.0 μg m−3 in Tianjin to 32.5 μg m−3 in Inner Mongolia and from −9.5 μg m−3 in the provinces with average per-capita incomes higher than US$5,000 to 6.5 μg m−3 in the provinces with average per-capita incomes lower than US$3,000, suggesting unequal changes in PWE by income group. The inequality in the PWE changes is further confirmed by the significant negative correlation between the PWE changes and provincial per-capita income levels (P < 0.001) and survives the assessments using county-level data or focusing on the rural population exclusively (Supplementary Fig. 4a,b).

Fig. 5: The spatial distribution of the changes in PWE among the Chinese population due to responses to COVID-19.
figure 5

Changes in ambient and indoor air quality, population migration and time spent indoors are considered. The PWE changes are illustrated by county. The white line marks China’s Qinling Mountains–Huai River Line (Qin–Huai Line) Qin–Huai Line divides China into two regions that differ in climate and is commonly used as a reference line in policy-making to determine the heating (northern) and non-heating (southern) regions56. The map of China map is from Harvard Dataverse and is publicly available under the Creative Commons CC0 Public Domain Dedication57.

Full size image

The urban population does not show significant inequality (Supplementary Fig. 4c) probably due to the much lower dependence on solid fuels and therefore being less affected by HAP than their rural counterparts. Regression analysis reveals a significant interaction (the regression coefficient of the interaction term is −0.69 (−0.84 to −0.54), P < 0.001) between the per-capita income and the epidemic severity in the model to predict the quarantine-induced changes in PWE (Methods and Supplementary Table 2) and suggests that regions with more severe epidemic situation are associated with greater inequality. In Hubei, every 20% reduction in income is estimated to be associated with an increase of 6.7 (5.9–7.5) μg m−3 in PWE (P < 0.001) due to the quarantine, which is almost twice the increase for the national average, 3.4 (3.2–3.5) μg m−3 (P < 0.001) (Supplementary Fig. 4a and Supplementary Table 2).

The effect of Clean Heating Plan on PWE changes

Despite the heterogeneity and inequality, the quarantine-induced increases in PWE in the heating (north) and non-heating (south) regions, including both rural and urban areas, were comparable—6.2 (1.6–10.6) μg m−3 and 5.9 (2.1–9.5) μg m−3 in the heating and non-heating regions, respectively (Figs. 4b and 5). We find that a recently implemented campaign called Clean Winter Heating Plan in Northern China (Clean Heating Plan for short), played an important role in balancing the PWE increases between heating and non-heating regions. Clean Heating Plan was launched by the Chinese central government in 2017 and set stringent and differentiated goals through 2021 toward a high rate of clean heating (the rate of clean energy used for heating) in the northern region, with the rates ranging from 40% in rural areas to 100% in some major cities30. This campaign, if successfully implemented, would reduce the amount of annual coal consumption by 150 Tg, and recent progress has shown much success in the implementation of this campaign such that it is expected to be achieved ahead of schedule31.

We estimate that Clean Heating Plan had phased out 44.4% of the solid fuels used in households in the northern provinces by the end of 2019. If there was no such campaign, we estimate that the increase in PWE induced by COVID-19 would be almost doubled in the heating region (12.0 μg m−3). In addition, the population inequality in the PWE increase, measured by the increase in PWE per 20% reduction in income (4.6 (4.3–4.9) μg m−3), would be 30.1% higher than is estimated in the real case (3.5 (3.2–3.8) μg m−3) in the heating region (Supplementary Fig. 5 and Supplementary Table 2). In an ideal case where Clean Heating Plan was fully phased in, the quarantine would only lead to an increase of 2.3 μg m−3 in PWE in the heating region, with the inequality decreased by 15.6% (Supplementary Fig. 5 and Supplementary Table 2). Our analysis thus reveals that Clean Heating Plan moderated the quarantine-induced increases in PWE in the heating region, reduced the inequality of the PWE increases among different income groups of people, and put the PWE increases of the heating and non-heating regions in the balance. Still, the PWE in the heating region (137 μg m−3) was 61% higher than it was in the non-heating region (85 μg m−3) and, in rural areas, the quarantine-induced increase in PWE in the heating region (24.4 μg m−3) was 31% higher than in the non-heating region (18.6 μg m−3).

Discussion

In this study, we integrate multiple data sources and modelling techniques to dynamically track the changes in PWE due to the national quarantine in China. We first show that the national population-weighted exposure to ambient PM2.5 reduced by 16.7 μg m−3 but the average PWE of the population considering both indoor and outdoor PM2.5 exposure is estimated to increase, which is mainly due to the worsened HAP and a higher opportunity for people to be exposed to indoor HAP during the pandemic. Changes to the actual dose of PM2.5 to the population will also depend on changes in use and effectiveness of facemasks during the period32,33, which deserves further studies.

With respect to the distribution of PWE, our assessment reveals an increase in the environmental inequality of air pollution exposure in response to the COVID-19 crisis. While the high-income group benefited from the reduction of PWE, the low-income group suffered a significant increase in PWE. Such inequality would be even higher if Clean Heating Plan that targets HAP in northern China was not implemented. In addition, given the reported association between short-term exposure to air pollution and the transmission of COVID-19 (ref. 23), this analysis shows how the COVID-19 pandemic itself, as well as the quarantine, may have deepened health inequalities. Our assessment highlights the importance of mitigating HAP for reducing the environmental inequality and protecting human health. If society is to confine people to their homes for their protection, it is far better that they have clean heating to start with.

Our study has several limitations. The emission changes in response to the COVID-19 quarantine are modelled implicitly by assuming constant reduction rates for individual sectors across the nation. Although the resulting changes in CMAQ-modelled outdoor PM2.5 concentrations show general agreement with observation-based PM2.5 changes (data fusion), the model simulation based on such implicit emission estimates is probably subject to large uncertainty. To better reproduce the spatial heterogeneity in the air quality response, subnational-level information on energy consumption, traffic intensities and industrial activities needs to be incorporated in the emission estimation. Note that compared to the CMAQ simulation, the PM2.5 fields derived from the data fusion approach show relatively small uncertainty (Supplementary Methods 1 and 2 and Supplementary Figs. 6 and 7).

Notably, our uncertainty analysis identifies the HAP-induced indoor PM2.5 assessment as the predominant source of uncertainty in PWE (Supplementary Methods 3 and Supplementary Figs. 8 and 9). This large uncertainty associated with HAP is introduced by uncertainties in rural household energy consumption and indoor PM2.5 concentrations, due primarily to the limited availability of such data in China. We show that HAP plays a central role in PWE during both normal periods (contributing >60% of the total PWE) and the quarantine period (80%). Nationwide surveys on household energy consumption and large-scale indoor air quality monitoring networks in rural China (Supplementary Figs. 8 and 10) may provide necessary information to narrow down the uncertainty in PWE. The importance of HAP in overall exposure, the differentiated responses to the quarantine between different population groups, and the large uncertainty in indoor PM2.5 assessment indicate that a better understanding of HAP is probably instrumental in improving air quality while maintaining environmental equity of the Chinese population.

Methods

Household energy consumption

Provincial-level household energy consumption data for 1992 to 2012 were collected and compiled on the basis of a representative national survey19 and China Statistical Yearbook34. The data were downscaled to county level and extrapolated to 2020 (the study year) on the basis of the fuel type-specific empirical models developed by Shen et al.35,36. Following a previous study10, the clean heating targets set by Clean Heating Plan were incorporated into the energy trends in the heating region.

Migration data

We derived detailed origin and destination information from the sixth National Census26 to characterize population migration on the county level36. The census data classified the migrants into four groups (rural-to-urban, urban-to-urban, rural-to-rural and urban-to-rural) and are representative of the migration pattern in 2010. The census data showed a total of 138 million migrant workers in 2010, noting that not all the migrants intended to return home during the Spring Festival holidays. The Ministry of Human Resources and Social Security reported 125 million migrant workers returning home in 202014. Therefore, the census data was scaled down by a factor of 0.9 to represent the migration pattern in 2020. The Spring Festival migration can be divided into two phases with opposite directions. In the first (back-home) phase, migrants go back home to reunite with their families and celebrate the upcoming New Year; in the second (returning-to-work) phase, migrants return to work. The first phase typically ends by New Year’s Eve (that is, the end of the old year and the start of the Spring Festival holidays), followed by the second phase. However, there is a small fraction of people who get back home after New Year’s Eve, which yields a temporal overlap between the first and second phases. We therefore assumed that the back-home migration was achieved before the second day of the Spring Festival holidays, and that the returning-to-work migration started from the first day of the Spring Festival holidays. The migration flows (the number of migrants) were temporally allocated using the daily cross-province mobility intensities reported by the Baidu real-time mobility monitoring platform as a surrogate27. Our assessment showed that the migrants who got back home after New Year’s Eve only occupied 1.6% of the total migrants in the first phase. For the 2019 Spring Festival, of which the detailed provincial-level Baidu mobility data were not available, the national-level mobility intensities were used to scale the 2020 migration pattern to 2019, assuming that the relative difference in the migration flows across provinces remained unchanged between 2019 and 2020.

Survey on human activity pattern

The information on the daily time spent indoors and in different indoor compartments (kitchen, living room and bedroom) in wintertime were derived from Exposure Factors Handbook of Chinese Population37, as summarized by Chen et al.12, and used in this study to represent the time-activity pattern when there was no COVID-19. The time-activity patterns during the pandemic were derived from an online questionnaire survey (https://www.wjx.cn/m/59666734.aspx) which collected information on the frequencies of going out during the quarantine together with demographic information, including gender, age, urban/rural residents, highest level of education, province/city and occupation38. To answer the question on frequencies of going out during the quarantine period, participants were asked to choose one of the following options: (1) more than twice a day, (2) once a day, (3) once per 2 days, (4) once per 3 days, (5) once per 4 days, (6) once per 5 days, (7) once per 6 days, (8) once a week, (9) less frequent than once a week or (10) never.

This survey adopted strict quality control measures during data processing and analysis38. The questionnaires with missing values, logical errors and data format errors were excluded. Two groups of personnel independently derived the data and completed the comparison to ensure the accuracy of the results. A total of 8,330 questionnaires were distributed with a recovery rate of 100%. A total of 7,784 valid questionnaires were obtained, covering 31 provinces in China38.

Among the participants, 3,364 were males and 4,420 were females. There were 183 participants at an age of less than 18 yr, 4,646 between 18 and 30 yr, 1,920 between 45 and 60 yr and 103 older than 60 yr. Information on going-out frequencies by gender and age was provided on Open Science Framework at https://osf.io/x46tb/. Overall going-out frequencies were adjusted by the gender and age distributions of the general Chinese population. The survey showed that the more severe the epidemic, the less frequently people went out each day. Going-out frequencies were translated into time lengths of outdoor stay by assuming the time length for each going-out event ranging from 200 min in the provinces that were the least affected by the COVID-19 outbreak (Qinghai and Tibet) to 120 min in Hubei where the outbreak was the most severe. The uncertainty induced by this assumption was considered in the uncertainty analysis specified in following section. The average time spent indoors by province before and during the pandemic is summarized in Supplementary Table 1.

Emissions and air quality modelling

We used AiMa emission inventory39,40 as the emission input to conduct the air quality modelling for ambient PM2.5 assessment. The emission inventory has been compiled by integrating a variety of sources (for example, statistical data and emission inventories)40 and has undergone continuous updates. This inventory is currently used by the online operational AiMa system (http://www.aimayubao.com/) that provides air quality forecast for the Chinese government and public. The base year of the latest version of the AiMa inventory is 2017.

The ambient PM2.5 concentrations were obtained by combining hourly ground-level observations reported by the China National Urban Air Quality Real-time Publishing Platform5 with model predictions by the Community Multiscale Air Quality (CMAQ) model v.5.0.2 (ref. 41) using an ensemble deep learning data fusion approach25. Meteorological variables were derived from the AiMa system, which were generated by the Weather Research Forecasting (WRF) model v.3.4.1 (ref. 42) driven by the 0.5-degree global weather forecast products produced by the National Centres for Environmental Prediction Global Forecast System43. The downscaled meteorology together with the AiMa emission inventory was used to drive CMAQ simulation which was conducted to cover the mainland China on a horizontal resolution of 12 km with 13 vertical layers extending up to ~16 km above ground. It should be noted that our CMAQ simulation was conducted on a 12-km horizontal resolution, primarily due to the resolution in the emissions inventory, and that PM2.5 and ozone fields are typically relatively homogeneous spatially44,45. Further, the computational costs of large-scale simulations on a finer resolution would be high. The model output was fused with observations from the China National Environmental Monitoring Centre (http://www.cnemc.cn) to get the final ambient PM2.5 concentration fields across China on a daily resolution over the study period (from 25 December 2019 to 25 March 2020). We compared the fused data with observations from an independent dataset and found a good agreement between the two sets of data in terms of the spatial distribution and temporal trends (Supplementary Figs. 6 and 7). This evaluation confirmed the good performance of the data fusion approach in reproducing ground-level PM2.5 concentrations. Details are provided in Supplementary Methods 1 and Supplementary Figs. 6 and 7. More information on the emission inventory, the model configuration, and the data fusion procedures can be found in a previous study25.

We conducted adjoint analysis to decompose the contributions of various emission sources to outdoor PM2.5 concentrations. The emission sources, as categorized in the AiMa inventory, included power generation, industry, residential (household), transportation, agriculture, solvent usage, fugitive dust and fires. CMAQ-Adjoint v.5.0 (ref. 46) was applied to calculate the adjoint sensitivities. The adjoint analysis provides location- and time-specific gradients (adjoint sensitivities) and can be used in applications such as backward sensitivity analysis, source attribution, optimal pollution control, data assimilation and inverse modelling46. The CMAQ-Adjoint v.5.0 is the most up-to-date version of CMAQ-Adjoint. Discrete adjoint is implemented for gas-phase chemistry, aerosol formation, cloud chemistry and dynamics, and diffusion. Continuous adjoint is implemented for advection. The model performance has been comprehensively evaluated in the previous study46, showing good agreements with the results given by forward sensitivity analysis.

In this study, the cost function of the adjoint analysis was defined as the ambient population-weighted average PM2.5 concentration over the study period across China. The adjoint model thus provided sensitivities of this cost function to per-unit emissions of various species in each model grid cell. Using the source-specific emission information, we evaluated the source contributions of household (residential) energy consumption and other sectors on ambient air pollution by province. Details about the principle equations, development and evaluation of CMAQ-Adjoint can be found in previous studies46,47.

Using the adjoint sensitivities, we further evaluated the changes in the population-weighted concentration in response to the emission reduction during the quarantine. Huang et al.28 estimated sector-specific emission reductions due to the COVID-19 quarantine in China on the basis of dynamic economic and industrial activity levels. The description on the emission reduction estimation was detailed in Huang’s study28. Following Huang’s study, we adopted a reduction rate of 10% in power plant emissions, 30% in industrial emissions and 70% in mobile emissions. The changes in residential emissions due to population migration were evaluated using the procedures as specified in our previous studies35,36.

Indoor exposure model

We used an indoor exposure model developed by Chen et al.12 to quantify the indoor PM2.5 levels. The model was modified to take into account the change in the amount of household energy consumption and outdoor infiltration and to unify the estimation approach for urban and rural household conditions as follows

$$C_{\mathrm{in}} = C_{\mathrm{in,add}} + C_{\mathrm{out,add}}$$
(1)

where Cin is the indoor PM2.5 concentration in μg m−3, Cin,add is the Cin component contributed by indoor sources and Cout,add is the Cin component contributed by outdoor infiltration. The value of Cin,add was calculated by the following equation:

$$C_{{\mathrm{in}},{\mathrm{add}}} = \frac{{{\sum} {E_{{f}} \times C_{{f,k}} \times T_{{k}}} }}{{\overline E \times {\sum} {T_{{k}}} }}$$
(2)

where subscripts f and k denote the type of fuel (wood, straw, coal and cleaner energy) and indoor compartment (kitchen, living room and bedroom), respectively; Ef is the per-household daily consumption of fuel type f in terms of thermal energy amount (the amount of energy consumption after thermal efficiency conversion); Ē is the average per-household daily thermal energy required for cooking and heating; Cf,k is the Cin,add in indoor compartment k when Ef = Ē and the household consumes fuel f solely; Tk is the time spent daily in indoor compartment k. Following a previous study48, the thermal efficiencies of biomass, coal, gas and electricity are 0.154, 0.244, 0.555 and 0.84, respectively. The Ē value 40 MJ day−1 household−1 was calculated as the national average daily household thermal energy consumption for cooking and heating in winter. Values for Cf,k were adopted from a previous study12 in which the means and variations of Cf,k were determined by a meta-analysis of 28 field measurement studies in China (see https://osf.io/x46tb/ for the complete dataset). The mean heating-season Cf,k in kitchen/living room are 283, 434 and 547 μg m−3 for coal, crop and wood, respectively, and in bedroom are 211, 267, 359 μg m−3 for coal, crop and wood, respectively. Cleaner energy was assumed to cause little addition to indoor PM2.5 and thus the Cf,k for cleaner energy was set to be 0. Equation (2) assumes that with all others being equal, Cin,add is proportional to the thermal amount of daily energy consumption of the household. This assumption was testified and supported by sensitivity tests using a single-box model49, as recommended in World Health Organization’s indoor air quality guidelines50, to predict Cin,add on the basis of varying amounts of energy consumption. The value Cout,add was calculated by multiplying ambient PM2.5 concentrations with region-specific infiltration factors following the method of Xiang et al.51. The modelled temporal trend in indoor PM2.5 concentrations captured observed temporal trends measured in 30 households in a rural area of Shandong province during the study period (Supplementary Fig. 10). The PM2.5 exposure of individuals at a specific location was calculated as the average of the indoor and outdoor PM2.5 concentrations weighted by the time fractions of indoor and outdoor stays. The calculated exposure was evaluated against personal exposure measurements reported in the literature (Supplementary Fig. 8), which showed general agreement. The PWE in a region was calculated as the population-weighted average of the individuals’ exposure within this region. The same approach to calculating PWE has been adopted in previous studies10,11.

Regression analysis

We conducted regression analysis to predict the county-level quarantine-induced changes in PWE. The regression showed significant interaction between per-capita income and the epidemic severity. The regression equation including the interaction term is as follows

$$\mathrm{dPWE} = - {\mathrm{31}}{\mathrm{.9}} \times \ln (\mathrm{INC}_{\rm{per}}) - {\mathrm{0}}{\mathrm{.69}} \times \mathrm{SEV} \times \ln (\mathrm{INC}_{\rm{per}}) + {\mathrm{124}}{\mathrm{.6}}$$
(3)

where dPWE denotes the change in PWE due to the COVID-19-induced quarantine in μg m−3; INCper is per-capita annual income in US$; SEV is the epidemic severity determined by the confirmed cases in the provinces (Supplementary Table 1), ranging from 1 in Qinghai and Tibet (the least severe) to 5 in Hubei (the most severe). In addition to equation (3), we developed single models (excluding the interaction term) to investigate the relationship between per-capita income and dPWE for specific regions (for example, Hubei, Northern China and urban/rural areas) and for different counterfactual cases (for example, assuming no Clean Heating Plan or assuming that Clean Heating Plan has been fully phased in). See Supplementary Table 2 for the full model results, including coefficients, standard errors and P values of individual models.

Uncertainty analysis

The uncertainty in the PWE estimates stemmed from various sources, including the uncertainties in the modelled ambient and indoor concentrations, population migration, time-activity patterns and infiltration factors. We conducted Monte Carlo simulations to propagate the uncertainties from these input sources to PWE. For some input variables (for example, concentration, migration intensity and time spent indoors), we assumed log-normal distributions to avoid negative values and used geometric coefficient of variation (GCV)52 instead of coefficient of variation (CV) to measure the uncertainty. GCV is defined as follows:

$${\mathrm{GCV}} = \sigma _{\mathrm{g}} - 1$$
(4)

where σg is the geometric standard deviation53. Variable uncertainties measured by the same value of CV and GCV have similar widths of uncertainty intervals. For example, the uncertainties with CV and GCV of 20% have semi-interquartile ranges (half of the difference between the upper quartile and the lower quartile) equal to 13.6% and 13.5% of the variable’s mean, respectively, and have 95% CI widths equal to 79% and 80% of the variable’s mean, respectively.

The GCVs of fused ambient PM2.5 concentrations were derived on the basis of the performance of the data fusion approach (Supplementary Methods 2). The GCV differed by region and urban/rural area, depending on the region-specific performance of the data fusion approach and the number of monitoring sites available within each region (Supplementary Methods 2). The GCV of the population migration intensity was assumed to be 20%. The GCV of the time spent indoors in normal days was set to be 5%, following Chen et al.12. We assigned a CV (instead of GCV) of 20% to the changes in the time spent indoors due to the quarantine, thus resulting in a larger uncertainty in the time spent indoors during the quarantine period than in normal days.

For Ē, we assumed a uniform distribution with a variation interval of 20% which is usually applied to reflect the uncertainty in the statistics of household solid fuel use35,54. The CV of infiltration factors in indoor/outdoor air exchange was set to be 12.5% following Shi et al.55. The uncertainties in indoor PM2.5 concentrations in households using solid fuels were derived by Chen et al. on the basis of over 2,000 observations collected from the literature12. Given the large uncertainty in the estimated emission reduction due to the responses to COVID-19, CVs of the magnitudes of the emission reduction were set to be 30% for all sectors. GCVs were converted to standard deviations on the logarithmic scale by the follow equation

$$\sigma _{\ln } = \ln \left( {{\mathrm{GCV}} + {\mathrm{1}}} \right)$$
(5)

where σln is the standard deviation of a given parameter on the natural logarithmic scale and was used to determine probability distributions of corresponding parameters in Monte Carlo simulations for uncertainty analysis. Monte Carlo simulations were performed 1,000,000 times to propagate the uncertainties in input variables into the uncertainty in PWE.

In addition, we investigated the contributions of various sources to the final uncertainty in PWE. The uncertainty sources we considered included the data fusion of ambient PM2.5, the HAP-related indoor PM2.5 estimation, the migration intensity, the time-activity pattern and the infiltration factors. The result showed that the indoor PM2.5 concentration is the most important source of uncertainty in the estimated national average PWE. Details can be found in Supplementary Methods 3 and Supplementary Fig. 9.

Reporting Summary

Further information on research design is available in the Nature Research Reporting Summary linked to this article.