Cohort-based long-term ozone exposure-associated mortality risks with adjusted metrics: A systematic review and meta-analysis.

Long-term ozone (O3) exposure may lead to non-communicable diseases and increase mortality risk. However, cohort-based studies are relatively rare, and inconsistent exposure metrics impair the credibility of epidemiological evidence synthetization. To provide more accurate meta-estimations, this study updates existing systematic reviews by including recent studies and summarizing the quantitative associations between O3 exposure and cause-specific mortality risks, based on unified exposure metrics. Cross-metric conversion factors were estimated linearly by decadal observations during 1990-2019. The Hunter-Schmidt random-effects estimator was applied to pool the relative risks. A total of 25 studies involving 226,453,067 participants (14 unique cohorts covering 99,855,611 participants) were included in the systematic review. After linearly unifying the inconsistent O3 exposure metrics , the pooled relative risks associated with every 10 nmol mol-1 (ppbV) incremental O3 exposure, by mean of the warm-season daily maximum 8-h average metric, were as follows: 1.014 with 95% confidence interval (CI) ranging 1.009-1.019 for all-cause mortality; 1.025 (95% CI: 1.010-1.040) for respiratory mortality; 1.056 (95% CI: 1.029-1.084) for COPD mortality; 1.019 (95% CI: 1.004-1.035) for cardiovascular mortality; and 1.074 (95% CI: 1.054-1.093) for congestive heart failure mortality. Insignificant mortality risk associations were found for ischemic heart disease, cerebrovascular diseases, and lung cancer. Adjustment for exposure metrics laid a solid foundation for multi-study meta-analysis, and widening coverage of surface O3 observations is expected to strengthen the cross-metric conversion in the future. Ever-growing numbers of epidemiological studies supported the evidence for considerable cardiopulmonary hazards and all-cause mortality risks from long-term O3 exposure. However, evidence of long-term O3 exposure-associated health effects was still scarce, so more relevant studies are needed to cover more populations with regional diversity.

Introduction

Atmospheric ozone (O3) is a short-lived climate forcer. Besides warming the global atmosphere, O3 in the stratosphere can abate the radiation hazards from UV rays on organisms, while O3 in the ambient air has detrimental effects on the ecosystem and human health,2, 3, 4 so health effects caused from exposure to surface O3 have become a serious public concern. Short-term (i.e., hours to days) exposure to high-level O3 can cause acute symptoms like asthma, respiratory tract infection, myocardial infarction, and cardiac arrest;5, 6, 7, 8 and long-term (i.e., over years) exposure can lead to chronic health conditions including but not limited to preterm delivery, stroke, chronic obstructive pulmonary diseases, and cerebrovascular diseases.9, 10, 11, 12 Long-term ambient O3 exposure was estimated to be responsible for over 0.36 million premature deaths globally in 2019, according to the Global Burden of Disease (GBD) report released by the Institute for Health Metrics and Evaluation (IHME).

Systematic reviews summarizing the associations between adverse health outcomes and both short-term and long-term O3 exposure have been performed in previous studies.14, 15, 16 Studies on short-term O3 exposure-induced morbidities are comparatively more abundant than the long-term O3 exposure studies where the epidemiological evidence is less congruous. Some deficiencies exist in the two reviews of long-term O3 exposure-associated mortality risk studies,, the primary issue being the inconsistent use of various O3 exposure metrics; however, no other reviews are found to remedy these flaws. As a secondary photolytic gaseous air pollutant, the warm-season and diurnal concentrations of surface O3 will be much higher than cool-season and nocturnal concentrations,, so the average and peak metrics of O3 concentrations have drastically different implications. Therefore, directly pooling the relative risks scaled with different metrics may lead to biases.

Atkinson et al. (2016) explored six types of mortality causes, but searched the literature only until 2015. Huangfu et al. (2019) updated the search to 2018, but only three types of mortality causes were considered. We update the review of the health effects of O3 to include more categories of mortality, together with covering the most recent studies. Additionally, GBD estimations ascribed long-term O3 exposure induced all-cause mortality for chronic obstructive pulmonary disease, which might lead to underestimations without considering other causes. It is reasonable to deduce that long-term O3 exposure will exacerbate the severity of certain diseases given that short-term exposure increases the morbidity risks of the same diseases, so scrutinizing epidemiological evidence for multiple causes of mortality will provide more credible support to fill in this gap.

The primary innovation of this updated review derives from taking full advantage of global stationary observations to explore the feasibility of adjusting the various exposure metrics, and pooling the multi-study risks with the unified exposure metric, the mean of warm-season daily maximum 8-h average, in response to the recent suggestions from the Lancet global environmental health collaboration. Through this updated systematic review and meta-analysis on long-term O3 exposure-associated cause-specific mortality risks, we aim to present and evaluate the epidemiological evidence for three major questions not fully addressed by the two previous reviews: (1) which mortality causes can be ascribed to long-term O3 exposure, (2) whether the risk associations have changed given the latest studies, and (3) how to estimate the risk association strengths with the suggested exposure metric. Both our methods and discoveries are expected to inspire future O3 health studies, and support relevant policy-making to benefit the global population.

Methods

Search strategy

We searched three research databases (MEDLINE, Embase, and Web of Science) from September 1, 2015 to February 1, 2022, to include the latest studies in our systematic review and meta-analysis, updating the studies included in two previous reviews on long-term O3 exposure-associated mortality., Search terms were similar to these two previous systematic reviews with modifications to enhance the inclusion of potentially relevant studies, as we combined the keywords related with the cause-specific mortalities (i.e., “mortality,” “death,” “premature death,” “all-cause,” “non-accidental,” “cardiopulmonary,” “respiratory,” “chronic obstructive pulmonary disease,” “pneumonia,” “cardiovascular,” “lung cancer,” “cerebrovascular,” “stroke,” “ischemic heart disease,” “congestive heart failure”), the pollutant of research interest (i.e., “ozone”), and qualified epidemiological study types (i.e., “long-term,” “cohort study,” “prospective,” “retrospective,” “longitudinal study”). The detailed search strategies are listed in Table S1. Health outcomes considered in the systematic review were as follows: mortality from (1) all causes (AC, ICD9: 001–799, ICD10: A00–R99); (2) all respiratory diseases (RESP, ICD9: 460–519, ICD10: J00–J98); (3) chronic obstructive pulmonary diseases and allied conditions (COPD, ICD9: 490–496, ICD10: J19–J46); (4) all cardiovascular diseases (CVD, ICD9: 390–459, ICD10:I00–I98); (5) all cerebrovascular diseases (CEVD, ICD9: 430–438, ICD10: I60–I69); (6) ischemic heart disease (IHD, ICD9: 410–414, ICD10: I20–I25); (7) congestive heart failure (CHF, ICD9: 428, ICD10: I50); and (8) lung cancer (LC, ICD9: 162, ICD10: C33–C34).

Study eligibility criteria

As an updated systematic review, studies identified in the previous two reviews were examined together with the newly retrieved ones. Studies were included during screening following the criteria as follows: (1) the epidemiological research should be conducted based on cohorts; (2) the exposure should include O3 as an individual risk factor; (3) the health outcomes should be all-cause or cause-specific deaths at individual level; (4) studies provided hazard ratio (HR), risk ratio (RR), or odds ratio (OR) and their 95% confidence intervals (CIs) clearly and reported per increased unit (e.g., 10-ppbV) of exposure concentration, assuming linear risk relationships with adjustment of key confounders; (5) the study was published as an original research article in peer-reviewed journals in English. For articles from the same cohort, only one study covering the widest populations and the longest follow-up period was reserved for meta-analysis, unless the subgroups of participants and study follow-up periods were clearly stated to be of mild overlap. We followed the Preferred Reporting Items for Systematic Review and Meta-Analyses (PRISMA) guidelines to process the included studies on ambient O3 exposure-induced mortality.

Study selection and scrutinization

All studies found from the literature search were archived in Clarivate Analytics Endnote X9.3.1 reference manager software. Two literature review investigators (H.Z.S. and C.L.) conducted title and abstract pre-screening independently for all web-searched records and reviewed the full text for the pre-screened studies. Disagreements were resolved by discussions with a third reviewer (P.Y.).

Data extraction

Details from each selected study were extracted and labeled for the purpose of meta-analysis, including (1) the authors with publication year as study labels of reference; (2) basic descriptive information of the study cohort, including the cohort name, country, follow-up periods, numbers of cases and total participants, population genders and ethnicities, exposure metrics, health outcomes, and major confounders; (3) the risk association strengths (exchangeable with effect size or risk values in terminology) preferably quantified in HR (and also RR/OR as substitute choices) per unit incremental exposure with 95% CI.

Study quality assessment

All selected studies underwent quality evaluation using the Quality Assessment Tool of Observational Cohort and Cross-Sectional Studies developed by the National Institute of Health (NIH) (https://www.nhlbi.nih.gov/health-topics/study-quality-assessment-tools), aiming to ensure that the studies considered for meta-analysis were adequately reliable. The assessments were cross-validated by two authors (H.Z.S. and C.L.) independently, with the third author (P.Y.) adjudicating on any disagreements. Table S2 lists 14 assessment items assigned with one score for each, and the tallied scores were translated into a literature-specific rating of quality. Studies scoring full marks, 14, were categorized to be “Good,” with 10–13 as “Fair” and <10 as “Poor.”

Besides applying the quality assessment tool to determine which reviewed studies should be included for meta-analysis, the epidemiological evidence quality of the included studies for each cause of mortality was evaluated with the Grading of Recommendations Assessment, Development, and Evaluation (GRADE) system, to yield rating bands “high,” “moderate,” “low,” and “very low.” This grading system by default rated “high-quality” for cohort studies as the starting point of evaluation, and the rate was downgraded by five limitations: the existence of (1) risk of bias examined by the quality assessment tool (Table S2), (2) imprecision (i.e., studies did not report the central risk estimations with CIs), (3) inconsistency (i.e., the directions of the estimated risks were controversial across studies), (4) indirectness (i.e., studies did not include the desired population, exposure, or health outcomes), and (5) publication bias (i.e., researchers tended to publish studies with positive results). Studies could be upgraded by three strengths: reporting (1) exposure-response trend, (2) residual confounding (i.e., adjusting the confounders highlighted the risks), and (3) strong associations. Publication biases were graphically presented by funnel plots, and statistically tested by the trim-and-fill method. The review was registered in PROSPERO (CRD42021270637).

Exposure adjustment

Unit unification

There were two major units used to quantify the surface O3 concentrations, nmol mol-1 (or parts per billion by volume mixing ratio, ppbV), more frequently used by atmospheric modeling researchers,,, and milligram per cubic meter by mass concentration (μg/m3) widely used by public health studies. These two units are interchangeable based on the ideal gas law PV = nRT, if the air temperatures (T) and pressures (P) are given, as presented in Equations 1, 2, 3, and 4.

Assuming T = 298.65 K (25.5°C) and P = 101.325 kPa, the ppbV to μg/m3 conversion factor could be approximated as 1 ppbV ∼1.96 μg/m3. Though the surface air temperatures and pressures would vary across seasons, this simplification was widely used in previous studies,,, being more credible for long-term surface O3 studies that average the surface air temperatures and pressures over longer periods. For example, even at the very low temperature of 270 K, the conversion factor was 2.17, which corroborated the stability of linear conversion.

Metric unification

Surface O3, as a secondary photochemical pollutant produced by photolysis of NO2 to trigger chains of radical reactions, has concentrations that vary significantly between daytime and nighttime, and between warm and cool seasons, as discussed by a number of studies.,28, 29, 30, 31 Various daily metrics to quantify the surface O3 concentrations have emerged due to a series of considerations, which has led to greater difficulty when assimilating various epidemiological evidence. The previous reviews simply pooled the reported risk association strengths without adjusting the diverse metrics., We have improved upon this approach.

We updated the meta-analysis by unifying the exposure metrics for pooled O3 exposure-associated risks. As suggested by the US EPA final report of Air Quality Criteria for Ozone and Other Photochemical Oxidants, linear relationships were assumed to estimate the cross-metric conversion factors using long-term reliable observations such as the Tropospheric Ozone Assessment Report archive (TOAR, https://b2share.fz-juelich.de/communities/TOAR) and China National Environmental Monitoring Center (CNEMC, http://www.cnemc.cn/en/) in our review, and correlation matrix was used to validate that the assumptions of linearity were not violated. Both TOAR and CNEMC sites measured the surface O3 by means of the UV absorption technique with strict quality control to ensure the comparability of the records across different countries and regions., We considered six complex metrics for mutual conversion as (1) annual mean of 24-h daily average (ADA24), (2) 6-month warm-season mean of 24-h daily average (6mDA24), (3) annual mean of daily maximum 8-h average (ADMA8), (4) 6-month warm-season mean of daily maximum 8-h average (6mDMA8), (5) annual mean of daily maximum 1-h average (ADMA1), and (6) 6-month warm-season mean of daily maximum 1-h average (6mDMA1). Long-term averaging-based metric conversion could smooth the temporal variations resulting from the seasonal and geographic solar radiation variabilities. The linear conversion factors (k) were mathematically defined by Equation 5, to adjust the original metric to the target one with irreducible regression errors є. In principle, such a metric conversion method only adjusts the regression-obtained RRs without altering the exposure-mortality associations.

Meta-analysis

We collectively named RR for HR/RR/OR throughout our meta-analysis. All literature-reported RRs were converted into adjusted incremental risk ratios with a 10-ppbV O3-exposure increase by target metric (i.e., 6mDMA8 in this study), following Equation 6 as shown below:where ln is the natural logarithm, RR is the originally reported risk estimates scaled into 10-ppbV incremental exposure, and is the conversion factor for metric unification. Multi-study pooled risks with 95% CI were calculated from the adjusted RRs by Hunter-Schmidt random-effects meta-regression estimator. The random-effects model used the inverse variances of the risk association CIs to weight each study, assuming the reported RRs conformed to an underlying distribution (Gaussian distribution by default). The Hunter-Schmidt approach could provide an effective correction to the potential systematic errors and biases caused from the diversity of study population and methodologies, by subtracting the relevant part of variances.

We applied the Higgins I2 to quantify the heterogeneity across studies. The Higgins statistics I2 is defined aswhere Q is the Cochran’s non-parametric heterogeneity statistic assessing whether there are any cross-study differences in risks based on χ2 distribution, and df is the corresponding degrees of freedom. Low I2 values indicate no important heterogeneity observed, and high I2 values, especially >75%, indicate considerable heterogeneity.

Subgroup analyses were conducted by grouping the selected studies upon the gender, regions, O3 exposure metrics, and methodological reliability of individual exposure assignment, together with the adjustment of ethnicity, body mass index (BMI), smoking history, lifestyle features, and exposure to PM2.5 and NO2. Subgroups had to contain at least three studies. Leave-one-out sensitivity analyses were also carried out to test the robustness of synthesized overall risks by meta-analysis. All meta-analyses were performed in R 4.1.1 with packages meta, metafor, and metainf.

The most widely recognized approach to construct the integrated exposure-response (IER) relationships required sufficient epidemiological studies to comprehensively sample the population exposure levels. However, studies on long-term O3 exposure health effects were relatively limited, under which circumstance we made methodological modifications to make better use of the variabilities in exposure levels by statistically imputing the exposure distributions for each study from the provided statistics (e.g., mean, standard deviation, and percentiles) for curve fitting as described in supplemental information S1. Supplemental information S2 describes the detailed procedures of exposure distribution imputations, with a demonstration provided in supplemental information S3, through which high uncertainties were still observed in the fitted IER curves due to insufficient epidemiological studies.

Results

Study characters

From the three databases searched between September 1, 2015 and February 1, 2022, a total of 339 studies (77 from MEDLINE, 102 from Embase, and 160 from Web of Science) were found; and together with 34 additional studies added manually from the two previous systematic reviews,, 373 studies underwent duplication censoring, deleting 101 duplicated studies. After detailed scrutinization of the remaining 272 studies, a total of 25 studies concerning long-term O3 exposure and multi-cause mortalities were included for quality evaluation, meta-analysis, and further discussions (Figure 1).38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62 Table 1 summarizes the basic information of the 25 included studies sorted by the publication year and surname of the first author, and Table S3 lists the key confounders adjusted in each study.

Figure 1

Schematic flowchart of study assessment and selection processes for literature review and meta-analysis

Table 1

Summary of cohort characteristics included for meta-analysis

Study	Cohort	Country	Follow-up duration	Population type	Sample size	Sex	Age	Mortality causes
Abbey et al. 199938	AHS	USA	1977–1992	occupational	6,182	FM	27–95	AC, RESP, LC
Lipfert et al. 200639	WU-EPRI	USA	1976–1996	general	67,108	M	51 (12)a	AC
Jerrett et al. 200940	ACS CPS II	USA	1977–2000	general	448,850	FM	≥30	AC, RESP, CVD, IHD
Krewski et al. 200941	ACS CPS II	USA	1982–2000	general	488,370	FM		AC, IHD, LC
Smith et al. 200942	ACS CPS II	USA	1982–2000	general	352,242	FM		AC, RESP, CVD
Lipsett et al. 201143	CTS	USA	1998–2005	occupational	124,614	F	≥20	AC, RESP, CVD, IHD, CEVD, LC
Zanobetti et al. 201144	Medicare	USA	1985–2006	general	8,894,473	FM	≥65	COPD, CHF
Carey et al. 201345	CPRD	UK	2003–2007	general	824,654	FM	40–89	AC, RESP, LC
Jerrett et al. 201346	ACS CPS II	USA	1982–2000	general	73,711	FM	57 (11)	AC, RESP, CVD, IHD, LC
Bentayeb et al. 201547	GAZEL	France	1989–2013	occupational	20,327	FM	44 (4)	AC, RESP, CVD
Crouse et al. 201548	CANCHEC	Canada	1991–2006	general	2,521,525	FM	≥25	AC, RESP, COPD, CVD, IHD, CEVD, LC
Tonne et al. 201649	MINAP	UK	2003–2010	MI survivorsb	18,138	FM	68 (14)	AC
Turner et al. 201650	ACS CPS II	USA	1982–2004	general	669,046	FM	≥30	AC, RESP, COPD, CVD, CHF, IHD, CEVD
Di et al. 201751	Medicare	USA	2000–2012	general	60,925,443	FM	≥65	AC
Weichenthal et al. 201752	CANCHEC	Canada	2001–2011	general	2,448,500	FM	25–89	AC, RESP, CVD
Cakmak et al. 201853	CANCHEC	Canada	1991–2011	general	2,291,250	FM	≥25	AC, COPD, IHD, LC
Hvidtfeldt et al. 201954	DDCH	Denmark	1993–1997	general	49,596	FM	50–64	AC, RESP, CVD
Kazemiparkouhi et al. 201955	Medicare	USA	2000–2008	general	22,159,190	FM	≥65	AC, RESP, COPD, CVD, IHD, CHF, CEVD, LC
Lim et al. 201956	NIH-AARP	USA	1995–2011	general	548,780	FM	50–71	AC, RESP, COPD, CVD, IHD, CHF, CEVD, LC
Paul et al. 202057	ONPHEC	Canada	1996–2015	diabetes	452,590	FM	35–85	CVD
Shi et al. 202158	Medicare	USA	2001–2017	general	44,684,756	FM	≥65	AC
Strak et al. 202159	ELAPSE	six countriesc	1985–2015	general	325,367	FM	49 (13)	AC, RESP, COPD, CVD, IHD, CEVD
Yazdi et al. 202160	Medicare	USA	2000–2016	general	44,430,747	FM	≥65	AC
Bauwelinck et al. 202261	BC2001	Belgium	2001–2011	general	5,474,470	FM	≥30	AC, RESP, CVD, LC
Stafoggia et al. 202262	ELAPSE	seven countriesd	2000–2017	general	28,153,138	FM	≥30	AC, RESP, CVD, LC

Cohort abbreviations: AHSMOG, Adventist Health Study of Smog; WU-EPRI, Washington University–Electric Power Research Institute; ACS CPS, American Cancer Society Cancer Prevention Study; CTS, California Teacher Study; CPRD, Clinical Practice Research Datalink; GAZEL, GAZ de France and ÉLectricité; CANCHEC, Canadian Census Health and Environment Cohort; MINAP, MyocardialIschaemia National AuditProject; DDCH, Danish Diet, Cancer and Health; NIH-AARP, National Institute of Health, American Association of Retired Persons; ONPHEC, Ontario Population Health and Environment Cohort; BC2001, Belgian 2001 Census.

Key confounding adjustments were listed in Table S3.

Population ages were reported by mean with standard deviation (in parenthesis).

MI, myocardial infarction.

Sweden, Denmark, France, Netherlands, Germany, and Austria.

Belgium, Denmark, England, Netherlands, Norway, Switzerland, and Italy.

Metrics and exposure assignments

Our updated systematic review focuses on the exposure metrics and methodologies used to obtain O3 exposure, as summarized in Table 2. Abbey et al. (1999), Jerrett et al. (2013), and Lipsett et al. (2011) did not clearly state the metric they used, but based on comparisons between the reported surface O3 concentrations and TOAR observational archives, we reasonably assumed ADA24 for the first study, and ADMA8 for the other two. Details of the metric matching are given in supplemental information S4. Lipfert et al. (2006) used the highest 95th percentile by hourly resolved O3 concentrations as the peak exposure metric, which was only used in this one study, so it was approximated to DMA1. Krewski et al. (2009) and Smith et al. (2009) were both studies on ACS CPS II, so the same exposure assignment methodologies and metrics were assumed as Jerrett et al. (2009). Likewise, Cakmak et al. (2018) and Weichenthal et al. (2017) were assumed to follow the methodology of Crouse et al. (2015) as all three studies were on CANCHEC. Warm season was defined as the 6 months from April to September for the northern hemisphere by default, except for three studies due to the limited number of studies included: Zanobetti et al. (2011) using May to September, and Crouse et al. (2015) and Paul et al. (2020) using May to October.

Table 2

Data sources and statistical methods of O3 exposure assignment

Study	Data sources	Methods	Resolution	Rating	Metrics	Level of incremental risk ratio
Abbey et al. 199938	monitoring station observations	IDW interpolation	N/Ra	low	ADMA8	12.03 ppbV
Lipfert et al. 200639	monitoring station observations	nearest matching (assumed)b	N/R	low	ADMA1	40 ppbV
Jerrett et al. 200940	monitoring station observations	nearest matching (assumed)	N/R	low	6mDMA1	10 ppbV
Krewski et al. 200941	monitoring station observations	ordinary kriging interpolation	N/R	low	6mDMA1	10 ppbV
Smith et al. 200942	monitoring station observations	nearest matching (assumed)	N/R	low	6mDMA1	1 μg/m³
Lipsett et al. 201143	monitoring station observations	IDW interpolation	250 m	low	ADA24	22.96 ppbV
Zanobetti et al. 201144	monitoring station observations	nearest matching (assumed)	N/R	low	6mDMA8	5 ppbV
Carey et al. 201345	monitoring station observations	interpolation (IDW assumed)	1 km	low	ADA24	3.0 μg/m3
Jerrett et al. 201346	monitoring station observations	IDW interpolation	N/R	low	ADA24	24.1782 ppbV
Bentayeb et al. 201547	monitoring station observations, model simulation, other auxiliary predictors	universal kriging-embedded land use regression	2 km	good	6mDMA8	12.3 μg/m3
Crouse et al. 201548	monitoring station observations, model simulation	linear data assimilation	21 km	good	6mDMA8	9.5 ppbV
Tonne et al. 201649	KCLurban air dispersion model simulation	N/Ac	20 m	moderate	ADA24	5.3 μg/m3
Turner et al. 201650	monitoring station observations, CMAQ model simulation	hierarchical Bayesian space-time data assimilation	12 km	high	ADMA8 6mDMA8	10 ppbV
Di et al. 201751	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	ensemble machine learning	1 km	high	6mDMA8	10 ppbV
Weichenthal et al. 201752	monitoring station observations, model simulation	linear data assimilation	21 km	good	6mDMA8	10.503 ppbV
Cakmak et al. 201853	monitoring station observations, model simulation	linear data assimilation	21 km	good	6mDMA8	10 ppbV
Hvidtfeldt et al. 201954	AirGIS dispersion model simulation	N/A	1 km	moderate	ADA24	10 μg/m3
Kazemiparkouhi et al. 201955	monitoring station observations	nearest matching (assumed)	6 km	low	6mDMA1 6mDMA8 6mDA24	10 ppbV
Lim et al. 201956	monitoring station observations, CMAQ model simulation	Bayesian space-time downscaling	12 km	high	6mDMA8	10 ppbV
Paul et al. 202057	monitoring station observations, model simulation	linear data assimilation	21 km	good	6mDMA8	6.4 ppbV
Shi et al. 202158	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	ensemble machine learning	1 km	high	6mDMA8	10 ppbV
Strak et al. 202159	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	universal kriging-embedded land use regression	100 m	high	6mDMA8	10 μg/m3
Yazdi et al. 202160	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	ensemble machine learning	1 km	high	6mDMA8	1 ppbV
Bauwelinck et al. 202261	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	land use regression	100 m	high	6mDMA8	10 μg/m3
Stafoggia et al. 202262	monitoring station observations, model simulation, satellite remote sensing observations, other auxiliary predictors	universal kriging-embedded land use regression	100 m	high	6mDMA8	10 μg/m3

Methodological ratings were based on spatial interpolation and multi-data assimilation approaches. Spatial resolutions, exposure metrics, and levels of incremental risk ratio were also listed.

N/R, not reported.

The statistical methods were not clearly stated in literature, so the most basic method was assumed. The nearest neighborhood matching shall be the simplest way to assign spatially sparse observations onto cohort participants, and the inverse distance weighting (IDW) is the simplest spatial interpolation approach.

N/A, not applicable. The chemical transport model simulations were directly used for individual exposure assignment without further statistical processing.

Across all included studies, multiple methods were applied to obtain gap-free surface O3 concentrations for individual-level exposure assignment. The most basic approach was to match nearest neighbors between participant residential locations and in situ observation sites, which were more frequently used in earlier studies., A comparatively more complicated approach was statistical spatial interpolation, by inverse distance weighting or ordinary kriging. Full spatial coverage products, such as satellite-based remote sensing and chemistry transport models, were used in some studies by supervised-learning-based data fusion techniques including but not limited to universal kriging-embedded land use regression, Bayesian hierarchical model, and ensemble learning to enhance spatial extrapolation accuracy, which were considered to be of higher credibility than the aforementioned basic ones. All basic interpolation methods using merely the observations were rated as “Low,” applying chemical transport model simulations without calibration from the observations as “Moderate,” linearly coupling the observations with simulations as “Good,” and multi-source data assimilation by means of more sophisticated approaches as “High.” In total, eight studies were rated “High,” five were “Good,” two were “Moderate,” and 10 were “Low.” Methodological progress with time was evident as shown in Table 2, indicating an explosion of population-based environmental health studies in the age of big data.

Based on the TOAR and CNEMC in situ observations, the cross-metric non-intercept linear conversion factors were estimated with regression accuracies given in Figure 2. As suggested by relevant recent studies, the 6mDMA8 metric was typically recommended to highlight peak exposure. Therefore, we chose to convert the originally reported RRs uniformly into the 6mDMA8 scale as standard. O3 exposure levels by the original and unified metrics are listed in supplemental information S1. Demonstrations for the conversion interpretation and procedures are presented in supplemental information S5.

Figure 2

Cross-metric linear relationships and conversion accuracies

The cross-metric linear relationships were quantified by Pearson’s correlation coefficients. The cross-metric conversion factors with 95% confidence intervals (95% CI) were estimated by non-intercept linear regression models, accompanied with fitting accuracies quantified by coefficient of determination (R2) and root-mean-square error (RMSE) in ppbV. The conversion factors were defined as multiples from the original metric by column into the target harmonized metric by row, e.g., ADMA8 = 1.671×ADA24, R2 = 0.9736, RMSE = 7.78 ppbV. Note that by non-intercept linear regression, the values of R2 should no longer be equal to the squared Pearson’s linear correlation coefficients. As the cross-metric conversion coefficients were estimated statistically, indirect conversions were not recommended, since regression noises restricted the validity of equation .

Meta-analysis results

We conducted meta-analyses for long-term O3 exposure-associated mortalities across 8 categories: (1) all causes (AC), (2) all respiratory diseases (RESP), (3) chronic obstructive pulmonary disease and allied conditions (COPD), (4) all cardiovascular diseases (CVD), (5) all cerebrovascular diseases (CEVD), (6) ischemic heart disease (IHD), (7) congestive heart failure (CHF), and (8) lung cancer (LC), with the exposure metrics adjusted to 6mDMA8. All selected studies measured effect sizes based on time-to-event survival outcomes. Therefore, the reported RRs (generalized definition) were equivalent to the HRs (narrow definition), thus negating the need to consider any RR/HR/OR adjustment throughout this meta-analysis study.

All-cause mortality

A total of 23 studies were included into the O3 exposure-associated all-cause mortality meta-analysis, pooling the overall risk into RR = 1.014 (95% CI: 1.009–1.019, I2: 97.8%) per 10-ppbV incremental exposure increase in 6mDMA8, as presented in Figure 3. Subgroup meta-analysis by originally reported metrics concluded that risk significance varies across metrics. The metric 6mDMA8, which is characterized by high concentrations, had the highest positive risk (RR = 1.022, 95% CI: 1.014–1.030), while the smoothed metric ADA24 reported negative association (RR = 0.980, 95% CI: 0.960–1.001), as shown in Figure S1. When grouped by study region, significant risk pattern discrepancies were found (Figure S2). The studies in North America revealed positive associations, as RR = 1.019 (95% CI: 1.014–1.024), while European populations showed reversed risks, as RR = 0.910 (95% CI: 0.827–1.001), although not statistically significant. The cross-region divergence did not necessarily indicate differences in population vulnerability, because the European studies involved (1) smaller and younger study populations, (2) shorter follow-up durations, and (3) use of smoothed exposure metrics, which could all potentially obscure the risk associations. Subgroup analysis indicated that high inter-study heterogeneities originated from metric inconsistency, methodological reliability of individual exposure assignment, and population variabilities (Table S4). The funnel plot was visually symmetrical (Figure S3), and studies reporting risks below the pooled value were slightly greater in number, indicating no severe publication biases.

Figure 3

Pooled estimates of all-cause mortality risk associated with every 10-ppbV incremental O3 exposure by 6mDMA8 metric

No significant inter-gender differences were observed, based on the limited studies reporting gender-specific risk association strengths. Further subgroup analyses were unfeasible due to the lack of reporting in the literature. Alternatively, grouped RRs were estimated based on whether the original research had adjusted for confounding effects from ethnicity, BMI, smoking history, lifestyle features, and exposure levels of PM2.5 and NO2. No inter-group divergences were observed (Table S4).

Respiratory mortality

Meta-analysis for O3 exposure-associated all respiratory mortality included 16 studies, which gave the pooled RR = 1.025 (95% CI: 1.010–1.040, I2: 83.9%) 10-ppbV incremental O3 exposure increase by 6mDMA8 (Figure 4). Based on subgroup meta-analysis for different metrics (Figure S4), peak metrics showed more significant increasing risks than ADA24, which exhibited most of the heterogeneities (I2 = 87.0%). Cross-metric divergences were generally higher than intra-metric discrepancies. Studies on North American populations showed better homogeneity in positive risks (RR = 1.029, 95% CI: 1.011–1.047, I2 = 71.1%, Figure S5) than the European cohorts, pooling from which the overall risks were congruously insignificant (RR = 0.941, 95% CI: 0.856–1.036, I2 = 91.2%). For O3-COPD mortality association, the pooled RR was 1.056 (95% CI: 1.029–1.084, I2 = 94.5%) per 10-ppbV incremental O3 exposure increase by 6mDMA8 from seven studies. No apparent positive publication biases were detected for both respiratory and COPD mortalities from the funnel plot (Figure S3).

Figure 4

Pooled estimates of respiratory mortality risks associated with every 10-ppbV incremental O3 exposure by 6mDMA8 metric

Cardiovascular mortality

A total of 15 studies were included to pool the overall O3 exposure-induced CVD mortality risks as RR = 1.019 (95% CI: 1.004–1.035, I2 = 97.6%) per 10-ppbV additional O3 exposure by 6mDMA8 (Figure 5). Given that the lower uncertainty bound was so close to the null hypothesis (i.e., RR = 1), the positive risk association found in this review could be controversial, so it would require more studies to support or refute the finding. Heterogeneities (I2 > 79.2%) were observed through all three metric-grouped studies, as presented in Figure S6. Positive risk associations were found on 10 North American cohorts (RR = 1.036, 95% CI: 1.017–1.056) but negatively for five European cohorts (RR = 0.934, 95% CI: 0.866–1.008), as shown in Figure S7. There were no concerns with publication bias, and no more inter-group divergences were spotted except for grouping by exposure assignment methodological credibility (Table S4). The pooled risk for congestive heart failure-induced mortality from four studies was RR = 1.074 (95% CI: 1.054–1.093, I2 = 85.8%) per 10-ppbV incremental O3 exposure increase by 6mDMA8.

Figure 5

Pooled estimates of cardiovascular mortality risks associated with every 10-ppbV incremental O3 exposure by 6mDMA8 metric

Other mortality causes

The other cause-specific mortality risks attributable to long-term O3 exposure were not statistically significant (Figure 6), as IHD mortality risk pooled from 10 studies was RR = 1.012 (95% CI: 0.987–1.039, I2 = 98.7%), CEVD mortality risk pooled from six studies was RR = 0.993 (95% CI: 0.979–1.008, I2 = 80.6%), and LC mortality risk pooled from 12 studies was RR = 0.966 (95% CI: 0.926–1.007, I2 = 84.2%). For all eight studied mortality causes, we also provide pooled risks by three more widely used metrics (6mDA24, ADMA8, and ADA24) besides 6mDMA8, as listed in Table 3 for reference.

Figure 6

Pooled estimates of ischemic heart disease, cerebrovascular diseases, and lung cancer mortality risks associated with every 10-ppbV incremental O3 exposure by 6mDMA8 metric

Table 3

Pooled RRs for long-term 10-ppbV incremental O3 exposure-associated mortalities by four major metrics and crude risks without harmonization

Mortality causes	6mDMA8	6mDA24	ADMA8	ADA24	Crude
All causes (n = 23)	1.014 (1.009, 1.019)	1.023 (1.014, 1.032)	1.016 (1.010, 1.022)	1.027 (1.017, 1.037)	1.017 (1.011, 1.023)
Respiratory diseases (n = 16)	1.025 (1.010, 1.040)	1.042 (1.016, 1.069)	1.029 (1.011, 1.047)	1.049 (1.019, 1.081)	1.031 (1.017, 1.046)
Chronic obstructive pulmonary disease (n = 7)	1.056 (1.029, 1.084)	1.098 (1.050, 1.149)	1.066 (1.034, 1.098)	1.116 (1.058, 1.176)	1.055 (1.032, 1.078)
Cardiovascular diseases (n = 15)	1.019 (1.004, 1.035)	1.033 (1.006, 1.061)	1.022 (1.004, 1.041)	1.038 (1.007, 1.071)	1.024 (1.009, 1.038)
Ischemic heart disease (n = 10)	1.012 (0.987, 1.039)	1.021 (0.977, 1.067)	1.014 (0.984, 1.045)	1.024 (0.973, 1.078)	1.017 (0.994, 1.041)
Congestive heart failure (n = 4)	1.074 (1.054, 1.093)	1.130 (1.094, 1.168)	1.086 (1.063, 1.110)	1.155 (1.110, 1.198)	1.083 (1.059, 1.107)
Cerebrovascular diseases (n = 6)	0.993 (0.979, 1.008)	0.988 (0.964, 1.013)	0.992 (0.976, 1.009)	0.986 (0.958, 1.015)	0.992 (0.979, 1.006)
Lung cancer (n = 12)	0.966 (0.926, 1.007)	0.943 (0.878, 1.012)	0.960 (0.915, 1.008)	0.933 (0.859, 1.014)	0.960 (0.909, 1.013)

Study assessment

All 25 studies included in our final meta-analysis were rated above “Fair” (14 “Fair” and 11 “Good”) by the Quality Assessment Tool for Observational Cohort Studies, as listed in Table S5. All studies met 10 out of 14 assessment items. Nine studies did not sufficiently clarify the participant exclusion criteria; two re-analysis study reports did not clearly state the O3 exposures;, two studies were of insufficient follow-up durations (e.g., less than 5 years) to observe the outcomes resulting from long-term exposure;, and 10 studies were methodologically deficient in individual exposure assignment,38, 39, 40, 41, 42, 43, 44, 45, 46, most of which were conducted before 2013 when data assimilation techniques were not comprehensively developed to fuse observations with other full spatial coverage products such as satellite-based remote sensing and atmospheric mechanistic simulations. The satisfactory assessment results overall indicated indiscernible risks of bias, laying the reliable foundation for meta-analyses.

Table S6 displays GRADE epidemiological evidence assessment results for each mortality cause from all involved corresponding studies. In brief, the overall judgements for all-cause, respiratory, cardiovascular, ischemic heart disease, congestive heart failure, and lung cancer mortality risks were “High,” while the rating for the remaining two cause-specific mortality risks (COPD and cerebrovascular diseases) were both “Moderate.” Inconsistency of the risk directions (i.e., positive or negative associations) was the most common reason for downgrading, except for CHF-induced mortality. There were six studies that reported O3-mortality exposure-response trends to support the additional risks, thus upgrading the pooled RRs of all-cause, respiratory, and cardiovascular mortality. Cakmak et al. (2018) reported higher RRs after adjusting for confounders compared to the crude values, which gave prominence to the positive risk associations and hence upgraded the rating for all-cause, ischemic heart disease, and lung cancer mortalities. No substantial positive publication biases were found based on the collected evidences.

Sensitivity analysis

Leave-one-out successive elimination sensitivity analyses showed stable risk estimates as summarized in Table S7, except for lung cancer mortality risk after eliminating Kazemiparkouhi et al. (2019), the only study reporting positive risk association, while the other 11 studies concluded insignificant risks or even protective effects. Since the metric harmonization in our study was an innovative approach, we provided both metric-adjusted and unadjusted crude results for reference as presented in Table 3. The crude results were pooled from the originally reported relative risk values unified into per 10-ppbV incremental exposure, without being transformed into any metrics for congruity. Along with the meta-analyses on all qualified studies, the RR were also pooled by retaining only the latest study with the largest population for each separate cohort, as summarized in Table S8. With this approach, the pooled unit incremental mortality risks per 10-ppbV O3 exposure increase by 6mMDA8 metric were RR = 1.008 (95% CI: 1.006–1.009, I2 = 82.6%) for all causes, RR = 1.034 (95% CI: 1.017–1.050, I2 = 81.7%) for all respiratory diseases, RR = 1.060 (95% CI: 1.040–1.080, I2 = 90.2%) for COPD, RR = 1.032 (95% CI: 1.010–1.055, I2 = 98.2%) for all cardiovascular diseases, RR = 1.008 (95% CI: 0.973–1.045, I2 = 99.2%) for ischemic heart disease, and RR = 0.966 (95% CI: 0.931–1.002, I2 = 83.8%) for lung cancer. Studies for mortality risks by cerebrovascular diseases and congestive heart failure were respectively conducted on different cohorts, so such supplementary analysis was unnecessary.

Discussion

Improvements as an updated review

This work improves on two previous high-quality reviews, by covering up-to-date peer-reviewed studies, and expanding the O3 exposure-associated causes of mortality into a wider range of categories. To the best of our knowledge, it is the first systematic review of the association between long-term O3 exposure and cause-specific mortality that highlights the issue of inconsistent use of exposure metrics. Since tropospheric O3 is a photochemical pollutant that largely depends on solar radiation, surface O3 concentrations can vary drastically between day and night, as well as warmer and cooler seasons. We point out that a 10-ppbV increase in annual daily 24-h average concentration (ADA24) is more constrained in magnitude than a 10-ppbV increase in warm-season daily 8-h maximum average concentration (6mDMA8) owing to the wider variability in the range of the latter metric. Taking the observations by TOAR and CNEMC in situ monitoring networks during 1990–2019 as an example, the surface O3 concentrations were 27.6 ± 6.1 (IQR: 24.1–31.0) ppbV by ADA24, while correspondingly 53.1 ± 10.6 (IQR: 47.7–61.4) ppbV by 6mDMA8, which indicated that a 10-ppbV change fell below the IQR by the 6mDMA8, but it could exceed the IQR using the ADA24 metric. This was why exposure metric adjustment was necessary for O3 exposure-attributable health risk meta-analysis.

We also put forward an approach to mutually convert the O3 exposure concentrations and corresponding risk strengths in various metrics by non-intercept linear projections, following the operational guidelines from US EPA, but updated the linear conversion factors using global in situ surface O3 observations during 1990–2019. This methodological innovation took advantage of multi-dimensional information from the original studies, which could inspire future data collection and research for corroboration and improvements.

Metric relevant issues

Although linear coefficients were applied to the cross-metric conversions, irreducible noises still existed given the high root mean squared errors (RMSE), as shown in Figure 2. This exposed the limitation of risk strength adjustment into the same exposure metric by simple linear conversion, as the actual cross-metric relationships are likely to be more complicated. However, the only appropriate approach was to use linear conversion coefficients to unify the RRs reported by different metrics in the original studies. Therefore, to avoid uncertainties brought about by the conversion of metrics, future long-term O3 exposure epidemiology studies should use a promissory consistent exposure metric or estimate the unit excess RRs in multiple metrics.

Such linear conversion of risk associations could be validated by Kazemiparkouhi et al. (2020), where multiple metrics were applied to estimate the mortality risks. For COPD mortality, the RR was 1.072 (95% CI: 1.067–1.077) by 6mDMA1 for every 10-ppbV additional exposure. After converting into 6mDMA8 using the linear coefficient 0.831 (Figure 2), the estimated RR was 1.087 (95% CI: 1.081–1.093), where Kazemiparkouhi et al. (2020) reported 1.084 (95% CI: 1.079–1.089), which justified our linear conversion method. Cross-metric linear conversions would not change the risk association direction, but using different exposure metrics when estimating the O3 exposure-attributable mortality risks could potentially cause discrepancies. For instance, Kazemiparkouhi et al. (2020) concluded excess hazards of long-term O3 exposure on all-cause mortality using 6mDMA1 and 6mDMA8 as quantitative metrics. However, 6mDA24 led to a specious prevention effect (RR = 0.990, 95% CI: 0.988–0.991), which should be attributed to the existence of a theoretical exposure safety level for O3 below which no negative health effects should occur. Under this circumstance, lower-level metrics (e.g., ADA24) that average peak O3 exposures may obscure effective doses above the threshold, and also reduce the signal-to-noise ratios, resulting in lower credibility for recognizing hazardous population exposures than higher-level metrics (e.g., 6mDMA8).

Data mining techniques are able to realize high-accuracy predictions of surface O3 concentrations, but errors were never avoidable. Carey et al. (2013) used the basic inverse distance weighting (IDW) spatial interpolation approach to obtain surface O3 concentrations with R2 = 0.24–0.76, while years later, Di et al. (2017) applied an ensemble learning approach, achieving R2 = 0.80, RMSE = 2.91 ppbV. Carey et al. (2013) reported the IQR of O3 exposure concentrations as 3.0 ppbV, which was comparable to the RMSE of Di et al. (2017). Lower R2 values are typically accompanied by higher prediction errors, which may conceal the highest and lowest quartiles, and lead to failures in distinguishing the population-level exposures. This concern is reflected in our subgroup meta-analysis by exposure metrics, where lower-level metrics were more inclined to report insignificant risks. This also casts doubt on the reliability of studies covering narrow exposure variabilities. We therefore support the Lancet suggestions to use peak metrics to quantify long-term O3 exposure, such as 6mDMA8, and the use of state-of-the-art data techniques to reduce errors in O3 prediction, so as to make a distinction between high- and low-exposure populations.

Rebuttal to a previous review

Numerous pathogenesis mechanisms have been at least partially ascertained by laboratory experiments. Inhaled O3 can constrict muscles in the airways, leading to shortness of breath, and damage the lining via inflammation. Long-term O3 exposure can increase oxidative stress in the cardiovascular system and cause progressive thickening of the carotid arteries to restrict cerebral blood supply. Short-term O3 exposure has also been strongly associated with a variety of cardiopulmonary symptoms, as reported by a number of observational epidemiological studies. This also supports the long-term effects, presuming incremental risks by long-term exposure given the verified significant short-term effects. We therefore approve of the opinion that long-term O3 exposure increases mortality risks, in agreement with the GBD report.

However, Atkinson et al. (2016) concluded insignificant pooled risks for long-term O3 exposure-associated all-cause and respiratory mortality, which contradicted with our results. It should mainly be ascribed to the heterogeneity between the more recent studies and earlier ones. The majority of studies collected in Atkinson et al. (2016) applied primitive statistical methods (i.e., nearest neighborhood matching, IDW, and ordinary kriging interpolation) for individual exposure assignment, which might have weakened the individual-level exposure distinguishment. In addition, some studies using ADA24 as the exposure metric could have also oversmoothed the peak exposures so obscured the significance of associations.,, In contrast, studies after 2016 more frequently applied advanced numerical simulation models and data assimilation techniques to increase the precision of population exposure assessment; most of these used the 6mDMA8 metric to foreground high exposures.,,,55, 56, 57 These recent studies highlighted significant O3-mortality associations.

To alleviate the population health loss resulting from O3 exposure, the US EPA appealed for optimizations in real-time accessibility of an air quality index, with which residents could avoid unnecessary high pollution exposure (https://www.epa.gov/ground-level-ozone-pollution/health-effects-ozone-pollution). Appropriate diets and supplements including carotenoids, vitamin D, and vitamin E were recognized to be preventive against air pollution-induced respiratory damage, as a practical protective measure for vulnerable individuals. These emergent announcements and studies demonstrated a growing number of researchers inclined to the view that long-term O3 exposure should be considered a health hazard.

Concentration-response relationship

Few studies have examined the concentration-response curves between long-term O3 exposure and mortality, so the threshold exposure level (also known as theoretical minimum risk exposure level, TMREL), below which no adverse health effects would be assumed to occur, remains controversial. For all-cause mortality, Di et al. (2017) reported a safe exposure level as 30 ppbV by the 6mDMA8 metric, while Shi et al. (2021) suggested a more conservative level at 40 ppbV by 6mDMA8, both estimated from the Medicare beneficiary cohort. For respiratory mortality, Jerrett et al. (2009) tested the concentration-response relationships and estimated the threshold level as 60 ppbV by 6mDMA1 (approximately 49.9 ppbV by 6mDMA8), while Lim et al. (2019) failed to identify a significant threshold level. For cardiovascular mortality, Lim et al. (2019) showed no apparent health hazards below 45 ppbV by 6mDMA8, and Paul et al. (2020) prescribed a threshold level around 35 ppbV by the 6mDMA8 metric for diabetic patients. These evidence-based threshold exposure levels were all within the current standards: 70 ppbV for daily maximum 8-h exposure under NAAQS (The National Ambient Air Quality Standards regulated by the US EPA) and 50 ppbV by warm-season DMA8 under WHO global air quality guidelines. However, further studies should determine whether stricter standard guidelines are necessary.

To synthesize epidemiological evidence, Burnett et al. (2014) developed an integrated exposure-response (IER) function-based curve-fitting method to pool the risk associations from multiple studies. We constructed the IER for long-term O3 exposure-associated mortalities in this review, with statistically reproduced exposure levels to enhance the curve fitting, as illustrated in supplemental information S1–S3. The exposure imputation revealed high reliability, but the high uncertainties of the estimated IER curves could not be addressed. This is attributed to the limited availability of effective epidemiological evidence. Empirically, this approach would require sufficient studies to cover a wide range of exposure levels, as has been frequently adopted for particulate matter exposure research,70, 71, 72 but seldom used for O3 health studies. The reason might be that the population long-term O3 exposure levels are not as comparably distinguishable as those of particulate matter. In addition, a reasonable prescribed TMREL would be necessary to establish the IER curves; the indeterminacy of the threshold level could exacerbate uncertainties in the estimated concentration-response trends. Therefore, we urgently appeal for more relevant studies on long-term O3 exposure-associated risks, and encourage in-depth discussions, optimizations, or corrections on our enhanced exposure-response relationship curve-fitting methodologies.

Hierarchical classification of diseases

The causes of mortality analyzed in our study followed hierarchical subordinate relationships, as the all-cause mortality consisted of cardiovascular diseases, respiratory diseases, cancer, and other causes; chronic obstructive pulmonary disease belonged to the respiratory category; and ischemic heart disease, congestive heart failure, and cerebrovascular diseases were all subordinated to cardiovascular symptoms. On this occasion, estimating all O3 exposure-induced mortalities could follow a bottom-up scheme by adding up subgroups of diseases. However, for historical O3-associated mortalities, GBD attributed all O3-associated mortalities onto COPD-induced premature deaths, which could be erroneous Long-term O3 exposure has shown significant association with excess cardiovascular mortalities, so mortality estimations in future studies should also include CVD.

Applications for mortality estimation

The widest applications of estimated risk association strengths have been to project how many people could be affected by long-term ambient O3 exposure. For example, Malley et al. (2017) estimated 1.23 (95% UI: 0.85–1.62) million respiratory deaths attributable to O3 exposure in 2010, using the risk strength by Turner et al. (2016) as HR = 1.12 (95% UI: 1.08–1.16). This estimation was much higher than the 2019 GBD reported figure of 0.31 (95% UI: 0.15–0.49) million, as highlighted in another recent study. This should be attributed to the use of high HR value among all included studies. We also found other studies using a singular HR value for population risk estimations,74, 75, 76, 77 but would encourage further relevant studies to consider multi-study pooled RRs. This could effectively reduce the potential biases from a single study. The adaptability of pooled RRs could be verified from the coverage of exposure levels, as the 25 studies identified in our review had embraced a wide range of exposure concentrations (supplemental information S1) to encompass global surface O3 variability. On the other hand, the leave-one-out sensitivity analyses (Table S7) revealed the robustness of the meta-analysis results when including sufficient numbers of studies. This reflects the representativeness of the synthesized risk association strengths. The annual GBD reports also presumed generalizability of the synthesized epidemiological evidences, but cohort-based research in unstudied regions remains a requirement for more convincing conclusions.

Limitations

Although the total number of study participants for risk pooling was adequately high to ensure statistical power, the cohort-based O3 health studies were factually rare according to our literature search, so long-term follow-up studies are urgently encouraged. Furthermore, few studies reported grouped RRs (gender, age, socio-economic status, smoking and alcohol history, etc.), which made meta-analyses by sub-categories unfeasible. High inter-study heterogeneities existed across all mortality causes, which were not substantially alleviated by subgroup clustering. This indicates that the causes for cross-study incongruity have not been fully identified based on current evidence, so future studies are strongly encouraged to help identify the causes of the heterogeneity. Scarcity of credible evidence also restricted the effects of conventional approaches to construct exposure-response curves, and our methodological innovation requires further relevant studies for substantiation. The cross-metric linear conversion factors were estimated using observations from available sites, which may not be sufficiently representative of global residential areas, since observational sites in India, Africa, and Latin America were still sparse. With ever-increasing deployment of in situ monitoring networks, the cross-metric conversion factors could be calibrated with more comprehensive observations; the pooled RRs should also be updated accordingly. In principle, the meta-analysis pooled RRs with exposure metric adjustment were reasonable estimations for the unknown true values of the effect size, requiring more relevant studies to support or correct. Prospective researchers should be aware of this limitation and use the meta-analysis results cautiously.

Further study suggestions

We suggest that further environmental epidemiology studies, especially long-term O3 exposure-related research, clearly report (1) the methodologies used to obtain ambient O3 concentrations, the spatiotemporal resolution, and prediction accuracy of the database; (2) the exposure metrics used for risk estimation; and (3) the statistical distribution of the O3 exposure concentrations. The data-oriented methodologies used to accomplish full spatial coverage ambient air O3 concentrations for individual-level exposure assignment should be transparent so that the construction credibility of air pollution concentration databases can be rigorously assessed, thus forming the foundation of epidemiological follow-up studies. We advocate the reporting of exposure metrics in future O3 health studies to avoid confusion when comparing the risks with literature and conducting meta-regression; according to recent consensus, warm-season average (6mDMA8) is preferred for epidemiological study metrics. We recommend that future studies estimate risks with multiple O3 metrics for reference and describe the statistical distribution of O3 exposure levels to assess the reliability of risk estimation models. This can be useful in exposure-response tendency exploration. We also propose that future cohort studies estimate subgroup-specific RRs, which can help to identify vulnerable populations.

Our review highlights a deficiency in current environmental health research literature: studies on long-term O3 exposure health effects are still rather rare compared to particulate matter-based studies. Additionally, the meta-analysis results may be geographically biased, since large-scale O3 exposure health risk studies up to 2022 did not cover Asian, African, or Latin American regions. However, there are some thriving cohort projects such as the Multi-Country Multi-City Collaborative Research Network covering over 22 countries or regions, and the China Kadoorie Biobank enrolling over 0.5 million people, enabling environmental exposure research. We are optimistic that more research will fill the literature gap of multi-region population-based studies.

Conclusions

Our updated systematic review has summarized cohort studies exploring the associations between long-term ambient O3 exposure and multi-cause mortality risks. Current studies support O3 exposure-attributable excess mortalities from all causes, respiratory diseases, chronic obstructive pulmonary disease, cardiovascular diseases, and congestive heart failure, but no significant mortality risks are found for ischemic heart diseases, all-type cerebrovascular diseases, and lung cancer. Exposure metrics are crucial for health risk estimations of long-term O3 exposure and meta-analysis to pool the multi-study risks, for which we develop a linear conversion approach to harmonize the different metrics. Further research on long-term O3 observations and exposure-induced mortalities are encouraged to corroborate or contradict our linear conversion factors and meta-analysis results by providing additional evidence to strengthen the O3 health literature.