Hepatitis C care cascade among patients with and without tuberculosis: Nationwide observational cohort study in the country of Georgia, 2015–2020Baliashvili, Davit;Blumberg, Henry M.;Gandhi, Neel R.;Averhoff, Francisco;Benkeser, David;Shadaker, Shaun;Gvinjilia, Lia;Turdziladze, Aleksandre;Tukvadze, Nestani;Chincharauli, Mamuka;Butsashvili, Maia;Sharvadze, Lali;Tsertsvadze, Tengiz;Zarkua, Jaba;Kempker, Russell R.
doi: 10.1371/journal.pmed.1004121pmid: 37141386
Background The Eastern European country of Georgia initiated a nationwide hepatitis C virus (HCV) elimination program in 2015 to address a high burden of infection. Screening for HCV infection through antibody testing was integrated into multiple existing programs, including the National Tuberculosis Program (NTP). We sought to compare the hepatitis C care cascade among patients with and without tuberculosis (TB) diagnosis in Georgia between 2015 and 2019 and to identify factors associated with loss to follow-up (LTFU) in hepatitis C care among patients with TB. Methods and findings Using national ID numbers, we merged databases of the HCV elimination program, NTP, and national death registry from January 1, 2015 to September 30, 2020. The study population included 11,985 adults (aged ≥18 years) diagnosed with active TB from January 1, 2015 through December 31, 2019, and 1,849,820 adults tested for HCV antibodies between January 1, 2015 and September 30, 2020, who were not diagnosed with TB during that time. We estimated the proportion of patients with and without TB who were LTFU at each step of the HCV care cascade and explored temporal changes. Among 11,985 patients with active TB, 9,065 (76%) patients without prior hepatitis C treatment were tested for HCV antibodies, of which 1,665 (18%) had a positive result; LTFU from hepatitis C care was common, with 316 of 1,557 (20%) patients with a positive antibody test not undergoing viremia testing and 443 of 1,025 (43%) patients with viremia not starting treatment for hepatitis C. Overall, among persons with confirmed viremic HCV infection, due to LTFU at various stages of the care cascade only 28% of patients with TB had a documented cure from HCV infection, compared to 55% among patients without TB. LTFU after positive antibody testing substantially decreased in the last 3 years, from 32% among patients diagnosed with TB in 2017 to 12% among those diagnosed in 2019. After a positive HCV antibody test, patients without TB had viremia testing sooner than patients with TB (hazards ratio [HR] = 1.46, 95% confidence intervals [CI] [1.39, 1.54], p < 0.001). After a positive viremia test, patients without TB started hepatitis C treatment sooner than patients with TB (HR = 2.05, 95% CI [1.87, 2.25], p < 0.001). In the risk factor analysis adjusted for age, sex, and case definition (new versus previously treated), multidrug-resistant (MDR) TB was associated with an increased risk of LTFU after a positive HCV antibody test (adjusted risk ratio [aRR] = 1.41, 95% CI [1.12, 1.76], p = 0.003). The main limitation of this study was that due to the reliance on existing electronic databases, we were unable to account for the impact of all confounding factors in some of the analyses. Conclusions LTFU from hepatitis C care after a positive antibody or viremia test was high and more common among patients with TB than in those without TB. Better integration of TB and hepatitis C care systems can potentially reduce LTFU and improve patient outcomes both in Georgia and other countries that are initiating or scaling up their nationwide hepatitis C control efforts and striving to provide personalized TB treatment. Why was this study done? There is ample evidence that hepatitis C prevalence is disproportionally high among patients with tuberculosis (TB). Highly effective new treatment options for hepatitis C allowed many countries, including Georgia, to implement large-scale hepatitis C programs. It has not been well characterized how often patients with current or past TB are offered and provided with hepatitis C testing and treatment services. What did the researchers do and find? We conducted an observational cohort study comparing the hepatitis C care cascade among patients with and without TB to explore if patients with tuberculosis receive hepatitis C treatment completely and timely. The proportion of patients with TB tested for hepatitis C virus (HCV) antibodies increased per year. Among patients diagnosed with TB in 2015, 60% were tested for HCV antibodies sometime during the study period. This proportion reached 90% among patients diagnosed with TB in 2019 Loss to follow-up (LTFU) from hepatitis C care was more common among patients with TB, with 20% of patients with a positive antibody test not undergoing viremia testing and 43% of patients with viremia not starting treatment for hepatitis C. For comparison, the respective numbers among patients without TB were 14% and 19%. What do these findings mean Our findings highlight the importance of improving integration and linkage to hepatitis C diagnostic and treatment services among patients with TB. Existing large-scale public health programs for both TB and hepatitis C in Georgia and other countries with nationwide programs create a unique opportunity for integrated care of these 2 infectious diseases, which could potentially reduce LTFU and improve overall health outcomes. Introduction Both tuberculosis (TB) and hepatitis C virus (HCV) infection cause substantial morbidity and mortality worldwide. In 2020, there were an estimated 10 million new cases of active TB globally and 1.5 million deaths due to TB [1]. There are an estimated 59 million people living with chronic HCV infection, 1.5 million new cases occur each year, and 290,000 people died due to hepatitis C globally in 2019 [2]. TB and HCV infection are often concentrated in the same high-risk population subgroups, such as homeless or incarcerated individuals and people who inject drugs [3–5]. However, it has not been well characterized how often patients with current or past TB are offered and provided with hepatitis C testing and treatment services. Hepatitis C diagnosis, treatment, and assessment of cure consist of multiple consecutive steps, collectively known as the hepatitis C care cascade [6]. The cascade starts with hepatitis C antibody testing and ends with the determination of cure based on sustained virologic response (SVR) at 12 weeks after treatment completion [7–16]. Available data has found that loss to follow-up (LTFU) at different stages of the hepatitis C care cascade is common [17]. In 2019, globally, only 26% of people living with hepatitis C knew their infection status, and only 16% were treated [2]. The hepatitis C care cascade among patients with both TB and hepatitis C has not been previously studied. However, given the increased rate of adverse events including hepatoxicity among patients with TB and HCV coinfection, screening and potentially early treatment should be a high priority in this group [18,19]. Historically, patients receiving treatment for TB have not been recommended to receive hepatitis C treatment until they finish TB treatment due to immunosuppressive effect of interferon or due to drug–drug interactions of more recently introduced directly-acting antivirals (DAAs) [9,20]. However, the combined effect of TB drugs and HCV infection can cause long-term impairment of liver function, and timely management of hepatitis C in patients with TB as soon as they become eligible should be an important research priority. Since patients with TB are already receiving care, there is a key opportunity to engage them in hepatitis C care without interruption, reducing the likelihood of LTFU. The Eastern European country of Georgia (population 3.7 million) is designated by the WHO as a high-priority country for TB control in European region, with an estimated incidence of 70 TB cases per 100,000 population in 2020 and a high percentage of multidrug-resistant (MDR) TB [1,21]. Chronic hepatitis C is also highly prevalent in Georgia, affecting an estimated 5.4% of the general adult population (approximately 150,000 individuals) based on a 2015 serosurvey—the highest prevalence among countries of Eastern Europe and Central Asia [22,23]. Historically, there has been a substantial overlap in the population affected by these 2 infectious diseases in Georgia; previous studies reported that >20% of newly diagnosed active TB cases had HCV antibodies, indicative of current or past infection [18,24]. To address the high burden of HCV infection, Georgia initiated the first nationwide hepatitis C elimination program in 2015, including free hepatitis C testing and treatment with new DAAs for all citizens [25–27]. The cost of confirmatory viremia testing was initially fully covered by the government only for persons with low income but became free for everyone beginning in March 2018 [28]. Free hepatitis C screening through antibody testing was gradually integrated into multiple existing programs and settings, including the National TB Program (NTP), which provides free diagnostic and treatment services for TB countrywide for residents of Georgia [26,28–31]. Before 2018, HCV antibody testing was performed routinely among patients with MDR TB and sporadically (by indication) among patients with drug-susceptible (DS) TB. According to the TB treatment guideline in Georgia adopted in July 2018, all patients newly diagnosed with TB should be routinely tested for HCV antibodies. If an antibody test is positive, blood samples are sent to the National Center for Disease Control and Public Health for HCV viremia testing. However, currently there is no formal referral system in place linking viremia–positive patients with TB to hepatitis C care. The objectives of this study were to (1) compare the hepatitis C care cascade among patients with and without TB diagnosis in Georgia between 2015 and 2019. We hypothesized that the proportion of HCV seropositive individuals who undergo HCV confirmatory testing, initiate and complete HCV treatment was lower among patients with TB compared to patients without TB; and (2) identify factors associated with LTFU in hepatitis C care among patients with TB. Findings from this analysis will be a valuable contribution to the Georgian Hepatitis C Elimination Program and other similar settings by providing data for future targeted interventions to improve linkage to and retention in hepatitis C care among patients with TB. Methods Study design and setting We conducted an observational cohort study and analyzed the hepatitis C care cascade among patients diagnosed with active TB disease and compared it to the care cascade among people with HCV without active TB disease. TB-related information was obtained from the NTP of Georgia. Hepatitis C screening, diagnostic, and treatment information was obtained from the Georgian Hepatitis C Elimination Program. Study population The study population consisted of 2 groups: (1) Adults (aged ≥18 years) diagnosed with active TB through the Georgian NTP from January 1, 2015 through December 31, 2019; and (2) adults tested for HCV antibodies in Georgia between January 1, 2015 and September 30, 2020, who were not diagnosed with TB during that time. We did not include patients diagnosed with TB in 2020 since they would have very limited follow-up time for HCV-related procedures. For both groups, hepatitis C testing and treatment information was obtained through September 30, 2020. Hepatitis C care cascade steps and definitions We assessed the following steps in the hepatitis C care cascade (Fig 1): (1) Screening for HCV antibodies using either rapid or laboratory-based tests; (2) viremia testing via ribonucleic acid (RNA) or core-antigen testing to confirm viremic HCV infection among those with HCV antibody positive result; (3) hepatitis C treatment initiation among those with positive viremia test; (4) hepatitis C treatment completion; and (5) SVR testing from 12 to 24 weeks after treatment completion. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. Steps of hepatitis C care cascade in Georgian Hepatitis C Elimination Program. DAA, directly-acting antiviral; HCV, hepatitis C virus; PCR, polymerase chain reaction; RNA, ribonucleic acid; SVR, sustained virologic response. https://doi.org/10.1371/journal.pmed.1004121.g001 Two separate care cascade analyses were performed among persons with HCV infection with and without TB disease. The hepatitis C care cascade among patients with TB included the following groups: (1) Persons tested for HCV antibodies at the time of or any time after their first TB diagnosis; and (2) persons tested for HCV antibodies before their TB diagnosis who did not initiate treatment for hepatitis C before their TB diagnosis, i.e., people who were eligible for hepatitis C care by the time they were enrolled in TB treatment. Hepatitis C care cascade among patients without TB included the following groups: (1) Persons tested for HCV antibodies who were not diagnosed with TB during 2015 to 2019; and (2) persons diagnosed with and treated for chronic HCV infection before they were first diagnosed with TB (thus, did not have TB at the time of their HCV diagnosis and treatment). Patients were defined as LTFU after a positive antibody test if they had a positive HCV antibody test but did not undergo HCV viremia testing within the HCV elimination program and were still alive at the time of data extraction (September 30, 2020). Patients were defined as LTFU after an HCV viremia test if they had a positive result on viremia testing (i.e., positive RNA or core antigen test) and did not start treatment for HCV infection within the HCV elimination program. Due to a lack of complete drug-susceptibility testing data, patients with TB were defined as having DR TB if they received TB treatment with second-line drugs which are used to treat MDR-TB. Sources of data Four data sources were utilized. We obtained all TB-related data from the Georgian NTP surveillance database. This database contains diagnostic and treatment-related information on every patient enrolled in the NTP, including penitentiary system. HCV screening information was obtained from the national HCV screening registry—a real-time, nationwide web-based system. The screening registry collects data from all sites providing HCV screening throughout the country, including TB facilities. HCV viremia testing and treatment information were obtained from the Hepatitis C Elimination Program clinical database, called “Elimination C” (ElimC) [32]. ElimC collects information related to viremia testing, other diagnostic procedures, and antiviral treatment from all service providers, except the medical facilities in penitentiary system. Patients included in any of these databases were cross-checked in the National Death Registry, and death dates were obtained for deceased patients. A unique national ID number was used as a linking variable between all data sources. Statistical analysis We conducted descriptive statistics of patients with and without TB to describe available demographic and clinical characteristics (Tables 1–3). Due to a limited number of variables in HCV screening registry, we could not conduct extensive comparison of the 2 groups. In the care cascade analysis, we first calculated the proportions of people reaching each step of the hepatitis C care cascade, as described above. Next, we used survival analysis methods to explore if the time from a positive HCV antibody test to viremia testing and time from a positive viremia test to hepatitis C treatment initiation differed between patients with and without TB. For patients with TB who had positive viremia test during their TB treatment, follow-up time started at the end of their TB treatment. Patients were censored at the end of the follow-up period (September 30, 2020). We used subdistribution hazards model to calculate hazards ratios (HR) and assess the differences between groups numerically. Death was treated as a competing risk. We also created cumulative incidence curves to examine the differences graphically, and Gray’s test for equality of cumulative incidence functions was used to test for differences [33,34]. Among those who underwent viremia testing, median time from antibody testing to viremia testing was compared using Wilcoxon rank-sum test. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 1. Demographic and clinical characteristics of patients diagnosed with TB, Georgia, January 1 2015–December 31, 2019. https://doi.org/10.1371/journal.pmed.1004121.t001 Download: PPT PowerPoint slide PNG larger image TIFF original image Table 2. Baseline demographic characteristics of persons tested for HCV antibodies, Georgia, January 1, 2015–September 30, 2020. https://doi.org/10.1371/journal.pmed.1004121.t002 Download: PPT PowerPoint slide PNG larger image TIFF original image Table 3. Baseline demographic characteristics of persons tested positive for HCV antibodies, Georgia, January 1, 2015–September 30, 2020. https://doi.org/10.1371/journal.pmed.1004121.t003 An additional analysis was performed among patients with TB to identify demographic or TB-related factors associated with LTFU in hepatitis C care after a positive antibody or viremia test. The analysis of factors associated with LTFU after a positive antibody test included all patients with TB (both DS and DR) who had a positive HCV antibody test result. The analysis of factors associated with LTFU after a positive viremia test was restricted to patients with DS TB who completed their TB treatment. Patients still on DS TB treatment were excluded because they are ineligible for hepatitis C treatment until after TB treatment is completed. Patients with DR TB were excluded due to heterogeneous TB treatment duration and lack of a standardized approach regarding the timing of hepatitis C treatment initiation. Log-binomial regression was used to calculate adjusted risk ratios (aRR) and 95% confidence intervals (CI). Covariate selection was conducted a priori based on directed acyclic graph (DAG) theory (S1–S4 Figs) [35]. Different sets of covariates were used for each of the patient characteristics. Significance level of 0.05 was used in all analyses that involved significance testing. To address the missing data in some of the covariates, we conducted multiple imputation with fully conditional specification using the same covariates as in the substantive analysis (20 imputations) [36]. All analyses were performed using SAS version 9.4. Ethics The study was approved by the institutional review boards (IRBs) at the National Center for Disease Control and Public Health (Tbilisi, Georgia), the National Center for Tuberculosis and Lung Disease (Tbilisi, Georgia), and Emory University (Atlanta, Georgia, United States of America). CDC’s Human Subjects Research Office determined CDC’s role to be provision of technical assistance that did not constitute engagement in the research. This study is reported as per the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) guideline (S1 Checklist). Inclusivity in global research Additional information regarding the ethical, cultural, and scientific considerations specific to inclusivity in global research is included in the Supporting information (S2 Checklist). Study design and setting We conducted an observational cohort study and analyzed the hepatitis C care cascade among patients diagnosed with active TB disease and compared it to the care cascade among people with HCV without active TB disease. TB-related information was obtained from the NTP of Georgia. Hepatitis C screening, diagnostic, and treatment information was obtained from the Georgian Hepatitis C Elimination Program. Study population The study population consisted of 2 groups: (1) Adults (aged ≥18 years) diagnosed with active TB through the Georgian NTP from January 1, 2015 through December 31, 2019; and (2) adults tested for HCV antibodies in Georgia between January 1, 2015 and September 30, 2020, who were not diagnosed with TB during that time. We did not include patients diagnosed with TB in 2020 since they would have very limited follow-up time for HCV-related procedures. For both groups, hepatitis C testing and treatment information was obtained through September 30, 2020. Hepatitis C care cascade steps and definitions We assessed the following steps in the hepatitis C care cascade (Fig 1): (1) Screening for HCV antibodies using either rapid or laboratory-based tests; (2) viremia testing via ribonucleic acid (RNA) or core-antigen testing to confirm viremic HCV infection among those with HCV antibody positive result; (3) hepatitis C treatment initiation among those with positive viremia test; (4) hepatitis C treatment completion; and (5) SVR testing from 12 to 24 weeks after treatment completion. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. Steps of hepatitis C care cascade in Georgian Hepatitis C Elimination Program. DAA, directly-acting antiviral; HCV, hepatitis C virus; PCR, polymerase chain reaction; RNA, ribonucleic acid; SVR, sustained virologic response. https://doi.org/10.1371/journal.pmed.1004121.g001 Two separate care cascade analyses were performed among persons with HCV infection with and without TB disease. The hepatitis C care cascade among patients with TB included the following groups: (1) Persons tested for HCV antibodies at the time of or any time after their first TB diagnosis; and (2) persons tested for HCV antibodies before their TB diagnosis who did not initiate treatment for hepatitis C before their TB diagnosis, i.e., people who were eligible for hepatitis C care by the time they were enrolled in TB treatment. Hepatitis C care cascade among patients without TB included the following groups: (1) Persons tested for HCV antibodies who were not diagnosed with TB during 2015 to 2019; and (2) persons diagnosed with and treated for chronic HCV infection before they were first diagnosed with TB (thus, did not have TB at the time of their HCV diagnosis and treatment). Patients were defined as LTFU after a positive antibody test if they had a positive HCV antibody test but did not undergo HCV viremia testing within the HCV elimination program and were still alive at the time of data extraction (September 30, 2020). Patients were defined as LTFU after an HCV viremia test if they had a positive result on viremia testing (i.e., positive RNA or core antigen test) and did not start treatment for HCV infection within the HCV elimination program. Due to a lack of complete drug-susceptibility testing data, patients with TB were defined as having DR TB if they received TB treatment with second-line drugs which are used to treat MDR-TB. Sources of data Four data sources were utilized. We obtained all TB-related data from the Georgian NTP surveillance database. This database contains diagnostic and treatment-related information on every patient enrolled in the NTP, including penitentiary system. HCV screening information was obtained from the national HCV screening registry—a real-time, nationwide web-based system. The screening registry collects data from all sites providing HCV screening throughout the country, including TB facilities. HCV viremia testing and treatment information were obtained from the Hepatitis C Elimination Program clinical database, called “Elimination C” (ElimC) [32]. ElimC collects information related to viremia testing, other diagnostic procedures, and antiviral treatment from all service providers, except the medical facilities in penitentiary system. Patients included in any of these databases were cross-checked in the National Death Registry, and death dates were obtained for deceased patients. A unique national ID number was used as a linking variable between all data sources. Statistical analysis We conducted descriptive statistics of patients with and without TB to describe available demographic and clinical characteristics (Tables 1–3). Due to a limited number of variables in HCV screening registry, we could not conduct extensive comparison of the 2 groups. In the care cascade analysis, we first calculated the proportions of people reaching each step of the hepatitis C care cascade, as described above. Next, we used survival analysis methods to explore if the time from a positive HCV antibody test to viremia testing and time from a positive viremia test to hepatitis C treatment initiation differed between patients with and without TB. For patients with TB who had positive viremia test during their TB treatment, follow-up time started at the end of their TB treatment. Patients were censored at the end of the follow-up period (September 30, 2020). We used subdistribution hazards model to calculate hazards ratios (HR) and assess the differences between groups numerically. Death was treated as a competing risk. We also created cumulative incidence curves to examine the differences graphically, and Gray’s test for equality of cumulative incidence functions was used to test for differences [33,34]. Among those who underwent viremia testing, median time from antibody testing to viremia testing was compared using Wilcoxon rank-sum test. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 1. Demographic and clinical characteristics of patients diagnosed with TB, Georgia, January 1 2015–December 31, 2019. https://doi.org/10.1371/journal.pmed.1004121.t001 Download: PPT PowerPoint slide PNG larger image TIFF original image Table 2. Baseline demographic characteristics of persons tested for HCV antibodies, Georgia, January 1, 2015–September 30, 2020. https://doi.org/10.1371/journal.pmed.1004121.t002 Download: PPT PowerPoint slide PNG larger image TIFF original image Table 3. Baseline demographic characteristics of persons tested positive for HCV antibodies, Georgia, January 1, 2015–September 30, 2020. https://doi.org/10.1371/journal.pmed.1004121.t003 An additional analysis was performed among patients with TB to identify demographic or TB-related factors associated with LTFU in hepatitis C care after a positive antibody or viremia test. The analysis of factors associated with LTFU after a positive antibody test included all patients with TB (both DS and DR) who had a positive HCV antibody test result. The analysis of factors associated with LTFU after a positive viremia test was restricted to patients with DS TB who completed their TB treatment. Patients still on DS TB treatment were excluded because they are ineligible for hepatitis C treatment until after TB treatment is completed. Patients with DR TB were excluded due to heterogeneous TB treatment duration and lack of a standardized approach regarding the timing of hepatitis C treatment initiation. Log-binomial regression was used to calculate adjusted risk ratios (aRR) and 95% confidence intervals (CI). Covariate selection was conducted a priori based on directed acyclic graph (DAG) theory (S1–S4 Figs) [35]. Different sets of covariates were used for each of the patient characteristics. Significance level of 0.05 was used in all analyses that involved significance testing. To address the missing data in some of the covariates, we conducted multiple imputation with fully conditional specification using the same covariates as in the substantive analysis (20 imputations) [36]. All analyses were performed using SAS version 9.4. Ethics The study was approved by the institutional review boards (IRBs) at the National Center for Disease Control and Public Health (Tbilisi, Georgia), the National Center for Tuberculosis and Lung Disease (Tbilisi, Georgia), and Emory University (Atlanta, Georgia, United States of America). CDC’s Human Subjects Research Office determined CDC’s role to be provision of technical assistance that did not constitute engagement in the research. This study is reported as per the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) guideline (S1 Checklist). Inclusivity in global research Additional information regarding the ethical, cultural, and scientific considerations specific to inclusivity in global research is included in the Supporting information (S2 Checklist). Results Participants A total of 14,993 records were obtained from the NTP database; 913 records (6.1%) were excluded due to a missing national ID number (Fig 2). The remaining 14,080 records corresponded to 12,767 individual patients, of whom 11,985 (94%) were aged ≥18 years and included in the analysis. Among those included (n = 11,985), 70% were males, 82% were unemployed at the time of TB diagnosis, 6.9% reported a history of incarceration, 11% had DR TB, 91% had 1 episode of TB treatment during the 2015 to 2019 period (i.e., diagnosed with TB only once), 7.9% had 2 episodes of TB treatment, and 1.6% had 3 to 8 episodes of TB treatment (Table 1). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. Flow chart of patients with TB and their HCV screening status—2015–2019. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g002 The hepatitis C screening registry contained data on 1,849,820 adults tested for HCV antibodies between January 1, 2015 and September 30, 2020, who were not diagnosed with TB in that period. The median age in this population was 46 years (interquartile range [IQR]: 31 to 62), and 54.8% were male (Tables 2 and 3). Hepatitis C care cascade among persons with and without TB After linking HCV testing information from the hepatitis C databases, 9,341 (78%) of 11,985 persons with TB were found to have been tested for HCV antibodies between January 1, 2015 and September 30, 2020 (Table 4). The proportion of patients with TB tested for HCV antibodies increased per year. Among patients diagnosed with TB in 2015, 60% were tested for HCV antibodies sometime during the study period. This proportion increased each subsequent year and reached 90% among patients diagnosed with TB in 2019. In the same period, HCV antibody-positivity in this population decreased from 25% among people diagnosed with TB in 2015 to 14% among those diagnosed in 2019. Among persons with TB included in our study, 276 (3.0%) had already started treatment for HCV infection before their first TB diagnosis and were included in the care cascade analysis of patients without TB. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 4. HCV testing by year of first TB diagnosis among adult patients diagnosed with TB in 2015–2019. https://doi.org/10.1371/journal.pmed.1004121.t004 Among the total of 9,065 patients with TB who were tested for hepatitis C, 1,665 (18%) patients (Fig 2) had a positive antibody test, of whom 108 (6%) died without undergoing viremia testing. Of the remaining 1,557 patients still alive, 316 (20%) were LTFU (i.e., never underwent viremia testing) (Figs 3A and S5); among 1,241 patients who underwent viremia results, viremic HCV infection was found in 1,025 (83%) patients, of whom 443 (43%) did not start hepatitis C treatment, including 80 persons that died (8% out of all viremic cases, 18% of those who did not start hepatitis C treatment). In a subset of patients with DR TB, the proportion of patients who were LTFU was even higher—27% and 52% after positive antibody and viremia tests, respectively (S6 Fig). For comparison, in the hepatitis C care cascade among patients without TB (Fig 3B), the proportion of patients who were LTFU from hepatitis C care after positive antibody test or after positive viremia test were 14% and 19%, respectively. Mortality was also lower among patients without TB: 3,842 (3%) of HCV antibody-positive persons died without viremia testing, and 3,127 (4%) viremic patients died without received hepatitis C treatment. Among patients initiating HCV treatment, those with TB were also less likely to complete hepatitis C treatment (89% versus 94%) and SVR assessment (59% versus 76%) than those without TB. However, among those who completed hepatitis C treatment and were tested for SVR, the cure rate was comparable to those without TB (98.3% versus 98.9%) (Figs 3A and 3B and S7 and S1 Table). Overall, among patients with confirmed viremic HCV infection, patients with TB were less likely to complete hepatitis C treatment (51%) and have a documented SVR achieved (28%), compared to patients without TB (76% and 55%, respectively; Fig 3C). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Hepatitis C care cascade among patients with and without tuberculosis, Georgia, 2015–2020. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. Note: Red lines represent the percent change between 2 consecutive steps in the care cascade, i.e., adjacent bars of the charts. https://doi.org/10.1371/journal.pmed.1004121.g003 Timeliness of HCV viremia testing and treatment initiation Among those who had information on follow-up duration available, we used survival analysis methods to compare time from positive antibody test to viremia testing or censoring between patients with TB (n = 1,499) and patients without TB (n = 118,396). After positive HCV antibody test, 3,840 (3.2%) persons died, and patients without TB had viremia testing sooner than patients with TB (HR = 1.46, 95% CI [1.39, 1.54], p < 0.001) (Figs 4 and S8). Among those who underwent viremia testing, the median time from a positive antibody test result to viremia testing was 17 days (IQR: 3 to 248 days) among patients with TB (n = 1,116) and 6 days (IQR: 1 to 24) among patients without TB (n = 97,554) (p < 0.001). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Comparison of HCV viremia testing cumulative incidence among HCV seropositive patients with and without TBa,b: January 1, 2015–September 30, 2020. a Blue line corresponds to the patients with confirmed HCV infection without evidence of TB diagnosis in 2015–2019. Red line corresponds to patients with confirmed HCV infection diagnosed with TB in 2015–2019. b Gray’s test p-value <0.001. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g004 We also compared time from positive viremia test to hepatitis C treatment initiation or censoring among 1,021 patients with TB and 87,964 patients without TB who had active viremia and information on follow-up duration available. After a positive viremia test, patients without TB started hepatitis C treatment sooner than patients with TB (HR = 2.05, 95% CI [1.87, 2.25], p < 0.001) (Figs 5 and S9). Among those who started hepatitis C treatment, the median time from confirmation to treatment initiation was 196 days (IQR: 76 to 523) among patients with MDR TB (n = 78), 99 days (IQR: 54 to 263) among patients with DS TB (n = 460), and 86 days (IQR: 56 to 203) among patients without TB (n = 71,128). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 5. Comparison of cumulative incidence of HCV treatment initiation among patients with confirmed HCV infection with and without TBa,b: January 1, 2015–September 30, 2020. a Blue line corresponds to patients with confirmed HCV infection without evidence of TB diagnosis in 2015–2019. Red line corresponds to patients with confirmed HCV infection diagnosed with TB in 2015–2019. b Gray’s test p-value <0.001. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g005 Factors associated with LTFU from hepatitis C care among patients with TB Among patients with TB who had a positive HCV antibody test result, LTFU before viremia testing was more common among females, patients enrolled in second-line TB treatment for drug-resistant TB, and patients with unsuccessful or unknown TB treatment outcomes (Table 5). LTFU at this step of the hepatitis C care cascade substantially decreased among patients diagnosed with TB in recent years, from 32% among patients diagnosed with TB in 2017 to 12% among patients diagnosed in 2019. In multivariable analysis, MDR TB was associated with increased risk of LTFU after positive HCV antibody test (aRR = 1.41, 95% CI [1.12, 1.76], p = 0.003). Download: PPT PowerPoint slide PNG larger image TIFF original image Table 5. LTFU from HCV care before HCV viremia testing among patients with TB. https://doi.org/10.1371/journal.pmed.1004121.t005 Among patients with DS TB who had confirmed viremic HCV infection, LTFU before hepatitis C treatment initiation was more common among females, and patients with successful TB treatment outcomes (Table 6). The proportion of patients LTFU before hepatitis C treatment initiation increased by year, from 21% among patients diagnosed with TB in 2015 to 56% among patients diagnosed with TB in 2019. In multivariable analysis, we could not identify any factors meaningfully associated with LTFU before HCV treatment initiation. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 6. LTFU from hepatitis C care before treatment initiation among patients with drug-susceptible TB. https://doi.org/10.1371/journal.pmed.1004121.t006 Participants A total of 14,993 records were obtained from the NTP database; 913 records (6.1%) were excluded due to a missing national ID number (Fig 2). The remaining 14,080 records corresponded to 12,767 individual patients, of whom 11,985 (94%) were aged ≥18 years and included in the analysis. Among those included (n = 11,985), 70% were males, 82% were unemployed at the time of TB diagnosis, 6.9% reported a history of incarceration, 11% had DR TB, 91% had 1 episode of TB treatment during the 2015 to 2019 period (i.e., diagnosed with TB only once), 7.9% had 2 episodes of TB treatment, and 1.6% had 3 to 8 episodes of TB treatment (Table 1). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. Flow chart of patients with TB and their HCV screening status—2015–2019. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g002 The hepatitis C screening registry contained data on 1,849,820 adults tested for HCV antibodies between January 1, 2015 and September 30, 2020, who were not diagnosed with TB in that period. The median age in this population was 46 years (interquartile range [IQR]: 31 to 62), and 54.8% were male (Tables 2 and 3). Hepatitis C care cascade among persons with and without TB After linking HCV testing information from the hepatitis C databases, 9,341 (78%) of 11,985 persons with TB were found to have been tested for HCV antibodies between January 1, 2015 and September 30, 2020 (Table 4). The proportion of patients with TB tested for HCV antibodies increased per year. Among patients diagnosed with TB in 2015, 60% were tested for HCV antibodies sometime during the study period. This proportion increased each subsequent year and reached 90% among patients diagnosed with TB in 2019. In the same period, HCV antibody-positivity in this population decreased from 25% among people diagnosed with TB in 2015 to 14% among those diagnosed in 2019. Among persons with TB included in our study, 276 (3.0%) had already started treatment for HCV infection before their first TB diagnosis and were included in the care cascade analysis of patients without TB. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 4. HCV testing by year of first TB diagnosis among adult patients diagnosed with TB in 2015–2019. https://doi.org/10.1371/journal.pmed.1004121.t004 Among the total of 9,065 patients with TB who were tested for hepatitis C, 1,665 (18%) patients (Fig 2) had a positive antibody test, of whom 108 (6%) died without undergoing viremia testing. Of the remaining 1,557 patients still alive, 316 (20%) were LTFU (i.e., never underwent viremia testing) (Figs 3A and S5); among 1,241 patients who underwent viremia results, viremic HCV infection was found in 1,025 (83%) patients, of whom 443 (43%) did not start hepatitis C treatment, including 80 persons that died (8% out of all viremic cases, 18% of those who did not start hepatitis C treatment). In a subset of patients with DR TB, the proportion of patients who were LTFU was even higher—27% and 52% after positive antibody and viremia tests, respectively (S6 Fig). For comparison, in the hepatitis C care cascade among patients without TB (Fig 3B), the proportion of patients who were LTFU from hepatitis C care after positive antibody test or after positive viremia test were 14% and 19%, respectively. Mortality was also lower among patients without TB: 3,842 (3%) of HCV antibody-positive persons died without viremia testing, and 3,127 (4%) viremic patients died without received hepatitis C treatment. Among patients initiating HCV treatment, those with TB were also less likely to complete hepatitis C treatment (89% versus 94%) and SVR assessment (59% versus 76%) than those without TB. However, among those who completed hepatitis C treatment and were tested for SVR, the cure rate was comparable to those without TB (98.3% versus 98.9%) (Figs 3A and 3B and S7 and S1 Table). Overall, among patients with confirmed viremic HCV infection, patients with TB were less likely to complete hepatitis C treatment (51%) and have a documented SVR achieved (28%), compared to patients without TB (76% and 55%, respectively; Fig 3C). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Hepatitis C care cascade among patients with and without tuberculosis, Georgia, 2015–2020. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. Note: Red lines represent the percent change between 2 consecutive steps in the care cascade, i.e., adjacent bars of the charts. https://doi.org/10.1371/journal.pmed.1004121.g003 Timeliness of HCV viremia testing and treatment initiation Among those who had information on follow-up duration available, we used survival analysis methods to compare time from positive antibody test to viremia testing or censoring between patients with TB (n = 1,499) and patients without TB (n = 118,396). After positive HCV antibody test, 3,840 (3.2%) persons died, and patients without TB had viremia testing sooner than patients with TB (HR = 1.46, 95% CI [1.39, 1.54], p < 0.001) (Figs 4 and S8). Among those who underwent viremia testing, the median time from a positive antibody test result to viremia testing was 17 days (IQR: 3 to 248 days) among patients with TB (n = 1,116) and 6 days (IQR: 1 to 24) among patients without TB (n = 97,554) (p < 0.001). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Comparison of HCV viremia testing cumulative incidence among HCV seropositive patients with and without TBa,b: January 1, 2015–September 30, 2020. a Blue line corresponds to the patients with confirmed HCV infection without evidence of TB diagnosis in 2015–2019. Red line corresponds to patients with confirmed HCV infection diagnosed with TB in 2015–2019. b Gray’s test p-value <0.001. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g004 We also compared time from positive viremia test to hepatitis C treatment initiation or censoring among 1,021 patients with TB and 87,964 patients without TB who had active viremia and information on follow-up duration available. After a positive viremia test, patients without TB started hepatitis C treatment sooner than patients with TB (HR = 2.05, 95% CI [1.87, 2.25], p < 0.001) (Figs 5 and S9). Among those who started hepatitis C treatment, the median time from confirmation to treatment initiation was 196 days (IQR: 76 to 523) among patients with MDR TB (n = 78), 99 days (IQR: 54 to 263) among patients with DS TB (n = 460), and 86 days (IQR: 56 to 203) among patients without TB (n = 71,128). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 5. Comparison of cumulative incidence of HCV treatment initiation among patients with confirmed HCV infection with and without TBa,b: January 1, 2015–September 30, 2020. a Blue line corresponds to patients with confirmed HCV infection without evidence of TB diagnosis in 2015–2019. Red line corresponds to patients with confirmed HCV infection diagnosed with TB in 2015–2019. b Gray’s test p-value <0.001. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.g005 Factors associated with LTFU from hepatitis C care among patients with TB Among patients with TB who had a positive HCV antibody test result, LTFU before viremia testing was more common among females, patients enrolled in second-line TB treatment for drug-resistant TB, and patients with unsuccessful or unknown TB treatment outcomes (Table 5). LTFU at this step of the hepatitis C care cascade substantially decreased among patients diagnosed with TB in recent years, from 32% among patients diagnosed with TB in 2017 to 12% among patients diagnosed in 2019. In multivariable analysis, MDR TB was associated with increased risk of LTFU after positive HCV antibody test (aRR = 1.41, 95% CI [1.12, 1.76], p = 0.003). Download: PPT PowerPoint slide PNG larger image TIFF original image Table 5. LTFU from HCV care before HCV viremia testing among patients with TB. https://doi.org/10.1371/journal.pmed.1004121.t005 Among patients with DS TB who had confirmed viremic HCV infection, LTFU before hepatitis C treatment initiation was more common among females, and patients with successful TB treatment outcomes (Table 6). The proportion of patients LTFU before hepatitis C treatment initiation increased by year, from 21% among patients diagnosed with TB in 2015 to 56% among patients diagnosed with TB in 2019. In multivariable analysis, we could not identify any factors meaningfully associated with LTFU before HCV treatment initiation. Download: PPT PowerPoint slide PNG larger image TIFF original image Table 6. LTFU from hepatitis C care before treatment initiation among patients with drug-susceptible TB. https://doi.org/10.1371/journal.pmed.1004121.t006 Discussion In this study of 2 large-scale public health programs in the country of Georgia, we found that LTFU from the hepatitis C care cascade was high overall, but substantially more common among persons with TB than those without TB. Persons diagnosed with TB between 2015 and 2019 were found to have a high prevalence (18%) of HCV antibodies, and approximately 20% of those with HCV antibodies did not receive viremia testing to evaluate for chronic HCV infection. Furthermore, almost 2 out of 5 patients with TB who had confirmed viremic HCV infection had not started hepatitis C treatment within the National Hepatitis C Elimination program in Georgia even after they completed TB treatment. We also found that proportion of patients with TB who get HCV antibody testing and subsequent viremia testing increased in recent years, while proportion of those with viremia who receive treatment decreased. Our findings highlight the importance of improving integration and linkage to hepatitis C diagnostic and treatment services among patients with TB. LTFU at different stages of hepatitis C care is a major barrier in controlling HCV infection globally [17]. Timely treatment of hepatitis C reduces the risk of cirrhosis, hepatocellular carcinoma, and death [37,38]. Therefore, LTFU from hepatitis C care can increase the risk of adverse health outcomes. Additionally, untreated hepatitis C contributes to further transmission of the infection, hindering progress towards achieving elimination goals globally [17]. Georgia faces the same problem despite the nationwide Hepatitis C Elimination Program, which offers free diagnostic and treatment services. Previous studies from Georgia have found that a high proportion of patients are LTFU before viremia testing and before hepatitis C treatment initiation [26,28]. In the first few years of the program, some patients had to co-pay for the viremia testing and other pretreatment diagnostic procedures, which was one of the main barriers for enrollment in the treatment [39,40]. Since then, this barrier has been removed by eliminating any copayment obligation for patients [28], but the LTFU still remains a challenge. Recent efforts that aimed to re-engage the persons who are LTFU in the care by actively contacting them achieved some success [41], but the gap is still wide. Similar interventions to promote linkage to care among specific target groups, such as patients with TB, have not been put in place. Our findings demonstrate LTFU from hepatitis C care is even more pronounced among patients with active TB disease coinfected with HCV. Notably, among patients with TB, LTFU before HCV viremia testing has decreased in recent years, which could be explained by a change in policy. According to Georgian 2018 TB management guideline, if an HCV antibody test is positive, blood samples are taken from the patient and sent for HCV viremia testing. Therefore, this policy change, coupled with complete financial coverage of viremia testing by the government [28], might have removed geographic, logistical, and financial barriers to viremia testing for patients with TB. This improvement highlights that systemic changes such as financial and logistic support could be key to improving retention of patients in the hepatitis C care. Furthermore, LTFU before hepatitis C treatment initiation is also very common among patients with TB. Higher LTFU before hepatitis C treatment initiation among patients with TB compared to those without TB could be explained by the fact that TB treatment is long (from 6 to 24 months) and can be associated with severe adverse reactions [42]. Additionally, HCV treatment is usually not initiated until TB treatment completion due to drug–drug interaction concerns [20]. Therefore, in the absence of a formal linkage program, patients with TB might be more reluctant to start another long treatment course for hepatitis C due to treatment fatigue and previous negative experiences [43]. This issue has not been explicitly studied among patients with both TB and HCV. Still, fatigue from diagnostic and treatment procedures as a risk factor for treatment discontinuation has been described among patients with TB and those coinfected with human immunodeficiency virus (HIV) [44–47]. This is even more likely in the case of HCV infection without advanced liver damage because patients might not experience any symptoms associated with HCV infection and might not feel the urgency to seek hepatitis C care [48]. We also found that the proportion of patients LTFU before hepatitis C treatment initiation is increasing by year among patients who completed treatment for DS TB. Several reasons could explain this finding: (1) Patients with TB diagnosed in recent years had less time to get enrolled in HCV elimination program after completing their TB treatment; (2) the COVID-19 pandemic and related restrictions in 2020 could have affected patients’ ability and willingness to start treatment for hepatitis C; and (3) the change in policy in 2018, as described above, removed a barrier to viremia testing for many patients, increasing the denominator of patients with TB diagnosed with HCV infection and eligible for treatment. This policy change that increased the viremia testing rates was not accompanied with a similar intervention to increase treatment initiation rates. Therefore, we can speculate that the number of patients who needed treatment (i.e., those with positive viremia result) increased at a higher rate than the number of patients who received treatment, reflected in increasing proportion of LTFU before treatment initiation in recent years. Additional surveys should be conducted among patients with TB to identify more specific reasons and challenges for not initiating hepatitis C treatment after a positive viremia test. Overall, hepatitis C care cascade analyses in other countries show a substantial heterogeneity in terms of viremia testing and treatment uptake. However, LTFU after positive antibody or viremia test is a universal challenge that many countries face, with Georgia usually having higher than average rates. For example, 1 review found that among countries reporting the nationwide data, Georgia had the highest rate of viremia testing among HCV antibody positive individuals [49]. Treatment uptake among HCV viremia–positive individuals ranged even more widely, from 5% to 95%. Georgia with 80% overall treatment uptake was ranked higher on this spectrum. Due to these similarities of challenges in hepatitis C care cascade across countries, results from this study have potential implications for other countries with national TB programs. In 2015, Georgia became the first country to formally initiate a nationwide hepatitis C elimination program, with significant progress achieved in scaling up HCV screening, diagnostic and treatment services throughout the country [26]. However, previously published literature and our analysis highlight that LTFU from different stages of hepatitis C care poses a major challenge to HCV elimination goals, and patients with TB are even more likely to be LTFU from hepatitis C care. To address this issue, program partners should consider introducing additional activities directed at more integrated care of hepatitis C and TB and interventions targeted at retaining patients in hepatitis C care. These interventions could include patient navigators that have been found effective in hepatitis C linkage to care [50] or providing hepatitis C treatment in the same facilities where they receive TB treatment. Additionally, patients with MDR TB could be treated for HCV infection and TB concomitantly. There is no documented drug–drug interaction between DAAs used in HCV treatment and second-line TB drugs [20], and concomitant treatment of TB and HCV infection have been successfully implemented in small samples of patients [51,52]. For patients with DS TB, replacing rifampin with rifabutin may allow for simultaneous TB and HCV treatment [20]. Therefore, the existing infrastructure of NTP can help the HCV elimination efforts in Georgia by achieving microelimination of HCV infection among patients with TB—a strategy that is proposed as a useful approach for achieving the overall elimination goals [53]. In addition to tuberculosis, there are several other conditions that have substantial overlap with hepatitis C, such as HIV infection, injection drug use, and incarceration. Georgia already has a good experience of integrating services for some of these different conditions. For example, Georgia is considered as an example of good practice in terms of integrating the diagnostic services for HIV, tuberculosis, and hepatitis C in primary healthcare [54]. Hepatitis C testing and treatment in corrections system started even before the launch of the nationwide hepatitis C elimination program [55]. Furthermore, hepatitis C treatment for people living with HIV is usually offered in the same facilities where they receive antiretroviral treatment. However, hepatitis C viremia testing and treatment is not systematically integrated in the harm reduction sites, even though it has been found feasible and preferable for people who inject drugs to receive these services in the same sites [56,57]. Our study has several limitations. First, due to missing national ID numbers, we had to exclude 6% of observations from the NTP database. Second, the limited number of variables in the hepatitis C screening registry and differences in variables available in hepatitis C and TB databases did not allow us to conduct a more in-depth analysis to adjust for potential confounders. For that reason, our time-to-event analyses are limited to the crude comparison of cumulative incidence curves and crude hazards ratios between patients with and without TB, rather than trying to explore any causal associations that would require confounding adjustment. However, even without exploring the causal association, the pronounced difference between patients with both TB and hepatitis C and those with only hepatitis C warrants the attention to reduce the gap in linkage to hepatitis C care and treatment among patients with TB. Third, we excluded patients with DR TB from the risk factor analysis of LTFU before HCV treatment initiation due to high heterogeneity in terms of when individual patients with DR TB become eligible for hepatitis C treatment. In conclusion, we found that LTFU from hepatitis C care after positive HCV antibody and viremia testing is more common among patients with TB than those without TB. Existing large-scale public health programs for both TB and hepatitis C in Georgia create a unique opportunity for integrated care of these 2 infectious diseases, which could potentially reduce LTFU. Though our study was not designed to identify effective interventions for reducing LTFU, integrated care should utilize a patient-centered approach and include the following suggested interventions: (1) merging care for TB and HCV infection for coinfected patients; (2) parallel treatment of TB and hepatitis C using drug regimens without substantial drug–drug interactions; (3) active support for patients with TB to navigate the referral process and enroll in the hepatitis C treatment program, whenever it is not feasible to receive HCV care at the same facility where they receive TB treatment. These interventions can positively impact both TB and hepatitis C programs by providing timely hepatitis C care to patients with TB, thereby reducing LTFU and improving overall patient health outcomes. Supporting information S1 Checklist. STROBE Statement—Checklist of items that should be included in reports of cohort studies. https://doi.org/10.1371/journal.pmed.1004121.s001 (DOCX) S2 Checklist. Inclusivity in global research. https://doi.org/10.1371/journal.pmed.1004121.s002 (DOCX) S1 Analysis Plan. Analysis plan. https://doi.org/10.1371/journal.pmed.1004121.s003 (DOCX) S1 Fig. Directed acyclic graph depicting the causal relations between exposure of interest (Employment), outcome of interest (loss to follow-up), and other covariates. IDP, internally displaced person; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s004 (DOCX) S2 Fig. Directed acyclic graph depicting the causal relations between exposure of interest (MDR TB), outcome of interest (loss to follow-up), and other covariates. MDR, multidrug-resistant; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s005 (DOCX) S3 Fig. Directed acyclic graph depicting the causal relations between exposure of interest (previously treated TB), outcome of interest (loss to follow-up), and other covariates. MDR, multidrug-resistant; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s006 (DOCX) S4 Fig. Directed acyclic graph depicting the causal relations between exposure of interest (TB treatment outcome), outcome of interest (loss to follow-up), and other covariates. MDR, multidrug-resistant; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s007 (DOCX) S5 Fig. Care cascade of hepatitis C among patients with tuberculosis, starting from the total number of patients with TB. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. https://doi.org/10.1371/journal.pmed.1004121.s008 (DOCX) S6 Fig. Hepatitis C care cascade among HCV seropositive patients with TB, stratified by their drug-resistance status. Note: Red lines represent the percent change between 2 consecutive steps in the care cascade, i.e., adjacent bars of the charts. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. https://doi.org/10.1371/journal.pmed.1004121.s009 (DOCX) S7 Fig. Comparison of hepatitis C care cascade among patients with and without TB in Georgia, 2015–2020. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. https://doi.org/10.1371/journal.pmed.1004121.s010 (DOCX) S8 Fig. Kaplan–Meier curves of time from positive antibody test to viremia testing, with 95% confidence bands. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s011 (DOCX) S9 Fig. Kaplan–Meier curves of time from positive viremia test or TB treatment completion (whichever occurred last) to hepatitis C treatment initiation, with 95% confidence bands. HCV, hepatitis C virus; TB, tuberculosis. https://doi.org/10.1371/journal.pmed.1004121.s012 (DOCX) S1 Table. Comparison of proportions at each step of care cascade between patients with and without TB, Georgia, 2015–2020. HCV, hepatitis C virus; SVR, sustained virologic response; TB, tuberculosis; Tx, treatment. https://doi.org/10.1371/journal.pmed.1004121.s013 (DOCX) Acknowledgments The authors would like to thank Amiran Gamkrelidze and Zaza Avaliani for their support in obtaining the data used in this study. Disclaimer The findings and conclusions in this report are those of the authors and do not necessarily represent the official position of the US Centers for Disease Control and Prevention. Disclaimer The findings and conclusions in this report are those of the authors and do not necessarily represent the official position of the US Centers for Disease Control and Prevention.
Global, regional, and national estimates of the impact of a maternal Klebsiella pneumoniae vaccine: A Bayesian modeling analysisKumar, Chirag K.;Sands, Kirsty;Walsh, Timothy R.;O’Brien, Seamus;Sharland, Mike;Lewnard, Joseph A.;Hu, Hao;Srikantiah, Padmini;Laxminarayan, Ramanan
doi: 10.1371/journal.pmed.1004239pmid: 37216371
Background Despite significant global progress in reducing neonatal mortality, bacterial sepsis remains a major cause of neonatal deaths. Klebsiella pneumoniae (K. pneumoniae) is the leading pathogen globally underlying cases of neonatal sepsis and is frequently resistant to antibiotic treatment regimens recommended by the World Health Organization (WHO), including first-line therapy with ampicillin and gentamicin, second-line therapy with amikacin and ceftazidime, and meropenem. Maternal vaccination to prevent neonatal infection could reduce the burden of K. pneumoniae neonatal sepsis in low- and middle-income countries (LMICs), but the potential impact of vaccination remains poorly quantified. We estimated the potential impact of such vaccination on cases and deaths of K. pneumoniae neonatal sepsis and project the global effects of routine immunization of pregnant women with the K. pneumoniae vaccine as antimicrobial resistance (AMR) increases. Methods and findings We developed a Bayesian mixture-modeling framework to estimate the effects of a hypothetical K. pneumoniae maternal vaccine with 70% efficacy administered with coverage equivalent to that of the maternal tetanus vaccine on neonatal sepsis infections and mortality. To parameterize our model, we used data from 3 global studies of neonatal sepsis and/or mortality—with 2,330 neonates who died with sepsis surveilled from 2016 to 2020 undertaken in 18 mainly LMICs across all WHO regions (Ethiopia, Kenya, Mali, Mozambique, Nigeria, Rwanda, Sierra Leone, South Africa, Uganda, Brazil, Italy, Greece, Pakistan, Bangladesh, India, Thailand, China, and Vietnam). Within these studies, 26.95% of fatal neonatal sepsis cases were culture-positive for K. pneumoniae. We analyzed 9,070 K. pneumoniae genomes from human isolates gathered globally from 2001 to 2020 to quantify the temporal rate of acquisition of AMR genes in K. pneumoniae isolates to predict the future number of drug-resistant cases and deaths that could be averted by vaccination. Resistance rates to carbapenems are increasing most rapidly and 22.43% [95th percentile Bayesian credible interval (CrI): 5.24 to 41.42] of neonatal sepsis deaths are caused by meropenem-resistant K. pneumoniae. Globally, we estimate that maternal vaccination could avert 80,258 [CrI: 18,084 to 189,040] neonatal deaths and 399,015 [CrI: 334,523 to 485,442] neonatal sepsis cases yearly worldwide, accounting for more than 3.40% [CrI: 0.75 to 8.01] of all neonatal deaths. The largest relative benefits are in Africa (Sierra Leone, Mali, Niger) and South-East Asia (Bangladesh) where vaccination could avert over 6% of all neonatal deaths. Nevertheless, our modeling only considers country-level trends in K. pneumoniae neonatal sepsis deaths and is unable to consider within-country variability in bacterial prevalence that may impact the projected burden of sepsis. Conclusions A K. pneumoniae maternal vaccine could have widespread, sustained global benefits as AMR in K. pneumoniae continues to increase. Why was this study done? Approximately 1 million newborns yearly die within the first 4 weeks of life due to bacteria infecting their bloodstream with Klebsiella pneumoniae (K. pneumoniae) as the leading cause of such infections. There have been numerous recent advancements in developing viable K. pneumoniae vaccine and antibody-based treatments in preclinical models, with some treatments reaching Phase 1 clinical trials. The impacts of vaccination must be quantified and may prove useful to prioritizing vaccine distribution and better understanding the burden of sepsis. What did the researchers do and find? Using a Bayesian mixture-model based on clinical surveillance of neonatal sepsis, we present country-specific estimates for the number of deaths and cases of antimicrobial-resistant neonatal sepsis caused by K. pneumoniae that would be averted if a vaccine with 70% efficacy were given to pregnant mothers. We find that most cases of K. pneumoniae neonatal sepsis are resistant to first-line treatments, such as ampicillin and gentamicin. We estimate that a vaccine with 70% efficacy could prevent 399,015 [95th percentile credible interval (CrI): 334,523 to 485,442] cases and 80,258 [CrI: 18,084 to 189,040] neonatal deaths. What do these findings mean? A maternal vaccine that confers newborns with protection from K. pneumoniae infection could reduce neonatal sepsis deaths in many low- and middle-income countries (LMICs) by nearly 15%. This would help to achieve targets set by the World Health Organization (WHO) for improved child health globally and to mitigate inequities in neonatal survival in low- and middle-income settings compared to high-income settings. Reducing cases of neonatal sepsis by vaccination could also contribute to reduced antibiotic use, subsequent improvements in antimicrobial resistance (AMR) rates, and a reduction in healthcare utilization and expenditure. Introduction Each year, an estimated 800,000 newborns die within the first 4 weeks of life due to sepsis [1–3]. Several studies, including most recently the Child Health and Mortality Prevention Surveillance (CHAMPS) [4] and the Burden of Antibiotic Resistance in Neonates from Developing Societies (BARNARDS) [5], indicate that Klebsiella pneumoniae is the leading cause of neonatal sepsis (Fig A in S1 Text). Increasing multidrug resistance (MDR) and lack of access to appropriate antibiotics contribute to high levels of morbidity and mortality among K. pneumoniae sepsis cases [6–8]. Though improved sanitation and infection control could limit K. pneumoniae transmission, they can be challenging to implement in many resource-poor settings [9]. A maternal vaccine against K. pneumoniae that confers protection to newborns via transplacental antibody transfer could reduce the burden of neonatal sepsis, but the extent of corresponding reductions in sepsis remains unknown. Primary approaches to developing a K. pneumoniae vaccine have targeted the K capsular polysaccharide antigens or the O lipopolysaccharide antigens—both of which demonstrate significant diversity across K. pneumoniae strains—or an inactivated whole cell [10]. Animal data from an early polysaccharide vaccine construct demonstrated a robust immune response and subsequent neutralization of K. pneumoniae bacterial isolates in mouse models when administered intravenously [11]. A 24-valent Klebsiella capsular polysaccharide vaccine tested in humans demonstrated high anti-polysaccharide immunoglobin titers [12], but a clinical trial to validate this vaccine was inconclusive because of a limited supply of materials [13]. More recently, there is growing interest in conjugate polysaccharide vaccines targeting the K and/or O antigens [14–17] with the first such vaccine entering a clinical phase [18,19]. Combined, these data suggest that a maternal vaccine to prevent K. pneumoniae demands further attention and investigation. As with previous [20] and future [21] proposed maternal vaccination programs, this vaccine would be administered to pregnant mothers, triggering an immune response that protects the mother and, critically, transfers maternal antibodies to the developing baby, providing protection to the neonate (or young infant) through the period of greatest risk [22]. Vaccination could reduce the overall burden of sepsis and decrease antibiotic use (fewer K. pneumoniae sepsis cases) [23] and K. pneumoniae nosocomial outbreaks (vaccination generates community-level immunity) [24], which are common in LMICs [5]. Additionally, by averting infection and reducing the need for treatment, vaccination may prevent lengthy hospital stays for newborn sepsis, thereby reducing the economic burden on LMIC families where healthcare costs are often deferred to the patient [7]. In this paper, we estimated the global, regional, and national impact of a potential maternal vaccine against K. pneumoniae deployed at the current levels of coverage of the maternal tetanus vaccine. We also examined the antimicrobial resistance (AMR) of sepsis-causing isolates. Finally, we estimated the continued benefits of such a vaccine in 2030 under future scenarios of increasing AMR. Methods We used data from 3 global clinical data sets of neonatal sepsis and a Bayesian modeling framework to estimate the annual neonatal infections and deaths by country caused by K. pneumoniae (although the methodology is applicable to any other sepsis-causing bacterial etiology; Figs O to U in S1 Text). We used target vaccine efficacy and coverage to project the total number of infections and deaths that could be averted using a maternal vaccine. We complemented these analyses with a mathematical model that estimated the fraction of infections and deaths that are resistant to the antibiotics recommended by World Health Organization (WHO) treatment guidelines for neonatal sepsis (ampicillin and gentamicin as a first-line treatment and amikacin and ceftazidime as a second-line treatment) [25] and those most used to treat neonatal sepsis (adding meropenem to the above list) [7]. Data Three data sets were used to quantify K. pneumoniae neonatal sepsis burden: the Child Health and Mortality Prevention Surveillance (CHAMPS) [4], the Burden of Antibiotic Resistance in Neonates from Developing Societies (BARNARDS) [5], and the Global Neonatal Sepsis Observational Study (NeoObs) [26]. CHAMPS, initiated in 2016, conducts minimally invasive tissue sampling (MITS) on deceased children under 5 to determine causal drivers of death at 8 study sites in 7 low- and middle-income countries (LMICs) across Africa and South-East Asia (Bangladesh, 2 sites in Ethiopia, Kenya, Mali, Mozambique, Sierra Leone, and South Africa). In contrast, BARNARDS and NeoObs are birth-cohort studies that tracked neonates prospectively through their first month of life, evaluating incidence of sepsis, the antibiotic susceptibility of associated isolates, and resulting mortality. BARNARDS surveilled neonates across 12 study sites in 7 countries in Africa and South-East Asia (2 sites in Bangladesh, India, 2 in Pakistan, Ethiopia, 3 in Nigeria, 2 in Rwanda, and South Africa) from 2015 to 2017. Finally, NeoObs reported neonatal sepsis cases, deaths, and antibiotic susceptibility of the isolated bacteria in 19 sites across 11 countries in Europe (Italy, Greece), Latin America (2 sites in Brazil), Africa (Kenya, Uganda, 3 sites in South Africa), South-East Asia (Bangladesh, 3 sites in India, 2 sites in Thailand), and the Western Pacific (3 sites in China, Vietnam) from 2018 to 2020. Table A in S1 Text provides the number of neonates surveilled over the study period by country. We restricted all analyses to cases within each study with culture-confirmed K. pneumoniae sepsis in neonates based on isolates obtained from blood or cerebrospinal fluids (CSF) (Fig B in S1 Text). Modeling techniques We used a Bayesian mixture-modeling approach to aggregate data across the 3 studies and produce consistent and robust estimates of neonatal sepsis infections, deaths, and drug resistance distributions. We estimated ps,l, the probability that K. pneumoniae was the cause of death for a neonate who died from neonatal sepsis, in each location, l, and for each study, s. We assumed a noninformative uniform prior on ps,l at all locations. We modeled whether K. pneumoniae was present in a clinical sample from a neonate who died of culture-confirmed sepsis as a binomial random variable and analytically solved for the posterior distribution of ps,l (Section 1 in S2 Text): where Nkps, l is the number of neonates who were culture-positive for K. pneumoniae and died, Nds, l is the number of neonates who had culture-confirmed sepsis (from any etiology) and died, both quantities at each location, l, in each study, s. β indicates the standard univariate beta distribution. In countries where multiple studies conducted surveillance (ex., both CHAMPS and BARNARDS had study sites in the country), we obtained multiple independent estimates of the distribution of ps,l, one for each study that had a site in the country. We used a mixture-modeling approach to aggregate these distributions so that each independent distribution for ps,l was weighted proportional to the number of neonates sampled in the study at that location: We used a mixture-modeling approach rather than aggregating the samples into 1 beta distribution to account for subpopulation variability across sites within the same country from different studies. To efficiently characterize the resulting distribution, we drew 10,000 samples using Latin hypercube sampling (LHS) [27]. We chose LHS over another distribution sampling technique such as Monte Carlo because LHS converges to the true distribution faster. We extended this analysis to estimate the absolute number of deaths attributable to K. pneumoniae using data from the United Nations Child Health Epidemiology Reference Group (CHERG) [1] and the University of Washington’s Institute for Health Metrics and Evaluation Global Burden of Disease (GBD) [28,29] on the total number of neonatal deaths and the total number of neonatal sepsis deaths. Estimates of pl were mapped into estimates of the number of deaths and the percentage of all neonatal deaths attributable to K. pneumoniae by scaling the estimates of pl by the estimates with uncertainty from CHERG and GBD for the total number of neonatal sepsis deaths. Unless otherwise mentioned, we report values derived using CHERG data. We determined the number of vaccine-avertable deaths associated with the deployment of a maternal K. pneumoniae vaccine at an effective coverage level equal to that of the maternal tetanus vaccine (median: 90%; range: 38.5% to 100% of pregnant women immunized) [30]. For modeling purposes, we assumed that the efficacy of the vaccine would be 70%: This estimate is based on a conjugate vaccine candidate that targets the 15 most common K. pneumoniae capsular serotypes that cause invasive infections in neonates. Current sero-epidemiology suggests that these serotypes account for approximately 70% of neonatal sepsis cases [14]. To determine the avertable cases from avertable deaths, we estimated the case fatality ratio (CFR) for K. pneumoniae sepsis using data from BARNARDS. We used the same Bayesian beta framework with a uniform prior described above (we treated whether the neonate died from K. pneumoniae or recovered as the binomial event for which we seek to estimate the probability parameter) to estimate the CFR. To correctly propagate uncertainty across combination products (i.e., when trying to determine the number of cases from the number of deaths using the CFR), we individually multiplied independent samples from the distributions for ps,l and the CFR. We summarized the resulting distribution using the median and 95% credible intervals (CrIs). We use this approach throughout for propagating uncertainty. Finally, we estimated antibiotic resistance profiles using data from BARNARDS and NeoObs at the WHO regional level (Section 2 in S2 Text). We used the Bayesian beta and mixture-modeling framework with a noninformative prior described earlier to estimate the probability that an isolate was resistant to a given antibiotic. We extrapolated our estimates from the selected countries for which we had neonatal sepsis surveillance data global estimates for all countries of avertable neonatal sepsis deaths and cases by fitting regression equations. We related the total number of neonatal deaths due to sepsis (for which estimates are available from CHERG) to the estimated number of avertable neonatal deaths if a K. pneumoniae vaccine were distributed as described above using standard ordinary least squares. We also calculated the number of avertable antibiotic-resistant neonatal deaths by country using our estimates for antibiotic resistance at the WHO regional level. To estimate the future benefits of such a vaccine, we analyzed the temporal trends in the presence of AMR genes in K. pneumoniae by using isolate genomes available on PathogenWatch [31] and in the BARNARDS data set. We limited our analysis to human isolates from clinical sources that were not taken as rectal or carriage samples (S1 Table and Fig C in S1 Text). We gathered metadata on the isolates, including their date, source, and country. In total, we analyzed 9,070 K. pneumoniae genomes from 68 countries from 2001 to 2020 (Figs D and E in S1 Text). To identify antimicrobial resistance genes (ARGs) in each genome, we used ABRicate (v1.1.0) with the reference database ResFinder (last updated on June 18, 2022). We fit linear probability model equations to the temporal trends in the number of aminoglycoside-resistant genes (AGRGs) and carbapenem-resistant genes (CRGs) and extrapolated those trends to 2030 (i.e., to predict the future number of AGRGs and CRGs). We used the fit linear probability model to determine the projected increase in ARGs and resistance by 2030 and used that to estimate the increase in the absolute number of aminoglycoside- and carbapenem-resistant neonatal sepsis cases and deaths. Fig F in S1 Text shows a schematic diagram of the model, and Table B in S1 Text denotes all model parameters. Data Three data sets were used to quantify K. pneumoniae neonatal sepsis burden: the Child Health and Mortality Prevention Surveillance (CHAMPS) [4], the Burden of Antibiotic Resistance in Neonates from Developing Societies (BARNARDS) [5], and the Global Neonatal Sepsis Observational Study (NeoObs) [26]. CHAMPS, initiated in 2016, conducts minimally invasive tissue sampling (MITS) on deceased children under 5 to determine causal drivers of death at 8 study sites in 7 low- and middle-income countries (LMICs) across Africa and South-East Asia (Bangladesh, 2 sites in Ethiopia, Kenya, Mali, Mozambique, Sierra Leone, and South Africa). In contrast, BARNARDS and NeoObs are birth-cohort studies that tracked neonates prospectively through their first month of life, evaluating incidence of sepsis, the antibiotic susceptibility of associated isolates, and resulting mortality. BARNARDS surveilled neonates across 12 study sites in 7 countries in Africa and South-East Asia (2 sites in Bangladesh, India, 2 in Pakistan, Ethiopia, 3 in Nigeria, 2 in Rwanda, and South Africa) from 2015 to 2017. Finally, NeoObs reported neonatal sepsis cases, deaths, and antibiotic susceptibility of the isolated bacteria in 19 sites across 11 countries in Europe (Italy, Greece), Latin America (2 sites in Brazil), Africa (Kenya, Uganda, 3 sites in South Africa), South-East Asia (Bangladesh, 3 sites in India, 2 sites in Thailand), and the Western Pacific (3 sites in China, Vietnam) from 2018 to 2020. Table A in S1 Text provides the number of neonates surveilled over the study period by country. We restricted all analyses to cases within each study with culture-confirmed K. pneumoniae sepsis in neonates based on isolates obtained from blood or cerebrospinal fluids (CSF) (Fig B in S1 Text). Modeling techniques We used a Bayesian mixture-modeling approach to aggregate data across the 3 studies and produce consistent and robust estimates of neonatal sepsis infections, deaths, and drug resistance distributions. We estimated ps,l, the probability that K. pneumoniae was the cause of death for a neonate who died from neonatal sepsis, in each location, l, and for each study, s. We assumed a noninformative uniform prior on ps,l at all locations. We modeled whether K. pneumoniae was present in a clinical sample from a neonate who died of culture-confirmed sepsis as a binomial random variable and analytically solved for the posterior distribution of ps,l (Section 1 in S2 Text): where Nkps, l is the number of neonates who were culture-positive for K. pneumoniae and died, Nds, l is the number of neonates who had culture-confirmed sepsis (from any etiology) and died, both quantities at each location, l, in each study, s. β indicates the standard univariate beta distribution. In countries where multiple studies conducted surveillance (ex., both CHAMPS and BARNARDS had study sites in the country), we obtained multiple independent estimates of the distribution of ps,l, one for each study that had a site in the country. We used a mixture-modeling approach to aggregate these distributions so that each independent distribution for ps,l was weighted proportional to the number of neonates sampled in the study at that location: We used a mixture-modeling approach rather than aggregating the samples into 1 beta distribution to account for subpopulation variability across sites within the same country from different studies. To efficiently characterize the resulting distribution, we drew 10,000 samples using Latin hypercube sampling (LHS) [27]. We chose LHS over another distribution sampling technique such as Monte Carlo because LHS converges to the true distribution faster. We extended this analysis to estimate the absolute number of deaths attributable to K. pneumoniae using data from the United Nations Child Health Epidemiology Reference Group (CHERG) [1] and the University of Washington’s Institute for Health Metrics and Evaluation Global Burden of Disease (GBD) [28,29] on the total number of neonatal deaths and the total number of neonatal sepsis deaths. Estimates of pl were mapped into estimates of the number of deaths and the percentage of all neonatal deaths attributable to K. pneumoniae by scaling the estimates of pl by the estimates with uncertainty from CHERG and GBD for the total number of neonatal sepsis deaths. Unless otherwise mentioned, we report values derived using CHERG data. We determined the number of vaccine-avertable deaths associated with the deployment of a maternal K. pneumoniae vaccine at an effective coverage level equal to that of the maternal tetanus vaccine (median: 90%; range: 38.5% to 100% of pregnant women immunized) [30]. For modeling purposes, we assumed that the efficacy of the vaccine would be 70%: This estimate is based on a conjugate vaccine candidate that targets the 15 most common K. pneumoniae capsular serotypes that cause invasive infections in neonates. Current sero-epidemiology suggests that these serotypes account for approximately 70% of neonatal sepsis cases [14]. To determine the avertable cases from avertable deaths, we estimated the case fatality ratio (CFR) for K. pneumoniae sepsis using data from BARNARDS. We used the same Bayesian beta framework with a uniform prior described above (we treated whether the neonate died from K. pneumoniae or recovered as the binomial event for which we seek to estimate the probability parameter) to estimate the CFR. To correctly propagate uncertainty across combination products (i.e., when trying to determine the number of cases from the number of deaths using the CFR), we individually multiplied independent samples from the distributions for ps,l and the CFR. We summarized the resulting distribution using the median and 95% credible intervals (CrIs). We use this approach throughout for propagating uncertainty. Finally, we estimated antibiotic resistance profiles using data from BARNARDS and NeoObs at the WHO regional level (Section 2 in S2 Text). We used the Bayesian beta and mixture-modeling framework with a noninformative prior described earlier to estimate the probability that an isolate was resistant to a given antibiotic. We extrapolated our estimates from the selected countries for which we had neonatal sepsis surveillance data global estimates for all countries of avertable neonatal sepsis deaths and cases by fitting regression equations. We related the total number of neonatal deaths due to sepsis (for which estimates are available from CHERG) to the estimated number of avertable neonatal deaths if a K. pneumoniae vaccine were distributed as described above using standard ordinary least squares. We also calculated the number of avertable antibiotic-resistant neonatal deaths by country using our estimates for antibiotic resistance at the WHO regional level. To estimate the future benefits of such a vaccine, we analyzed the temporal trends in the presence of AMR genes in K. pneumoniae by using isolate genomes available on PathogenWatch [31] and in the BARNARDS data set. We limited our analysis to human isolates from clinical sources that were not taken as rectal or carriage samples (S1 Table and Fig C in S1 Text). We gathered metadata on the isolates, including their date, source, and country. In total, we analyzed 9,070 K. pneumoniae genomes from 68 countries from 2001 to 2020 (Figs D and E in S1 Text). To identify antimicrobial resistance genes (ARGs) in each genome, we used ABRicate (v1.1.0) with the reference database ResFinder (last updated on June 18, 2022). We fit linear probability model equations to the temporal trends in the number of aminoglycoside-resistant genes (AGRGs) and carbapenem-resistant genes (CRGs) and extrapolated those trends to 2030 (i.e., to predict the future number of AGRGs and CRGs). We used the fit linear probability model to determine the projected increase in ARGs and resistance by 2030 and used that to estimate the increase in the absolute number of aminoglycoside- and carbapenem-resistant neonatal sepsis cases and deaths. Fig F in S1 Text shows a schematic diagram of the model, and Table B in S1 Text denotes all model parameters. Results Effect of a K. pneumoniae maternal vaccine on neonatal deaths We initially focused on the 18 countries for which we had surveillance data on bacterial neonatal sepsis (Fig A in S1 Text). A potential vaccine could reduce neonatal sepsis and improve newborn survival (Fig 1). Our estimates of the proportion of neonatal deaths averted with maternal K. pneumoniae vaccination are consistent whether we use CHERG or GBD data as a reference source. Under the assumption that K. pneumoniae vaccine coverage matches current effective maternal tetanus vaccination coverage levels [30], a median estimated 57,392 [95th percentile Bayesian CrI: 22,248 to 125,803] deaths, accounting for 3.83% [CrI: 1.48 to 8.40] of newborn fatalities, could be averted each year in the initial 18 countries analyzed. We find a similar benefit in reducing neonatal sepsis cases, with an estimated 286,499 [CrI: 239,259 to 346,496] cases averted yearly. The benefits of a maternal K. pneumoniae vaccine vary across countries. While the relative impact of vaccination on neonatal mortality ranges a few percentage points (Fig 1A), the estimates for the reduction in total neonatal sepsis deaths and cases span orders of magnitude across countries, reflecting the variability in the burden of neonatal deaths across countries. We project that vaccination will have the greatest relative effect on neonatal mortality in Ethiopia, Rwanda, and Sierra Leone and a smaller effect in Greece, China, and Pakistan (Fig 1A), countries where a lower fraction of neonatal sepsis deaths were associated with K. pneumoniae (Fig A in S1 Text). However, although the relative impact of vaccination on neonatal mortality (i.e., percent of deaths averted with respect to all neonatal deaths) in Pakistan is smaller compared to other countries considered, we estimate that Pakistan will have one of the greatest absolute reductions in neonatal mortality (i.e., a high number of absolute deaths/cases averted; Fig 1B and 1C), consistent with previous work that Pakistan experiences a significant K. pneumoniae burden from nosocomial outbreaks [5]. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. A maternal vaccine reduces global K. pneumoniae burden. (A) Median estimated percent of neonatal deaths that are avertable with a maternal vaccine with respect to all neonatal deaths in a given country. (B) Median estimated number of neonatal deaths that are avertable with a maternal vaccine in a given country. (C) Median estimated number of neonatal cases of sepsis that are avertable with a maternal vaccine in a given country. GBD indicates data from the Global Burden of Disease study and CHERG indicates data from the United Nations Child Health Epidemiology Reference Group. Error bars indicate 95th percentile CrIs. A pseudo log transformation is done for values between 0 and 1 in panel B. CHERG, Child Health and Epidemiology Reference Group; CrI, credible interval; GBD, Global Burden of Disease. https://doi.org/10.1371/journal.pmed.1004239.g001 Across all WHO regions, K. pneumoniae isolates are highly resistant to most antibiotics (Fig 2). We focused on those antibiotics often prescribed to treat neonatal sepsis: ampicillin, gentamicin, ceftazidime, amikacin, and meropenem (Fig G in S1 Text for all antibiotics used to treat sepsis) [7,25]. Prevalence of resistance to ampicillin is highest: 89.82% [CrI: 64.97 to 99.49] of K. pneumoniae isolates are resistant to ampicillin worldwide (Fig 2B). The average rate of resistance to gentamicin, which is usually prescribed in combination with ampicillin as the first-line treatment for sepsis per current WHO guidelines [32], is 57.22% [CrI: 31.540 to 80.42] (Fig 2B). Comparatively, fewer isolates were resistant to ceftazidime and amikacin—frequently employed as a second-line combination to treat neonatal sepsis—with 44.63% [26.84 to 72.20] and 35.80% [14.69 to 58.53], respectively, isolates identified as resistant (Fig 2B). Approximately 22.43% [5.24 to 41.42] of isolates from neonates who died were resistant to meropenem. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. AMR distribution of vaccine-avertable deaths and cases. (A, B) Median percent of isolates that are resistant to specific classes of antibiotics by WHO region. Error bars indicate 95th percentile CrIs. (A) Of those deaths that are averted, the percent of deaths that are caused by K. pneumoniae resistant to each class of antibiotics. (B) Of those cases that are averted, the percent of illnesses that are caused by K. pneumoniae resistant to each class of antibiotics. (C, D) Median percent of isolates that are resistant to multiple classes of antibiotics. Error bars indicate 95th percentile CrIs. (C) Of those deaths that are averted, the percent of isolates that have varying degrees of MDR. (D) Of those cases that are averted, the percent of isolates that have varying degrees of MDR. (E, F) Joint distribution of ampicillin and gentamicin resistance. Error bars indicate 95th percentile CrIs. (E) Of those deaths that are averted, the percent of isolates that are resistant to combinations of aminoglycoside and beta lactam drugs. (F) Of those cases that are averted, the percent of isolates that are resistant to combinations of aminoglycoside and beta lactam drugs. WHO regions with no data indicate that no isolates were gathered in those regions, so no informative modeling can be done. AMR, antimicrobial resistance; CrI, credible interval; MDR, multidrug resistance; WHO, World Health Organization. https://doi.org/10.1371/journal.pmed.1004239.g002 Because ampicillin and gentamicin are prescribed together as the first-line treatment for neonatal sepsis in many countries, we considered the joint distribution of sepsis isolates that were resistant to neither, one, or both antibiotics (Fig 2C and 2D). Resistance is high against both drugs, with an average of 86.38% [CrI: 73.09 to 94.74] of isolates resistant to both drugs across all regions (Fig 2D). Notably, rates of resistance to ampicillin alone were substantially higher than to gentamicin alone (Fig 2D). In fact, less than 5% of isolates were resistant to just gentamicin (Fig 2D). Consequently, resistance to gentamicin almost guaranteed resistance to ampicillin (Fig 2D). Different antibiotics within the same antibiotic class display varying geographic trends. For instance, resistance to ampicillin is approximately constant across all geographic regions. Comparatively, rates of resistance to ceftazidime and meropenem are variable across geographic regions: Rates of resistance to ceftazidime in Africa are high and close to those of ampicillin, approximately 60 percentage points greater than rates of resistance to meropenem. Resistance to meropenem is low (7.01% [CrI: 3.21 to 12.64]) in Africa and substantially higher in South-East Asia (66.58% [CrI: 52.68 to 78.82]), greater than resistance to ceftazidime (Fig 2B). Among aminoglycosides, resistance rates to amikacin and gentamicin are similar in all regions except Africa where gentamicin resistance is 60 percentage points greater than amikacin resistance (Fig 2B). Finally, we consider trends across antibiotics within the same treatment regimen; in particular, we compare ampicillin and gentamicin (first-line treatment) against ceftazidime and amikacin (second-line treatment). While rates of resistance to ampicillin are consistently greater than to gentamicin, rates of resistance to ceftazidime are greater than amikacin only in Africa and the Western Pacific (Fig 2B). Across all regions, we generally did not find an increase in drug-resistant K. pneumoniae isolates from neonates who died compared with all K. pneumoniae isolates in general: The prevalence of resistance to a specific drug class is approximately equal whether the isolate was taken from a septic neonate who died or not (Fig 2A and 2B). However, there are 2 key exceptions: In South-East Asia, rates of resistance to amikacin, gentamicin, and ceftazidime are higher in neonates who died, and rates of resistance to ampicillin are lower among all neonates. Nevertheless, we do not find that isolates from neonates who succumbed to sepsis were more likely to exhibit AMR (Fig 2A and 2B), resistance to ampicillin or gentamicin (Fig 2C and 2D), or higher MDR (Fig 2E and 2F) in comparison to all isolates. Most K. pneumoniae have become resistant to multiple antibiotics (Fig 2E and 2F). On average, more than 99% of isolates were resistant to at least 1 drug. Many isolates were resistant to multiple antibiotics, and isolates that were resistant to none were rare (0.15% [CrI: 2.35 × 10−9 to 18.17]). The highest MDR rates are observed in the Eastern Mediterranean and South-East Asia where over 98% of K. pneumoniae are resistant to either 4 or 5 different antibiotics (Fig 2F); however, in Africa, most isolates are resistant to 3 antibiotics only with fewer isolates resistant to 4 or 5 antibiotics (Fig 2F). Projected global and future benefits of a K. pneumoniae maternal vaccine We expanded our estimates of the effects of a maternal K. pneumoniae vaccine on childhood mortality and neonatal sepsis to all countries using regression models as described above (Table C in S1 Text, R2 of 84.0%, p < 0.001). We project that 14.91% [CrI: 5.26 to 21.24] of all neonatal sepsis deaths in countries of interest (all countries except for high-income nations)—that is, 80,258 [CrI: 18,084 to 189,040] neonatal deaths or 3.40% [CrI: 0.75 to 8.01] of all neonatal mortality—and 399,015 [CrI: 334,523 to 485,442] neonatal sepsis cases could be averted yearly once vaccination coverage reaches that of the maternal tetanus vaccine (Figs 3A and 3B, and Figs H and I in S1 Text for the 95th percentile CrIs, and S2 Table for the numeric values). We project the greatest reduction in overall neonatal mortality (i.e., fraction of deaths reduced due to vaccination) to be in Africa, specifically in sub-Saharan Africa, with a significant but lesser reduction in South-East Asia, the Eastern Mediterranean, and parts of Latin America. The greatest decreases in absolute deaths from neonatal sepsis are in India followed by sub-Saharan Africa; there were approximately equal decreases in the number of deaths in the other regions. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Projected global impact of a maternal K. pneumoniae vaccine. (A) Projected percent of neonatal deaths with respect to all neonatal deaths that could be averted due to a maternal K. pneumoniae vaccine. (B) Projected median number of neonatal sepsis deaths averted due to maternal vaccination. (C) Projected number of ceftazidime-resistant neonatal sepsis deaths averted due to maternal vaccination. (D) Projected number of meropenem-resistant neonatal sepsis deaths averted due to maternal vaccination. The maps are reprinted from pygal_maps_world under GNU GPL. https://doi.org/10.1371/journal.pmed.1004239.g003 In general, we project regions with the greatest reduction in K. pneumoniae neonatal sepsis are likely to experience the greatest reduction in neonatal deaths due to antimicrobial-resistant K. pneumoniae sepsis (Figs 3C and 3D and Figs J–L in S1 Text). However, whereas overall neonatal mortality would be considerably reduced across Africa and South Asia, the greatest reduction in resistant K. pneumoniae deaths would occur in Central Africa, East Africa, and India (Fig 3). Antibiotic resistance remains a significant issue in the treatment of K. pneumoniae because the projected number of both ceftazidime-resistant deaths (Fig 3C) and meropenem-resistant deaths (Fig 3D) is only marginally less than the projected number of total sepsis deaths (Fig 3B and S2 Table). There are likely to be more ceftazidime-resistant deaths than meropenem-resistant deaths averted globally. However, the relative burden of ceftazidime resistance is greater in North Africa while the burden of meropenem resistance is greater in South Asia. Regardless, for both drugs, we estimate high resistance rates worldwide and a significant number of deaths caused by drug-resistant infections could be averted through maternal vaccination. Finally, we project the benefits of a maternal vaccine in future scenarios of increased AMR prevalence. Available K. pneumoniae genomes indicate that the prevalence of carbapenemase ARGs (CRGs) is increasing significantly worldwide (Fig 4A) at a rate of 0.050 [95th percentile linear probability model confidence interval: 0.042 to 0.058, p < 0.001] ARGs being acquired per isolate yearly (Table D in S1 Text; R2 = 91.1%). Comparatively, aminoglycoside ARGs (AGRGs) increased significantly in the early 2000s (Fig 4A), but subsequently have not increased significantly (p = 0.436) and are constant at an average of 2.69 AGRGs per isolate. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Projected benefits of a maternal vaccine under future AMR scenarios. (A) Average number of ARGs that confer resistance to aminoglycosides (N = 24,359) and carbapenems (N = 5,189) per isolate by year of collection of K. pneumoniae isolates. (B) Projected median number of meropenem-resistant neonatal deaths averted by 2030 due to maternal vaccination based on the trends shown in Fig 4A. The maps are reprinted from pygal_maps_world under GNU GPL. AMR, antimicrobial resistance; ARG, antimicrobial resistance gene. https://doi.org/10.1371/journal.pmed.1004239.g004 Based on the increasing rates of CRGs, we extrapolated the trends from historic data to estimate the number of CRGs present in 2030 and calculated the associated future burden of meropenem-resistant K. pneumoniae deaths that could be averted by a potential maternal vaccine (Figs 4B and Fig M in S1 Text for the 95th percentile CrIs). The map both illustrates neonatal survival trends similar to those described above and the new regional benefits from vaccination as AMR increases. Countries that may experience only marginal benefits from vaccine distribution—such as those in the Americas including Venezuela and Colombia, those in the Western Pacific including Indonesia, those in the Eastern Mediterranean including Afghanistan, and those in South-East Asia including Bangladesh—are likely to experience greater reductions in neonatal mortality when accounting for increases in AMR. Effect of a K. pneumoniae maternal vaccine on neonatal deaths We initially focused on the 18 countries for which we had surveillance data on bacterial neonatal sepsis (Fig A in S1 Text). A potential vaccine could reduce neonatal sepsis and improve newborn survival (Fig 1). Our estimates of the proportion of neonatal deaths averted with maternal K. pneumoniae vaccination are consistent whether we use CHERG or GBD data as a reference source. Under the assumption that K. pneumoniae vaccine coverage matches current effective maternal tetanus vaccination coverage levels [30], a median estimated 57,392 [95th percentile Bayesian CrI: 22,248 to 125,803] deaths, accounting for 3.83% [CrI: 1.48 to 8.40] of newborn fatalities, could be averted each year in the initial 18 countries analyzed. We find a similar benefit in reducing neonatal sepsis cases, with an estimated 286,499 [CrI: 239,259 to 346,496] cases averted yearly. The benefits of a maternal K. pneumoniae vaccine vary across countries. While the relative impact of vaccination on neonatal mortality ranges a few percentage points (Fig 1A), the estimates for the reduction in total neonatal sepsis deaths and cases span orders of magnitude across countries, reflecting the variability in the burden of neonatal deaths across countries. We project that vaccination will have the greatest relative effect on neonatal mortality in Ethiopia, Rwanda, and Sierra Leone and a smaller effect in Greece, China, and Pakistan (Fig 1A), countries where a lower fraction of neonatal sepsis deaths were associated with K. pneumoniae (Fig A in S1 Text). However, although the relative impact of vaccination on neonatal mortality (i.e., percent of deaths averted with respect to all neonatal deaths) in Pakistan is smaller compared to other countries considered, we estimate that Pakistan will have one of the greatest absolute reductions in neonatal mortality (i.e., a high number of absolute deaths/cases averted; Fig 1B and 1C), consistent with previous work that Pakistan experiences a significant K. pneumoniae burden from nosocomial outbreaks [5]. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. A maternal vaccine reduces global K. pneumoniae burden. (A) Median estimated percent of neonatal deaths that are avertable with a maternal vaccine with respect to all neonatal deaths in a given country. (B) Median estimated number of neonatal deaths that are avertable with a maternal vaccine in a given country. (C) Median estimated number of neonatal cases of sepsis that are avertable with a maternal vaccine in a given country. GBD indicates data from the Global Burden of Disease study and CHERG indicates data from the United Nations Child Health Epidemiology Reference Group. Error bars indicate 95th percentile CrIs. A pseudo log transformation is done for values between 0 and 1 in panel B. CHERG, Child Health and Epidemiology Reference Group; CrI, credible interval; GBD, Global Burden of Disease. https://doi.org/10.1371/journal.pmed.1004239.g001 Across all WHO regions, K. pneumoniae isolates are highly resistant to most antibiotics (Fig 2). We focused on those antibiotics often prescribed to treat neonatal sepsis: ampicillin, gentamicin, ceftazidime, amikacin, and meropenem (Fig G in S1 Text for all antibiotics used to treat sepsis) [7,25]. Prevalence of resistance to ampicillin is highest: 89.82% [CrI: 64.97 to 99.49] of K. pneumoniae isolates are resistant to ampicillin worldwide (Fig 2B). The average rate of resistance to gentamicin, which is usually prescribed in combination with ampicillin as the first-line treatment for sepsis per current WHO guidelines [32], is 57.22% [CrI: 31.540 to 80.42] (Fig 2B). Comparatively, fewer isolates were resistant to ceftazidime and amikacin—frequently employed as a second-line combination to treat neonatal sepsis—with 44.63% [26.84 to 72.20] and 35.80% [14.69 to 58.53], respectively, isolates identified as resistant (Fig 2B). Approximately 22.43% [5.24 to 41.42] of isolates from neonates who died were resistant to meropenem. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. AMR distribution of vaccine-avertable deaths and cases. (A, B) Median percent of isolates that are resistant to specific classes of antibiotics by WHO region. Error bars indicate 95th percentile CrIs. (A) Of those deaths that are averted, the percent of deaths that are caused by K. pneumoniae resistant to each class of antibiotics. (B) Of those cases that are averted, the percent of illnesses that are caused by K. pneumoniae resistant to each class of antibiotics. (C, D) Median percent of isolates that are resistant to multiple classes of antibiotics. Error bars indicate 95th percentile CrIs. (C) Of those deaths that are averted, the percent of isolates that have varying degrees of MDR. (D) Of those cases that are averted, the percent of isolates that have varying degrees of MDR. (E, F) Joint distribution of ampicillin and gentamicin resistance. Error bars indicate 95th percentile CrIs. (E) Of those deaths that are averted, the percent of isolates that are resistant to combinations of aminoglycoside and beta lactam drugs. (F) Of those cases that are averted, the percent of isolates that are resistant to combinations of aminoglycoside and beta lactam drugs. WHO regions with no data indicate that no isolates were gathered in those regions, so no informative modeling can be done. AMR, antimicrobial resistance; CrI, credible interval; MDR, multidrug resistance; WHO, World Health Organization. https://doi.org/10.1371/journal.pmed.1004239.g002 Because ampicillin and gentamicin are prescribed together as the first-line treatment for neonatal sepsis in many countries, we considered the joint distribution of sepsis isolates that were resistant to neither, one, or both antibiotics (Fig 2C and 2D). Resistance is high against both drugs, with an average of 86.38% [CrI: 73.09 to 94.74] of isolates resistant to both drugs across all regions (Fig 2D). Notably, rates of resistance to ampicillin alone were substantially higher than to gentamicin alone (Fig 2D). In fact, less than 5% of isolates were resistant to just gentamicin (Fig 2D). Consequently, resistance to gentamicin almost guaranteed resistance to ampicillin (Fig 2D). Different antibiotics within the same antibiotic class display varying geographic trends. For instance, resistance to ampicillin is approximately constant across all geographic regions. Comparatively, rates of resistance to ceftazidime and meropenem are variable across geographic regions: Rates of resistance to ceftazidime in Africa are high and close to those of ampicillin, approximately 60 percentage points greater than rates of resistance to meropenem. Resistance to meropenem is low (7.01% [CrI: 3.21 to 12.64]) in Africa and substantially higher in South-East Asia (66.58% [CrI: 52.68 to 78.82]), greater than resistance to ceftazidime (Fig 2B). Among aminoglycosides, resistance rates to amikacin and gentamicin are similar in all regions except Africa where gentamicin resistance is 60 percentage points greater than amikacin resistance (Fig 2B). Finally, we consider trends across antibiotics within the same treatment regimen; in particular, we compare ampicillin and gentamicin (first-line treatment) against ceftazidime and amikacin (second-line treatment). While rates of resistance to ampicillin are consistently greater than to gentamicin, rates of resistance to ceftazidime are greater than amikacin only in Africa and the Western Pacific (Fig 2B). Across all regions, we generally did not find an increase in drug-resistant K. pneumoniae isolates from neonates who died compared with all K. pneumoniae isolates in general: The prevalence of resistance to a specific drug class is approximately equal whether the isolate was taken from a septic neonate who died or not (Fig 2A and 2B). However, there are 2 key exceptions: In South-East Asia, rates of resistance to amikacin, gentamicin, and ceftazidime are higher in neonates who died, and rates of resistance to ampicillin are lower among all neonates. Nevertheless, we do not find that isolates from neonates who succumbed to sepsis were more likely to exhibit AMR (Fig 2A and 2B), resistance to ampicillin or gentamicin (Fig 2C and 2D), or higher MDR (Fig 2E and 2F) in comparison to all isolates. Most K. pneumoniae have become resistant to multiple antibiotics (Fig 2E and 2F). On average, more than 99% of isolates were resistant to at least 1 drug. Many isolates were resistant to multiple antibiotics, and isolates that were resistant to none were rare (0.15% [CrI: 2.35 × 10−9 to 18.17]). The highest MDR rates are observed in the Eastern Mediterranean and South-East Asia where over 98% of K. pneumoniae are resistant to either 4 or 5 different antibiotics (Fig 2F); however, in Africa, most isolates are resistant to 3 antibiotics only with fewer isolates resistant to 4 or 5 antibiotics (Fig 2F). Projected global and future benefits of a K. pneumoniae maternal vaccine We expanded our estimates of the effects of a maternal K. pneumoniae vaccine on childhood mortality and neonatal sepsis to all countries using regression models as described above (Table C in S1 Text, R2 of 84.0%, p < 0.001). We project that 14.91% [CrI: 5.26 to 21.24] of all neonatal sepsis deaths in countries of interest (all countries except for high-income nations)—that is, 80,258 [CrI: 18,084 to 189,040] neonatal deaths or 3.40% [CrI: 0.75 to 8.01] of all neonatal mortality—and 399,015 [CrI: 334,523 to 485,442] neonatal sepsis cases could be averted yearly once vaccination coverage reaches that of the maternal tetanus vaccine (Figs 3A and 3B, and Figs H and I in S1 Text for the 95th percentile CrIs, and S2 Table for the numeric values). We project the greatest reduction in overall neonatal mortality (i.e., fraction of deaths reduced due to vaccination) to be in Africa, specifically in sub-Saharan Africa, with a significant but lesser reduction in South-East Asia, the Eastern Mediterranean, and parts of Latin America. The greatest decreases in absolute deaths from neonatal sepsis are in India followed by sub-Saharan Africa; there were approximately equal decreases in the number of deaths in the other regions. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Projected global impact of a maternal K. pneumoniae vaccine. (A) Projected percent of neonatal deaths with respect to all neonatal deaths that could be averted due to a maternal K. pneumoniae vaccine. (B) Projected median number of neonatal sepsis deaths averted due to maternal vaccination. (C) Projected number of ceftazidime-resistant neonatal sepsis deaths averted due to maternal vaccination. (D) Projected number of meropenem-resistant neonatal sepsis deaths averted due to maternal vaccination. The maps are reprinted from pygal_maps_world under GNU GPL. https://doi.org/10.1371/journal.pmed.1004239.g003 In general, we project regions with the greatest reduction in K. pneumoniae neonatal sepsis are likely to experience the greatest reduction in neonatal deaths due to antimicrobial-resistant K. pneumoniae sepsis (Figs 3C and 3D and Figs J–L in S1 Text). However, whereas overall neonatal mortality would be considerably reduced across Africa and South Asia, the greatest reduction in resistant K. pneumoniae deaths would occur in Central Africa, East Africa, and India (Fig 3). Antibiotic resistance remains a significant issue in the treatment of K. pneumoniae because the projected number of both ceftazidime-resistant deaths (Fig 3C) and meropenem-resistant deaths (Fig 3D) is only marginally less than the projected number of total sepsis deaths (Fig 3B and S2 Table). There are likely to be more ceftazidime-resistant deaths than meropenem-resistant deaths averted globally. However, the relative burden of ceftazidime resistance is greater in North Africa while the burden of meropenem resistance is greater in South Asia. Regardless, for both drugs, we estimate high resistance rates worldwide and a significant number of deaths caused by drug-resistant infections could be averted through maternal vaccination. Finally, we project the benefits of a maternal vaccine in future scenarios of increased AMR prevalence. Available K. pneumoniae genomes indicate that the prevalence of carbapenemase ARGs (CRGs) is increasing significantly worldwide (Fig 4A) at a rate of 0.050 [95th percentile linear probability model confidence interval: 0.042 to 0.058, p < 0.001] ARGs being acquired per isolate yearly (Table D in S1 Text; R2 = 91.1%). Comparatively, aminoglycoside ARGs (AGRGs) increased significantly in the early 2000s (Fig 4A), but subsequently have not increased significantly (p = 0.436) and are constant at an average of 2.69 AGRGs per isolate. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Projected benefits of a maternal vaccine under future AMR scenarios. (A) Average number of ARGs that confer resistance to aminoglycosides (N = 24,359) and carbapenems (N = 5,189) per isolate by year of collection of K. pneumoniae isolates. (B) Projected median number of meropenem-resistant neonatal deaths averted by 2030 due to maternal vaccination based on the trends shown in Fig 4A. The maps are reprinted from pygal_maps_world under GNU GPL. AMR, antimicrobial resistance; ARG, antimicrobial resistance gene. https://doi.org/10.1371/journal.pmed.1004239.g004 Based on the increasing rates of CRGs, we extrapolated the trends from historic data to estimate the number of CRGs present in 2030 and calculated the associated future burden of meropenem-resistant K. pneumoniae deaths that could be averted by a potential maternal vaccine (Figs 4B and Fig M in S1 Text for the 95th percentile CrIs). The map both illustrates neonatal survival trends similar to those described above and the new regional benefits from vaccination as AMR increases. Countries that may experience only marginal benefits from vaccine distribution—such as those in the Americas including Venezuela and Colombia, those in the Western Pacific including Indonesia, those in the Eastern Mediterranean including Afghanistan, and those in South-East Asia including Bangladesh—are likely to experience greater reductions in neonatal mortality when accounting for increases in AMR. Discussion K. pneumoniae is the leading cause of newborn sepsis deaths worldwide (Figs O–U in S1 Text); however, no K. pneumoniae vaccine is currently approved. We estimate the global benefits of a potential maternal K. pneumoniae vaccine in reducing neonatal mortality and AMR. A K. pneumoniae vaccine has the potential to significantly improve child wellbeing, reduce a large fraction (15%) of neonatal deaths and sepsis infections, and will benefit additional geographic areas as AMR is projected to increase. Most neonatal sepsis cases are multidrug resistant, further underscoring the benefits of vaccination. Nevertheless, our modeling suggests that the bacteria causing neonatal sepsis in neonates that eventually die do not have higher antimicrobial or MDR rates. Instead, resistance is high in all isolates and against all drugs, underscoring neonatal fragility and suggesting that other factors such as access to effective drugs and antibiotic dosing regiments or standard of care may be driving differential mortality [7]. Regardless, our estimates indicate that maternal vaccination against K. pneumoniae would have global benefits and improve efforts toward reaching child survival targets in LMICs worldwide. Our methodology represents an advantage over using just clinical data or a single study to draw conclusions and make estimates: We applied a Bayesian framework to aggregate data from various studies and generate increased predictive power. This approach more accurately reflects the underlying uncertainty and allows meaningful inference and projections, even when there are few reported and annotated cases of culture-confirmed neonatal sepsis. However, further work is required to develop a vaccine that is approved, widely distributable, and affordable. Most K. pneumoniae neonatal sepsis cases occur in premature infants [26], but the success of a maternal vaccine will depend on robust transplacental antibody transfer from mother to a developing fetus, which peaks late in the third trimester of pregnancy. Additional evaluation of the adequacy of maternal antibody transfer for premature infants is important to inform maternal vaccine development efforts. There are several potential caveats in our study. Our estimates of averted cases and deaths assume that vaccination rates instantly reach their maximum in all countries—that is, we do not model the reduction in cases of K. pneumoniae sepsis as vaccine rollout progresses. To address this, we have conducted sensitivity analyses where we vary the vaccine efficacy/coverage (S2 Table). Our estimates for the number of avertable deaths respond linearly to vaccine efficacy and coverage. Other studies that estimate the effect of a hypothetical intervention against a communicable disease have used a similar assumption [33,34]. As more details of the potential K. pneumoniae vaccine are determined, we can refine our model accordingly. Since little is known about the factors that drive nosocomial infections and hospital-acquired neonatal sepsis, the regression models we used to produce estimates of averted cases in all countries are simple and consider only health variables. Nevertheless, previous work that has attempted to extrapolate disease burden and rates of AMR across countries has used a similar econometric-like approach [29]. Finally, our estimates of the CFR of K. pneumoniae neonatal sepsis may be biased due to loss-to-follow-up: In particular, we are unable to account for whether there is a sampling bias of a disproportionate number of deaths being recorded. Further data is required to address this issue. The uncertainty in our estimates largely comes from uncertainty in the underlying values for annual neonatal sepsis deaths from CHERG and GBD, which is considerable. Although estimates for the percentage of averted deaths are consistent (Fig 1A), estimates for the number of avertable deaths using data from CHERG as opposed to GBD often differ by an order of magnitude (Fig 1B). This is because of underlying differences in the number of neonatal sepsis deaths from the 2 sources. Estimates using CHERG data consistently indicate a higher number, though estimates derived from GBD data are still substantial. In many cases, the CrIs from the estimates using different data sources overlap. To address this issue, better country-level data on the burden of neonatal sepsis are needed. Moreover, the sample sizes from all 3 studies were limited and we were unable to delineate hospital outbreaks of K. pneumoniae caused by a single strain that may have differential vaccine aversion. We have aimed to address this limitation using the Bayesian approach and noninformative prior, as discussed above, to enable meaningful inference while considering latent noise in the data. Finally, because of the limited data available on the etiologies causing neonatal sepsis, our analysis only considers country-level trends for K. pneumoniae neonatal sepsis cases and WHO regional trends for AMR. However, there is heterogeneity within a country and region; our modeling is unable to capture that. We do not know whether heterogeneity within a country is due to external factors or rather indicates study-level variability: Given our limited data and few countries for which we have multiple study sites, we do not model study-level variability, though this is an area of future work. Additional work is required to better understand within-country and within-region trends: Preliminary data indicates variability in sepsis levels within a country and across private and public hospitals [5,35,36]. Moreover, we may be subject to a sampling bias because the data from the 3 studies we used are point estimates of the impact of sepsis at specific healthcare facilities within a country rather than surveillance over the whole country. In particular, we must limit our analysis to cases of culture-confirmed K. pneumoniae, but doing so makes our estimates a probable lower bound on the actual number of averted deaths and cases: it is likely that not all cases of neonatal sepsis are identified as such and cultured. Even if they are cultured, they may not test positive due to low blood volume, a bias whose impact across study sites likely varies appreciably but is unknown [37]. Additional work is required to better understand the factors that impact neonatal sepsis rates within a region. We also lacked data from many countries. We tried to minimize the effect of this by projecting the number of averted cases by country of maternal vaccination on antimicrobial-resistant K. pneumoniae neonatal sepsis on a regional basis and limiting our extrapolations to LMICs. Moreover, we used a uniform prior over the probability of mortality from K. pneumoniae to consider the variability in K. pneumoniae mortality across all countries in a region: This ensures that the CrIs for our projections to other countries include other mortality rates that may be the true underlying mortality rate for that country. Nonetheless, additional data on K. pneumoniae neonatal sepsis across more countries is required to reduce the uncertainty in our estimates. Future work includes a cost-effectiveness analysis to estimate the economic and social implications of maternal vaccination. A K. pneumoniae vaccine would reduce the economic burden of infection for both the patient and the hospital: quantifying these impacts is critical to understanding the greater societal impacts of this vaccine, especially as the burden of sepsis is highest in LMICs that have larger health expenditure (Fig N in S1 Text). Additionally, K. pneumoniae is known to be the leading driver of sepsis, especially in the late neonatal stage: reducing sepsis during this timeframe would be significantly advantageous to maternal admissions and mitigating healthcare burden as it would reduce the length of prolonged hospital stays. Beyond reducing cases and deaths, a K. pneumoniae vaccine would likely also reduce antibiotic usage and thus may help improve antibiotic stewardship in LMICs and reduce AMR rates. It is noteworthy that in many LMICs, antibiotic availability is challenging (particularly for the more potent drugs) and the cost is in many cases deferred to the family; therefore, a maternal vaccine will also alleviate drug demand and provide local financial benefits. Moreover, vaccination may reduce the severity of disease and have positive social ramifications by reducing stigma surrounding infection and helping mothers who would otherwise leave their jobs to care for sick newborns [38]. We have highlighted the main benefits of the rollout of a potential K. pneumoniae maternal vaccine: reducing global cases and deaths of resistance and susceptible neonatal sepsis, reducing overall neonatal mortality, and improving childhood health. Supporting information S1 Text. Supplementary materials. Fig A. Raw data. Calculated percent of neonatal sepsis deaths that are associated (i.e., an isolate from the neonate who died was culture-positive) with various etiologies across each study by location. Table A. Number of neonates who died of neonatal sepsis divided by number of neonates surveilled by location. Fig B. Flow diagram summarizing cases of culture-confirmed sepsis used in the main analysis of vaccine-avertable sepsis and AMR. Fig C. Flow diagram summarizing data collection and cleaning of the Klebsiella pneumoniae genomes used in the antimicrobial resistance genes (ARG) prevalence analysis. BARNARDS refers to data gathered from the Burden of Antimicrobial Resistance in Neonates in Developing Societies study. Fig D. Distribution of available genomic data from PathogenWatch and the Burden of Antimicrobial Resistance in Neonates from Developing Societies study by year used in the antimicrobial resistance gene prevalence analysis. Fig E. Tree map of the distribution of available K. pneumoniae isolates across countries for use in the prevalence of antimicrobial resistance genes analysis. Colors have no meaning and are used to create contrast between countries. Fig F. Schematic diagram of the modeling framework. Table B. Model parameters. Note that this refers to values that are used as inputs to various modeling stages, not quantities that are predicted through the modeling process described in Fig F. Fig G. Raw resistance data by study. Table C. Regression analysis results for the model used to extrapolate the number of averted deaths from those countries for which we have data to all countries. Fig H. 2.5th percentile (i.e., lower bound) of estimates represented as maps shown in Fig 3. The maps are reprinted from pygal_maps_world under GNU GPL. Fig I. 97.5th percentile (i.e., lower bound) of estimates represented as maps shown in Fig 3. The maps are reprinted from pygal_maps_world under GNU GPL. Fig J. As Fig 3C and 3D but for ampicillin. Median estimates shown on top. 2.5th percentile shown in middle. 97.5th percentile shown on bottom. The maps are reprinted from pygal_maps_world under GNU GPL. Fig K. As Fig 3C and 3D but for gentamicin. Median estimates shown on top. 2.5th percentile shown in middle. 97.5th percentile shown on bottom. The maps are reprinted from pygal_maps_world under GNU GPL. Fig L. As Fig 3C and 3D but for amikacin. Median estimates shown on top. 2.5th percentile shown in middle. 97.5th percentile shown on bottom. The maps are reprinted from pygal_maps_world under GNU GPL. Table D. Regression analysis results for the model used to estimate the yearly rate of increase in antimicrobial resistance genes. Fig M. Credible interval of map shown in Fig 4B. 2.5th percentile shown on top and 97.5th percentile shown on bottom. The maps are reprinted from pygal_maps_world under GNU GPL. Fig N. Health expenditure as a fraction of the country’s GDP. Data courtesy of the WHO, Global Health Observatory (2022). Map reprinted from OurWorldInData under a CC-BY license. Original: https://ourworldindata.org/grapher/total-healthcare-expenditure-gdp. Fig O. As Fig 1A but for other etiologies of interest. Median estimated fraction of neonatal deaths averted given maternal vaccination against a specific pathogen at 70% efficacy and coverage equivalent to that of the maternal tetanus vaccine. Median shown; error bars indicate 95th percentile Bayesian credible intervals. GBD refers to data from the Global Burden of Disease study, and CHERG refers to data from the Child Health and Epidemiology Reference Group. Fig P. As Fig 1B but for other etiologies of interest. Median estimated number of avertable neonatal sepsis deaths given maternal vaccination against a specific pathogen at 70% efficacy and coverage equivalent to that of the maternal tetanus vaccine. Median shown; error bars indicate 95th percentile Bayesian credible intervals. A pseudo log transform is done for values between 0 and 1. GBD refers to data from the Global Burden of Disease study, and CHERG refers to data from the Child Health and Epidemiology Reference Group. Fig Q. As Fig 1C but for other etiologies of interest. Median estimated number of avertable neonatal sepsis cases given maternal vaccination against a specific pathogen at 70% efficacy and coverage equivalent to that of the maternal tetanus vaccine. Median shown; error bars indicate 95th percentile Bayesian credible intervals. GBD refers to data from the Global Burden of Disease study, and CHERG refers to data from the Child Health and Epidemiology Reference Group. Fig R. As Fig 2A but for other etiologies and relevant antibiotics of interest. Estimated median fraction of isolates from neonates who died with culture-confirmed sepsis that are resistant to various drugs across WHO regions. Median shown; error bars indicate 95th percentile Bayesian credible intervals. Fig S. As Fig 2B but for other etiologies and relevant antibiotics of interest. Estimated median fraction of isolates from neonates with culture-confirmed sepsis that are resistant to various drugs across WHO regions. Median shown; error bars indicate 95th percentile Bayesian credible intervals. Fig T. As Fig 2E but for other etiologies of interest. Antibiotics considered are shown in Figs R and S. Median shown; error bars indicate 95th percentile Bayesian credible intervals. Fig U. As Fig 2F but for other etiologies of interest. Antibiotics considered are shown in Figs R and S. Median shown; error bars indicate 95th percentile Bayesian credible intervals. https://doi.org/10.1371/journal.pmed.1004239.s001 (PDF) S2 Text. Extended methods. Derivation of the mixture model posterior distribution (Section 1) and details of antimicrobial susceptibility testing done (Section 2). https://doi.org/10.1371/journal.pmed.1004239.s002 (PDF) S1 Table. Details of the analyzed genomes used in this study including FASTA file ID, year, source, country of origin, and number of identified antimicrobial resistance genes. https://doi.org/10.1371/journal.pmed.1004239.s003 (XLSX) S2 Table. Impact of maternal vaccination by country. Country-specific estimates of overall averted cases, deaths, and result on neonatal mortality as well as antibiotic resistant deaths. https://doi.org/10.1371/journal.pmed.1004239.s004 (XLSX) Acknowledgments We thank the BARNARDS and NeoObs groups for sharing their data. In particular, TRW would like to personally thank Kathryn Thomson and Rebecca Milton, both part of the BARNARDS group. We thank Dianna Blau from the CHAMPS group for valuable feedback. We would also like to thank each CHAMPS site team. We thank Nicole Benson, at the Bill & Melinda Gates Foundation, Mateusz Hasso-Agopsowicz and Isabel Frost, both at the World Health Organization, and Sonya Davey, at Brigham and Women’s Hospital for valuable feedback. We thank Sally Atwater for valuable editorial suggestions on drafts.
SARS-CoV-2 transmission with and without mask wearing or air cleaners in schools in Switzerland: A modeling study of epidemiological, environmental, and molecular dataBanholzer, Nicolas;Zürcher, Kathrin;Jent, Philipp;Bittel, Pascal;Furrer, Lavinia;Egger, Matthias;Hascher, Tina;Fenner, Lukas
doi: 10.1371/journal.pmed.1004226pmid: 37200241
Background Growing evidence suggests an important contribution of airborne transmission to the overall spread of Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2), in particular via smaller particles called aerosols. However, the contribution of school children to SARS-CoV-2 transmission remains uncertain. The aim of this study was to assess transmission of airborne respiratory infections and the association with infection control measures in schools using a multiple-measurement approach. Methods and findings We collected epidemiological (cases of Coronavirus Disease 2019 (COVID-19)), environmental (CO2, aerosol and particle concentrations), and molecular data (bioaerosol and saliva samples) over 7 weeks from January to March 2022 (Omicron wave) in 2 secondary schools (n = 90, average 18 students/classroom) in Switzerland. We analyzed changes in environmental and molecular characteristics between different study conditions (no intervention, mask wearing, air cleaners). Analyses of environmental changes were adjusted for different ventilation, the number of students in class, school and weekday effects. We modeled disease transmission using a semi-mechanistic Bayesian hierarchical model, adjusting for absent students and community transmission. Molecular analysis of saliva (21/262 positive) and airborne samples (10/130) detected SARS-CoV-2 throughout the study (weekly average viral concentration 0.6 copies/L) and occasionally other respiratory viruses. Overall daily average CO2 levels were 1,064 ± 232 ppm (± standard deviation). Daily average aerosol number concentrations without interventions were 177 ± 109 1/cm3 and decreased by 69% (95% CrI 42% to 86%) with mask mandates and 39% (95% CrI 4% to 69%) with air cleaners. Compared to no intervention, the transmission risk was lower with mask mandates (adjusted odds ratio 0.19, 95% CrI 0.09 to 0.38) and comparable with air cleaners (1.00, 95% CrI 0.15 to 6.51). Study limitations include possible confounding by period as the number of susceptible students declined over time. Furthermore, airborne detection of pathogens document exposure but not necessarily transmission. Conclusions Molecular detection of airborne and human SARS-CoV-2 indicated sustained transmission in schools. Mask mandates were associated with greater reductions in aerosol concentrations than air cleaners and with lower transmission. Our multiple-measurement approach could be used to continuously monitor transmission risk of respiratory infections and the effectiveness of infection control measures in schools and other congregate settings. Why was this study done? Public health authorities worldwide closed businesses and schools during the Coronavirus Disease 2019 (COVID-19) pandemic. The closure of schools has been most intensely debated. The contribution of school children to Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) transmission and the role of school rooms remain unknown. What did the researchers do and find? We used molecular, environmental, and epidemiological data to understand the transmission of the virus causing COVID-19 (SARS-CoV-2) in 2 secondary schools (90 students) in Switzerland in the presence and absence of mask wearing and air cleaners. We detected SARS-CoV-2 in aerosols in the air and saliva samples from the students throughout the study. Aerosol and particle concentrations were on average 70% lower with mask mandates and 40% lower with air cleaners. The transmission model estimated that between 2 and 19 infections could be avoided during the study period with mask wearing. What do these findings mean? Molecular analyses indicated sustained airborne SARS-CoV-2 transmission. Mask wearing may be more effective than air cleaners in reducing aerosol concentrations and transmission of SARS-CoV-2. This approach can be used to assess transmission dynamics and the effectiveness of infection control measures in reducing transmission of respiratory infections during future epidemics. 1 Introduction The spread of respiratory infections such as Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) and influenza is difficult to control [1]. Person-to-person transmission mainly occurs through exhaled respiratory particles containing the respective pathogen, particularly via aerosols (defined as respiratory particles <100 μm [2,3]), rather than through larger droplets >100 μm. Multiple reports have provided evidence that a considerable part of SARS-CoV-2 transmission is likely to happen through small respiratory particles (<5 μm, also called fine aerosols), allowing for longer suspension times and airborne transmission at short (1 to 2 m) and long ranges (>2 m) [4–6]. Growing evidence suggests they contribute importantly to the overall spread of SARS-CoV-2 in indoor congregate settings such as clinics, workplaces, and schools [3,6–8]. Public authorities worldwide closed businesses and schools during the Coronavirus Disease 2019 (COVID-19) pandemic [9,10], but the closure of schools was particularly contentious. While the negative impact of school closures on student well-being and mental health is well documented [11], the role of children and adolescents in transmitting SARS-CoV-2 is less clear [12]. A study in Germany estimated that school contacts contributed between 2% and 20% to the overall transmission of SARS-CoV-2 in the population [13]. Studies on the effectiveness of government interventions at the population level are inconclusive regarding the effects of school closures in different epidemic waves [10,14]. The introduction of compulsory face mask wearing [13,15] and improved ventilation [15] in schools was associated with a lower incidence of COVID-19. In addition, the installation of portable HEPA-air filtration devices (air cleaners) was shown to remove SARS-CoV-2 bioaerosols in a hospital ward [16]. Finally, a simulation study of exhaled SARS-CoV-2 bioaerosols in an indoor space demonstrated the efficacy of mask wearing and portable air cleaners in reducing aerosols [17]. We used a multiple-measurement approach to study the transmission of SARS-CoV-2 and other respiratory viruses in school rooms. We hypothesized that our approach could demonstrate transmission (indicated by exposure to bioaerosols and epidemiological data) and that infection control measures (mask wearing and air cleaners) would reduce the concentration of aerosols and respiratory particles and lower the risk of infection for students in school rooms. We collected epidemiological (cases of respiratory diseases), environmental (CO2, aerosol, and particle concentrations), and molecular data (detection of respiratory viruses in bioaerosol and human saliva samples) over a seven-week study period from January to March 2022 in 2 secondary schools in Switzerland. We analyzed changes in environmental and molecular characteristics and estimated the probability of infection with SARS-CoV-2 during 3 different study conditions with and without infection control measures (general mask wearing and air cleaners). 2 Methods 2.1 Study setting and design We collected data in 2 secondary schools (age of students 13 to 15 years) over a seven-week study period from January 24 to March 26 (School 1) and March 18 (School 2), 2022. Both schools are located in the Canton of Solothurn, Switzerland, and have 1,500 (School 1) and 700 (School 2) students. Epidemiological data were collected in the 5 classes (School 1: classes A/B, C; School 2: classes D, E), and environmental and molecular data were collected in 2 classrooms (School 1: A/B, School 2: D). In School 1, the same classroom was used by 2 classes A/B due to half-class teaching. Fig 1 shows the schematic study setup. This study is reported as per the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) guideline (S1 Checklist). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. Study setting. Schematic study setup of classrooms where environmental data was collected in each school. One air cleaner was placed in the front and the other in the back of classrooms. All devices were placed at the head level of students when they were seated. Both classrooms were not equipped with an active HVAC (heating, ventilation, air conditioning) system, but were ventilated using passive window ventilation. For School 1, ventilation was additionally assisted by a CO2 guided opener of a small window at the top. https://doi.org/10.1371/journal.pmed.1004226.g001 2.2 Study interventions We distinguished 3 conditions (Fig 2): (i) wearing face masks as mandated by the public health authorities at that time (Mask mandate; typically Type II and Type IIR masks, although community masks were also allowed); (ii) standard condition following the lifting of mask mandates (No intervention); and (iii) environmental intervention using commercially available portable HEPA- filtration devices (Air cleaner; Xiaomi Mi Air Pro 70m2, Shenzhen, China; approx. USD 250 per device, running at 2 × 600 m3/h clean air delivery rate). Mask mandates applied to all classes (including teachers) and were generally well adhered to. In School 2, mask wearing continued for 1 week after the mandate was lifted (week 4). Air cleaners were only installed in 2 classrooms with bioaerosol and environmental sampling. Passive window ventilation occurred per recommendations of the national public health authorities during all study conditions. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. Study design and community transmission during the study period. (a) Study conditions over the seven-week study period. (b) Number of new cases of COVID-19 across all age groups (7-day moving average) and the reproduction number (average of the median published estimate) in the Canton of Solothurn [18] over the study period. https://doi.org/10.1371/journal.pmed.1004226.g002 2.3 Data collection Epidemiological data. At baseline, we collected aggregated data on age, sex, and COVID-19 vaccination status in the participating classes. We inferred data on the number of suspected and confirmed cases of COVID-19 based on the reported number of students absent from school due to sickness probably related to COVID-19 (Text C and Tables A–E in S1 Appendix). Reports about absences were entered electronically into a REDCap database [19,20]. Both schools participated in repetitive weekly testing for SARS-CoV-2. The saliva samples were transported to the laboratory and stored at −80°C until further processed [21–23]. Environmental data. CO2 and particle measurements. An air quality device (AQ Guard, Palas GmbH, Karlsruhe, Germany) continuously measured indoor CO2 levels, aerosol number concentrations (particle diameter between 175 nm to 20 μm) and particle mass concentrations (PM in μgm−3; PM1, PM2.5, PM4, PM10, i.e., particles of sizes <1 to <10 μm). This device has been used in previous work [24,25]. Data were filtered according to the times students were in the classroom, which were monitored together with the time that windows were opened. Bioaerosol sampling. We collected airborne respiratory viruses in the classroom with a bioaerosol sampling device (BioSpot-VIVAS, Aerosol Devices, Ft. Collins, Colorado, United States of America). This device samples airborne virus particles into a viral transport medium (VTM) using a laminar flow water-based condensation method. In parallel, we also used the Coriolis Micro Air (Bertin Instruments Montigny-le-Bretonneux, France) sampler, running at 200 l/min and collecting into 15 mL PBS as previously reported in clinic settings [26]. BioSpot-VIVAS operated throughout lessons while Coriolis Micro Air only operated shortly before and during break times to reduce noise exposure (approximately 60 min/day). The removable parts of both sampling devices were regularly autoclaved to avoid contamination. At the end of the day, samples were transported to the Institute for Infectious Diseases (IFIK) and stored at −80°C. At the end of the study period, the Xiaomi HEPA filters were carefully removed and 20 swabs were taken from each filter and stored at −80°C until further processed. 2.4 Laboratory and molecular analyses Prior to the real-time (RT)-PCR analysis, daily bioaerosol samples were combined for each sampling device and filtered using Amicon Ultra-15 Centrifugal Filters with Ultracel 10,000 Dalton molecular weight cutoffs filters (UFC9010; MilliporeSigma, Burlington, USA) to a volume of 1 mL. The human saliva samples were directly analyzed without prior filtration. The Allplex RV Master Assay (Seegene, Seoul, South Korea) was used to detect a panel of 19 respiratory viruses (Text A in S1 Appendix), including SARS-CoV-2. Viral genome load (VGL) of specimens was quantified using standardized dilution series and reported in genomic equivalent copies/L. For positive samples, we performed targeted sequencing of viral RNA to compare genetic relatedness between SARS-CoV-2 strains detected in the air and human samples [27]. However, we were unable to amplify and sequence any of the gene targets in the bioaerosol samples due to low RNA concentrations. 2.5 Statistical analyses and modeling Number of new cases of COVID-19 by date of symptom onset. The daily number of new cases of COVID-19 was inferred based on the number of students absent from school. Confirmed cases referred to absences due to a positive laboratory test result or isolation. Suspected cases referred to absences due to sickness probably related to COVID-19. Absences due to other known reasons were excluded. We used a probabilistic simulation approach (Text D in S1 Appendix) to estimate the number of suspected cases that were cases of COVID-19 and the date of symptom onset, which was not always reported by students. We generated 100 datasets for the daily number of new cases of COVID-19. Subsequent analyses were performed on each of these datasets and we report the mean of estimation results if not stated otherwise. Cases in teachers were excluded as teachers taught multiple classes and had different exposure. Aerosol and particle concentrations. We computed the average concentrations per day and compared these between study conditions. We report the mean ± standard deviation. We further estimated the reduction in concentrations using Bayesian log-linear regression models (Text H in S1 Appendix), adjusting for the ventilation rate (computed from indoor CO2 levels; see Text I in S1 Appendix), the daily number of students in class, school, and weekday effects. Risk of transmission. We estimated daily transmission of SARS-CoV-2 with a Bayesian semi-mechanistic hierarchical model [9] (Fig 3): (i) We modeled the number of new infections as a function of susceptible students in each class and day, where the probability of infection can vary by study condition. (ii) We adjusted the estimated effects of interventions for the daily proportion of all absent students and the effective reproduction number in the community. (iii) The number of new cases was computed as the weighted sum of the number of new infections in the previous days. (iv) The number of susceptible students was computed by removing the number of students who have already been infected. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Visual summary of the structure of the Bayesian semi-mechanistic hierarchical model. (1) The number of new infections was modeled as a function of the number of susceptible students, the class-specific daily rate of infections, and the reductions from active infection control measures; (2) the effects of control measures were adjusted for transmission in the community and the proportion of students in class; (3) the observed number of new cases was computed as the weighted sum of the number of new infections in the previous days; and (4) the number of susceptible students was computed as the total number of students minus the cumulative sum of infections in the previous days. https://doi.org/10.1371/journal.pmed.1004226.g003 A detailed description of the transmission model and choice of priors for all model parameters are provided in Text E in S1 Appendix. Model parameters were estimated with a Bayesian approach (Text F in S1 Appendix). Specifically, Markov chain Monte Carlo (MCMC) sampling was used as implemented by the Hamiltonian Monte Carlo algorithm with the No-U-Turn Sampler (NUTS) [28]. If not stated otherwise, we report posterior means and credible intervals (CrIs) based on the 50%, 80%, and 95% quantiles of the posterior samples, respectively. We estimated the total number of avoided infections for each intervention by computing the difference between the estimated number of infections in the presence and absence of interventions (counterfactual scenario). Software. All analyses were performed in R software (version 4.2.0) [29] and modeling in Stan (version 2.21.0) [30]. The code is available from https://github.com/nbanho/mcid. 2.6 Ethics statement The Ethics Committee of the Canton of Bern, Switzerland, approved the study (reference no. 2021–02377). For the saliva samples, we included all students who were willing to participate and obtained written informed consent from their caregivers. 2.1 Study setting and design We collected data in 2 secondary schools (age of students 13 to 15 years) over a seven-week study period from January 24 to March 26 (School 1) and March 18 (School 2), 2022. Both schools are located in the Canton of Solothurn, Switzerland, and have 1,500 (School 1) and 700 (School 2) students. Epidemiological data were collected in the 5 classes (School 1: classes A/B, C; School 2: classes D, E), and environmental and molecular data were collected in 2 classrooms (School 1: A/B, School 2: D). In School 1, the same classroom was used by 2 classes A/B due to half-class teaching. Fig 1 shows the schematic study setup. This study is reported as per the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) guideline (S1 Checklist). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 1. Study setting. Schematic study setup of classrooms where environmental data was collected in each school. One air cleaner was placed in the front and the other in the back of classrooms. All devices were placed at the head level of students when they were seated. Both classrooms were not equipped with an active HVAC (heating, ventilation, air conditioning) system, but were ventilated using passive window ventilation. For School 1, ventilation was additionally assisted by a CO2 guided opener of a small window at the top. https://doi.org/10.1371/journal.pmed.1004226.g001 2.2 Study interventions We distinguished 3 conditions (Fig 2): (i) wearing face masks as mandated by the public health authorities at that time (Mask mandate; typically Type II and Type IIR masks, although community masks were also allowed); (ii) standard condition following the lifting of mask mandates (No intervention); and (iii) environmental intervention using commercially available portable HEPA- filtration devices (Air cleaner; Xiaomi Mi Air Pro 70m2, Shenzhen, China; approx. USD 250 per device, running at 2 × 600 m3/h clean air delivery rate). Mask mandates applied to all classes (including teachers) and were generally well adhered to. In School 2, mask wearing continued for 1 week after the mandate was lifted (week 4). Air cleaners were only installed in 2 classrooms with bioaerosol and environmental sampling. Passive window ventilation occurred per recommendations of the national public health authorities during all study conditions. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 2. Study design and community transmission during the study period. (a) Study conditions over the seven-week study period. (b) Number of new cases of COVID-19 across all age groups (7-day moving average) and the reproduction number (average of the median published estimate) in the Canton of Solothurn [18] over the study period. https://doi.org/10.1371/journal.pmed.1004226.g002 2.3 Data collection Epidemiological data. At baseline, we collected aggregated data on age, sex, and COVID-19 vaccination status in the participating classes. We inferred data on the number of suspected and confirmed cases of COVID-19 based on the reported number of students absent from school due to sickness probably related to COVID-19 (Text C and Tables A–E in S1 Appendix). Reports about absences were entered electronically into a REDCap database [19,20]. Both schools participated in repetitive weekly testing for SARS-CoV-2. The saliva samples were transported to the laboratory and stored at −80°C until further processed [21–23]. Environmental data. CO2 and particle measurements. An air quality device (AQ Guard, Palas GmbH, Karlsruhe, Germany) continuously measured indoor CO2 levels, aerosol number concentrations (particle diameter between 175 nm to 20 μm) and particle mass concentrations (PM in μgm−3; PM1, PM2.5, PM4, PM10, i.e., particles of sizes <1 to <10 μm). This device has been used in previous work [24,25]. Data were filtered according to the times students were in the classroom, which were monitored together with the time that windows were opened. Bioaerosol sampling. We collected airborne respiratory viruses in the classroom with a bioaerosol sampling device (BioSpot-VIVAS, Aerosol Devices, Ft. Collins, Colorado, United States of America). This device samples airborne virus particles into a viral transport medium (VTM) using a laminar flow water-based condensation method. In parallel, we also used the Coriolis Micro Air (Bertin Instruments Montigny-le-Bretonneux, France) sampler, running at 200 l/min and collecting into 15 mL PBS as previously reported in clinic settings [26]. BioSpot-VIVAS operated throughout lessons while Coriolis Micro Air only operated shortly before and during break times to reduce noise exposure (approximately 60 min/day). The removable parts of both sampling devices were regularly autoclaved to avoid contamination. At the end of the day, samples were transported to the Institute for Infectious Diseases (IFIK) and stored at −80°C. At the end of the study period, the Xiaomi HEPA filters were carefully removed and 20 swabs were taken from each filter and stored at −80°C until further processed. Epidemiological data. At baseline, we collected aggregated data on age, sex, and COVID-19 vaccination status in the participating classes. We inferred data on the number of suspected and confirmed cases of COVID-19 based on the reported number of students absent from school due to sickness probably related to COVID-19 (Text C and Tables A–E in S1 Appendix). Reports about absences were entered electronically into a REDCap database [19,20]. Both schools participated in repetitive weekly testing for SARS-CoV-2. The saliva samples were transported to the laboratory and stored at −80°C until further processed [21–23]. Environmental data. CO2 and particle measurements. An air quality device (AQ Guard, Palas GmbH, Karlsruhe, Germany) continuously measured indoor CO2 levels, aerosol number concentrations (particle diameter between 175 nm to 20 μm) and particle mass concentrations (PM in μgm−3; PM1, PM2.5, PM4, PM10, i.e., particles of sizes <1 to <10 μm). This device has been used in previous work [24,25]. Data were filtered according to the times students were in the classroom, which were monitored together with the time that windows were opened. Bioaerosol sampling. We collected airborne respiratory viruses in the classroom with a bioaerosol sampling device (BioSpot-VIVAS, Aerosol Devices, Ft. Collins, Colorado, United States of America). This device samples airborne virus particles into a viral transport medium (VTM) using a laminar flow water-based condensation method. In parallel, we also used the Coriolis Micro Air (Bertin Instruments Montigny-le-Bretonneux, France) sampler, running at 200 l/min and collecting into 15 mL PBS as previously reported in clinic settings [26]. BioSpot-VIVAS operated throughout lessons while Coriolis Micro Air only operated shortly before and during break times to reduce noise exposure (approximately 60 min/day). The removable parts of both sampling devices were regularly autoclaved to avoid contamination. At the end of the day, samples were transported to the Institute for Infectious Diseases (IFIK) and stored at −80°C. At the end of the study period, the Xiaomi HEPA filters were carefully removed and 20 swabs were taken from each filter and stored at −80°C until further processed. 2.4 Laboratory and molecular analyses Prior to the real-time (RT)-PCR analysis, daily bioaerosol samples were combined for each sampling device and filtered using Amicon Ultra-15 Centrifugal Filters with Ultracel 10,000 Dalton molecular weight cutoffs filters (UFC9010; MilliporeSigma, Burlington, USA) to a volume of 1 mL. The human saliva samples were directly analyzed without prior filtration. The Allplex RV Master Assay (Seegene, Seoul, South Korea) was used to detect a panel of 19 respiratory viruses (Text A in S1 Appendix), including SARS-CoV-2. Viral genome load (VGL) of specimens was quantified using standardized dilution series and reported in genomic equivalent copies/L. For positive samples, we performed targeted sequencing of viral RNA to compare genetic relatedness between SARS-CoV-2 strains detected in the air and human samples [27]. However, we were unable to amplify and sequence any of the gene targets in the bioaerosol samples due to low RNA concentrations. 2.5 Statistical analyses and modeling Number of new cases of COVID-19 by date of symptom onset. The daily number of new cases of COVID-19 was inferred based on the number of students absent from school. Confirmed cases referred to absences due to a positive laboratory test result or isolation. Suspected cases referred to absences due to sickness probably related to COVID-19. Absences due to other known reasons were excluded. We used a probabilistic simulation approach (Text D in S1 Appendix) to estimate the number of suspected cases that were cases of COVID-19 and the date of symptom onset, which was not always reported by students. We generated 100 datasets for the daily number of new cases of COVID-19. Subsequent analyses were performed on each of these datasets and we report the mean of estimation results if not stated otherwise. Cases in teachers were excluded as teachers taught multiple classes and had different exposure. Aerosol and particle concentrations. We computed the average concentrations per day and compared these between study conditions. We report the mean ± standard deviation. We further estimated the reduction in concentrations using Bayesian log-linear regression models (Text H in S1 Appendix), adjusting for the ventilation rate (computed from indoor CO2 levels; see Text I in S1 Appendix), the daily number of students in class, school, and weekday effects. Risk of transmission. We estimated daily transmission of SARS-CoV-2 with a Bayesian semi-mechanistic hierarchical model [9] (Fig 3): (i) We modeled the number of new infections as a function of susceptible students in each class and day, where the probability of infection can vary by study condition. (ii) We adjusted the estimated effects of interventions for the daily proportion of all absent students and the effective reproduction number in the community. (iii) The number of new cases was computed as the weighted sum of the number of new infections in the previous days. (iv) The number of susceptible students was computed by removing the number of students who have already been infected. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Visual summary of the structure of the Bayesian semi-mechanistic hierarchical model. (1) The number of new infections was modeled as a function of the number of susceptible students, the class-specific daily rate of infections, and the reductions from active infection control measures; (2) the effects of control measures were adjusted for transmission in the community and the proportion of students in class; (3) the observed number of new cases was computed as the weighted sum of the number of new infections in the previous days; and (4) the number of susceptible students was computed as the total number of students minus the cumulative sum of infections in the previous days. https://doi.org/10.1371/journal.pmed.1004226.g003 A detailed description of the transmission model and choice of priors for all model parameters are provided in Text E in S1 Appendix. Model parameters were estimated with a Bayesian approach (Text F in S1 Appendix). Specifically, Markov chain Monte Carlo (MCMC) sampling was used as implemented by the Hamiltonian Monte Carlo algorithm with the No-U-Turn Sampler (NUTS) [28]. If not stated otherwise, we report posterior means and credible intervals (CrIs) based on the 50%, 80%, and 95% quantiles of the posterior samples, respectively. We estimated the total number of avoided infections for each intervention by computing the difference between the estimated number of infections in the presence and absence of interventions (counterfactual scenario). Software. All analyses were performed in R software (version 4.2.0) [29] and modeling in Stan (version 2.21.0) [30]. The code is available from https://github.com/nbanho/mcid. Number of new cases of COVID-19 by date of symptom onset. The daily number of new cases of COVID-19 was inferred based on the number of students absent from school. Confirmed cases referred to absences due to a positive laboratory test result or isolation. Suspected cases referred to absences due to sickness probably related to COVID-19. Absences due to other known reasons were excluded. We used a probabilistic simulation approach (Text D in S1 Appendix) to estimate the number of suspected cases that were cases of COVID-19 and the date of symptom onset, which was not always reported by students. We generated 100 datasets for the daily number of new cases of COVID-19. Subsequent analyses were performed on each of these datasets and we report the mean of estimation results if not stated otherwise. Cases in teachers were excluded as teachers taught multiple classes and had different exposure. Aerosol and particle concentrations. We computed the average concentrations per day and compared these between study conditions. We report the mean ± standard deviation. We further estimated the reduction in concentrations using Bayesian log-linear regression models (Text H in S1 Appendix), adjusting for the ventilation rate (computed from indoor CO2 levels; see Text I in S1 Appendix), the daily number of students in class, school, and weekday effects. Risk of transmission. We estimated daily transmission of SARS-CoV-2 with a Bayesian semi-mechanistic hierarchical model [9] (Fig 3): (i) We modeled the number of new infections as a function of susceptible students in each class and day, where the probability of infection can vary by study condition. (ii) We adjusted the estimated effects of interventions for the daily proportion of all absent students and the effective reproduction number in the community. (iii) The number of new cases was computed as the weighted sum of the number of new infections in the previous days. (iv) The number of susceptible students was computed by removing the number of students who have already been infected. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 3. Visual summary of the structure of the Bayesian semi-mechanistic hierarchical model. (1) The number of new infections was modeled as a function of the number of susceptible students, the class-specific daily rate of infections, and the reductions from active infection control measures; (2) the effects of control measures were adjusted for transmission in the community and the proportion of students in class; (3) the observed number of new cases was computed as the weighted sum of the number of new infections in the previous days; and (4) the number of susceptible students was computed as the total number of students minus the cumulative sum of infections in the previous days. https://doi.org/10.1371/journal.pmed.1004226.g003 A detailed description of the transmission model and choice of priors for all model parameters are provided in Text E in S1 Appendix. Model parameters were estimated with a Bayesian approach (Text F in S1 Appendix). Specifically, Markov chain Monte Carlo (MCMC) sampling was used as implemented by the Hamiltonian Monte Carlo algorithm with the No-U-Turn Sampler (NUTS) [28]. If not stated otherwise, we report posterior means and credible intervals (CrIs) based on the 50%, 80%, and 95% quantiles of the posterior samples, respectively. We estimated the total number of avoided infections for each intervention by computing the difference between the estimated number of infections in the presence and absence of interventions (counterfactual scenario). Software. All analyses were performed in R software (version 4.2.0) [29] and modeling in Stan (version 2.21.0) [30]. The code is available from https://github.com/nbanho/mcid. 2.6 Ethics statement The Ethics Committee of the Canton of Bern, Switzerland, approved the study (reference no. 2021–02377). For the saliva samples, we included all students who were willing to participate and obtained written informed consent from their caregivers. 3 Results The study population consisted of 90 students (39 female, 51 male; Table 1). Of these, 27 students were fully vaccinated and 34 students had recovered from an infection within the last year. Over the seven-week study period (3,150 student-days in total), students were absent from school for 644 days (20% of total) of which 147 days (23% of absences) were due to isolation related to COVID-19, and 247 (38% of absences) were due to sickness. Overall, there were 35 confirmed and 73 suspected cases of COVID-19, exceeding the number of students in some classes (Table C in S1 Appendix). This suggests that only a proportion of suspected cases were actual cases of COVID-19 as it is unlikely that students were infected twice. Accordingly, we estimated the actual number of cases of COVID-19 across schools to be 55 (95% CrI 35 to 76). Download: PPT PowerPoint slide PNG larger image TIFF original image Table 1. Overview of the study population, number of COVID-19 cases, and person-days of absences. https://doi.org/10.1371/journal.pmed.1004226.t001 3.1 Molecular analyses We analyzed 262 saliva, 130 bioaerosol samples and swabs from the filters of air cleaners (20 swabs per filter) in 2 classrooms. Overall, there were 21 positive saliva and 10 positive airborne samples. We detected SARS-CoV-2, adenovirus, and influenza virus (Table 2 and Tables A-B in S1 Appendix). SARS-CoV-2 made up the vast majority of positive saliva (19 out of 21) and bioaerosol samples (9 out of 10). We found 4 positive air-saliva samples in the same respective week (3 SARS-CoV-2 and 1 adenovirus), suggesting they were paired samples. We also detected SARS-CoV-2 and influenza viruses from the HEPA filters of the air cleaners. The number of positive saliva and airborne SARS-CoV-2 samples per week varied by study condition (Fig 4A). Without interventions, the proportion of positive samples per week was 11.5% for saliva and 8.1% for airborne samples. These proportions were lower with mask mandates (saliva: 5.7%, air: 7.1%) and air cleaners (saliva: 7.7%, air: 5.0%). The weekly average viral concentration of positive samples was 0.6 copies/L. There were also differences in viral concentration between study conditions (Fig 4B). Without interventions, it was 1.1 copies per liter per week, which was higher than with mask mandates (0.7 copies/L per week) and air cleaners (0.1 copies/L per week). Download: PPT PowerPoint slide PNG larger image TIFF original image Table 2. Analysis of molecular data and the number of positive and negative saliva and airborne samples in each school. https://doi.org/10.1371/journal.pmed.1004226.t002 Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Analysis of molecular data and comparison by study condition. (a) Proportion of positive airborne (solid line) and saliva (dashed line) samples (average per week). (b) Viral concentration in positive airborne samples from BioSpot-VIVAS (average per week). https://doi.org/10.1371/journal.pmed.1004226.g004 3.2 Analysis of aerosol and particle concentrations Particle concentrations differed by study condition (Fig 5A). Aerosol number concentrations were lower with mask mandates (mean 49 ± 52 1/cm3 standard deviation) and air cleaners (84 ± 56 1/cm3) than without interventions (177 ± 109 1/cm3). Similarly, particle mass concentrations (e.g., PM1) were lower with mask mandates (2.0 ± 1.6 μgm−3) and air cleaners (3.8 ± 2.9 μgm−3) than without interventions (6.9 ± 4.1 μgm−3). Overall daily average CO2 levels were 1,064 ± 232 ppm. CO2 levels without interventions (953 ± 198 ppm) were slightly lower than with mask mandates (1,155 ± 237 ppm) and air cleaners (1,088 ± 224 ppm), indicating increased ventilation through outdoor air exchange (Fig I in S1 Appendix), although the differences in the daily proportion of time that windows were opened between no intervention and mask mandates (0.03, 95% CrI −0.07 to 0.12) and between no intervention and air cleaners (0.00, 95% CrI −0.11 to 0.11) were indistinguishable from zero. When adjusting for different ventilation rates, the number of students in class, school effects, and weekday effects, the aerosol number concentration decreased by 69% (95% CrI 42% to 86%) with mask mandates and by 39% (95% CrI 4% to 69%) with air cleaners (Fig 5B and Table H in S1 Appendix). The concentration of smaller particles (PM1 and PM2.5) was more reduced during mask mandates, and the concentration of larger particles (PM4 and PM10) was more reduced during air cleaners. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 5. Analysis of particle concentrations and comparison by study condition. (a) Boxplot of the daily average values for aerosol number concentration (CN in 1/cm3) and particle mass concentration (PM for particles of sizes <1 to <10 μm, respectively in μgm−3). Results for CO2 and other environmental variables are provided in Text J in S1 Appendix. (b) Estimated reduction in aerosol number and particle mass concentrations with interventions (posterior mean as dot and 50%, 80%, and 95% CrI as lines, respectively). https://doi.org/10.1371/journal.pmed.1004226.g005 3.3 Estimating transmission risk of SARS-CoV-2 The cumulative number of cases of COVID-19 increased considerably in all classes without intervention and the large majority of students has been infected in School 1 by the time air cleaners were installed, except for class D in School 2 (Fig 6A). Based on our Bayesian transmission model, the probability of getting infected was lower with mask mandates than without interventions (adjusted odds ratio 0.19, 95% CrI 0.09 to 0.38), and comparable with air cleaners (1.00, 95% CrI 0.15 to 6.51). Excluding suspected cases from the model, these probabilities were similar for both mask mandates (0.21, 95% CrI 0.08 to 0.50) and air cleaners (0.98, 95% CrI 0.14 to 6.74), although with greater uncertainty. Considering both confirmed and suspected cases, mask mandates were associated with a considerable number of avoided infections (9.98, 95% CrI 2.16 to 19.00), but not air cleaners (Fig 6B). As an additional analysis, we used a modified Wells–Riley equation [31] and assumed that the change in the rate of emitted infectious quanta was proportional to the estimated reduction in the aerosol number concentration, while other parameters were kept constant under all study conditions (Text K in S1 Appendix). Accordingly, the daily risk of infection for a 6 h school day was 1.0% (95% CrI 0.4% to 1.9%) with mask mandates and 1.9% (95% CrI 1.0% to 3.0%) with air cleaners, compared to a 3.1% risk without interventions (Fig J in S1 Appendix). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 6. Analysis of epidemiological data and estimated transmission by study condition using the Bayesian hierarchical transmission model. (a) Estimated mean cumulative number of cases for each school class after probabilistic simulation accounting for uncertainty in suspected cases and the date of symptom onset (Text D in S1 Appendix). Dotted red lines indicate the number of students per class. A comparison of the estimated number of new cases after probabilistic simulation with the observed number of new confirmed and suspected cases is shown in Fig C in S1 Appendix. (b) Estimated number of avoided infections associated with interventions across schools (posterior mean as dot and 50%, 80%, and 95% CrI as shaded areas) based on the Bayesian hierarchical transmission model. Detailed estimation results are provided in Text G in S1 Appendix. https://doi.org/10.1371/journal.pmed.1004226.g006 3.1 Molecular analyses We analyzed 262 saliva, 130 bioaerosol samples and swabs from the filters of air cleaners (20 swabs per filter) in 2 classrooms. Overall, there were 21 positive saliva and 10 positive airborne samples. We detected SARS-CoV-2, adenovirus, and influenza virus (Table 2 and Tables A-B in S1 Appendix). SARS-CoV-2 made up the vast majority of positive saliva (19 out of 21) and bioaerosol samples (9 out of 10). We found 4 positive air-saliva samples in the same respective week (3 SARS-CoV-2 and 1 adenovirus), suggesting they were paired samples. We also detected SARS-CoV-2 and influenza viruses from the HEPA filters of the air cleaners. The number of positive saliva and airborne SARS-CoV-2 samples per week varied by study condition (Fig 4A). Without interventions, the proportion of positive samples per week was 11.5% for saliva and 8.1% for airborne samples. These proportions were lower with mask mandates (saliva: 5.7%, air: 7.1%) and air cleaners (saliva: 7.7%, air: 5.0%). The weekly average viral concentration of positive samples was 0.6 copies/L. There were also differences in viral concentration between study conditions (Fig 4B). Without interventions, it was 1.1 copies per liter per week, which was higher than with mask mandates (0.7 copies/L per week) and air cleaners (0.1 copies/L per week). Download: PPT PowerPoint slide PNG larger image TIFF original image Table 2. Analysis of molecular data and the number of positive and negative saliva and airborne samples in each school. https://doi.org/10.1371/journal.pmed.1004226.t002 Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 4. Analysis of molecular data and comparison by study condition. (a) Proportion of positive airborne (solid line) and saliva (dashed line) samples (average per week). (b) Viral concentration in positive airborne samples from BioSpot-VIVAS (average per week). https://doi.org/10.1371/journal.pmed.1004226.g004 3.2 Analysis of aerosol and particle concentrations Particle concentrations differed by study condition (Fig 5A). Aerosol number concentrations were lower with mask mandates (mean 49 ± 52 1/cm3 standard deviation) and air cleaners (84 ± 56 1/cm3) than without interventions (177 ± 109 1/cm3). Similarly, particle mass concentrations (e.g., PM1) were lower with mask mandates (2.0 ± 1.6 μgm−3) and air cleaners (3.8 ± 2.9 μgm−3) than without interventions (6.9 ± 4.1 μgm−3). Overall daily average CO2 levels were 1,064 ± 232 ppm. CO2 levels without interventions (953 ± 198 ppm) were slightly lower than with mask mandates (1,155 ± 237 ppm) and air cleaners (1,088 ± 224 ppm), indicating increased ventilation through outdoor air exchange (Fig I in S1 Appendix), although the differences in the daily proportion of time that windows were opened between no intervention and mask mandates (0.03, 95% CrI −0.07 to 0.12) and between no intervention and air cleaners (0.00, 95% CrI −0.11 to 0.11) were indistinguishable from zero. When adjusting for different ventilation rates, the number of students in class, school effects, and weekday effects, the aerosol number concentration decreased by 69% (95% CrI 42% to 86%) with mask mandates and by 39% (95% CrI 4% to 69%) with air cleaners (Fig 5B and Table H in S1 Appendix). The concentration of smaller particles (PM1 and PM2.5) was more reduced during mask mandates, and the concentration of larger particles (PM4 and PM10) was more reduced during air cleaners. Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 5. Analysis of particle concentrations and comparison by study condition. (a) Boxplot of the daily average values for aerosol number concentration (CN in 1/cm3) and particle mass concentration (PM for particles of sizes <1 to <10 μm, respectively in μgm−3). Results for CO2 and other environmental variables are provided in Text J in S1 Appendix. (b) Estimated reduction in aerosol number and particle mass concentrations with interventions (posterior mean as dot and 50%, 80%, and 95% CrI as lines, respectively). https://doi.org/10.1371/journal.pmed.1004226.g005 3.3 Estimating transmission risk of SARS-CoV-2 The cumulative number of cases of COVID-19 increased considerably in all classes without intervention and the large majority of students has been infected in School 1 by the time air cleaners were installed, except for class D in School 2 (Fig 6A). Based on our Bayesian transmission model, the probability of getting infected was lower with mask mandates than without interventions (adjusted odds ratio 0.19, 95% CrI 0.09 to 0.38), and comparable with air cleaners (1.00, 95% CrI 0.15 to 6.51). Excluding suspected cases from the model, these probabilities were similar for both mask mandates (0.21, 95% CrI 0.08 to 0.50) and air cleaners (0.98, 95% CrI 0.14 to 6.74), although with greater uncertainty. Considering both confirmed and suspected cases, mask mandates were associated with a considerable number of avoided infections (9.98, 95% CrI 2.16 to 19.00), but not air cleaners (Fig 6B). As an additional analysis, we used a modified Wells–Riley equation [31] and assumed that the change in the rate of emitted infectious quanta was proportional to the estimated reduction in the aerosol number concentration, while other parameters were kept constant under all study conditions (Text K in S1 Appendix). Accordingly, the daily risk of infection for a 6 h school day was 1.0% (95% CrI 0.4% to 1.9%) with mask mandates and 1.9% (95% CrI 1.0% to 3.0%) with air cleaners, compared to a 3.1% risk without interventions (Fig J in S1 Appendix). Download: PPT PowerPoint slide PNG larger image TIFF original image Fig 6. Analysis of epidemiological data and estimated transmission by study condition using the Bayesian hierarchical transmission model. (a) Estimated mean cumulative number of cases for each school class after probabilistic simulation accounting for uncertainty in suspected cases and the date of symptom onset (Text D in S1 Appendix). Dotted red lines indicate the number of students per class. A comparison of the estimated number of new cases after probabilistic simulation with the observed number of new confirmed and suspected cases is shown in Fig C in S1 Appendix. (b) Estimated number of avoided infections associated with interventions across schools (posterior mean as dot and 50%, 80%, and 95% CrI as shaded areas) based on the Bayesian hierarchical transmission model. Detailed estimation results are provided in Text G in S1 Appendix. https://doi.org/10.1371/journal.pmed.1004226.g006 4 Discussion We collected epidemiological, environmental, and molecular data to estimate transmission of respiratory infections in schools and assessed the association with infection control measures. Airborne detection of SARS-CoV-2 documented sustained SARS-CoV-2 transmission. Bioaerosol SARS-CoV-2 concentrations with general mask mandates were lower, and SARS-CoV-2 was detected on the filters of air cleaners. Both interventions were associated with significantly lower aerosol number and particle mass concentrations. The Bayesian transmission model using epidemiological data estimated that mask mandates avoided SARS-CoV-2 infections, but not air cleaners. Our study demonstrated airborne detection of SARS-CoV-2 in schools. Although sampling and molecular detection of infectious bioaerosols are challenging and there is no agreed standard [16,32,33], it provides important evidence on the airborne transmission of respiratory pathogens. So far, viral RNA in airborne samples of SARS-CoV-2 was mainly found in hospitals and healthcare facilities [34]. A related study in 2 South African schools detected airborne Mycobacterium tuberculosis [35]. Tuberculosis was the leading cause of death worldwide due to an infectious disease prior to the COVID-19 pandemic. An improved understanding of airborne transmission and the effectiveness of interventions can benefit the control of both infectious diseases [36]. Our study provides evidence on airborne SARS-CoV-2 viral transmission and the potential effects of interventions in schools based on airborne and saliva samples from students. In our study, positive samples mostly pertained to SARS-CoV-2, but we occasionally detected other respiratory viruses such as adenovirus and influenza. General public health measures during the study likely suppressed the spread of other respiratory viruses. It also must be noted that the detected viral concentrations were low, as shown by high cycle thresholds (CTs) in the RT-qPCR results. The molecular presence of other viral pathogens cannot be excluded. The low sensitivity for detection of airborne pathogens by molecular assays is a well-documented challenge [16,33]. An experimental study demonstrated the effectiveness of infection control measures (universal face mask wearing and air cleaners) using simulated exhaled SARS-CoV-2 bioaerosols in a closed indoor space [17]. In contrast, we studied infection control measures in a real-life setting and demonstrated their effectiveness using a multiple-measurement approach and obtain similar results. Our findings also align with existing evidence from population-level studies showing that the incidence of COVID-19 was lower in schools with mask use [13,15] and improved ventilation [15]. Similarly, a field study showed that adequate ventilation was associated with reduced incidence of tuberculosis, a strictly airborne disease, in a university building [37]. Altogether, these findings support arguments in favor of multiple intervention strategies to address airborne transmission of respiratory infections in crowded indoor settings [38]. Indoor ventilation is one of the key determinants of airborne transmission [2,3], but schools are often poorly ventilated [24,39]. We showed that the concentration of both larger and smaller particles were lower with air cleaners, in line with findings about their effectiveness in hospitals [16] and simulated indoor environments [17]. The detection of viruses on the filters of air cleaners further supports the evidence on the effective removal of bioaerosols. However, it was difficult to estimate changes in transmission following air cleaners because they were installed at the end of the study period when presumably a large proportion of students were already infected with SARS-CoV-2. Furthermore, air cleaners could not prevent transmission outside the classrooms (e.g., during breaks and outside classes), which limits their effectiveness compared to masks that had to be worn in all indoor settings following the general mask mandate. In addition, physicochemical properties of aerosols, environmental factors, and the distance to infectious people determine the survival of airborne particles and the loss of infectivity over time [3]. Thus, a predominance of close range high particle density aerosol transmission of SARS-CoV-2 in school settings could further explain why air cleaners were not associated with reduced transmission. Our study used portable, affordable air cleaners that could be implemented at scale, rather than large cleaners [25]. Noise exposure and a lack of acceptance by teachers [40] may still prevent their widespread use. Although not specifically assessed, we perceived good acceptance of our air cleaners during the short study period. Nevertheless, investments in professional building ventilation systems should be preferred to air cleaners in the long term [41]. School closures during the COVID-19 pandemic have been intensely debated as children and adolescents are particularly vulnerable to the negative impact of such interventions on their well-being and mental health [11]. Furthermore, numerous studies have examined the role of children in transmitting SARS-CoV-2 [12] and it remains unclear to what extent transmission of SARS-CoV-2 occurs in schools [42]. These findings contrast with studies of influenza viruses showing that school children may drive the seasonal influenza epidemic. Community studies in the US demonstrated that influenza transmission rates in children and adolescents were high in schools and that they easily transmit influenza viruses to household members and into their communities [43–45]. Our study suggests that also transmission of SARS-CoV-2 occurs to a considerable extent in schools. Our study has several limitations. First, aerosol measurements and molecular detection of pathogens in bioaerosol samples document exposure, but not transmission and the direction of transmission (human to air, air to human). Nevertheless, paired saliva and airborne samples suggest that infected students exhaled infectious particles into the air, indicating a considerable transmission risk in the school rooms. Second, a comparison of viral concentration between study conditions should be interpreted with care due to the technical limitations of molecular detection in bioaerosol samples, and because the number of possibly infectious (and thus susceptible) students decreased towards the end of the study. Third, our epidemiological analysis is based on observational data, thus subject to potential confounding, e.g., the incidence of COVID-19 (cases per week) in the community varied over the study period. However, levels were high throughout and included 2 Omicron subvariant peaks. Community transmission was also considered in our Bayesian transmission model. CO2 levels were not considered in the model but the levels were slightly lower without interventions, suggesting that lower ventilation during intervention phases may have actually reduced the estimated effectiveness of mask mandates and air cleaners. Fourth, epidemiological data may not always be complete due to the underreporting of COVID-19 among absent students. We thus used a probabilistic approach to estimate the proportion of suspected cases being actual cases of COVID-19 and the date of symptom onset where it was not reported. While this allows us to consider uncertainty in the observed data, the estimated effects of interventions will be less precise as reflected in larger uncertainty intervals. Fifth, the ordering of study conditions was the same across classes, mainly because our study period coincided with the lifting of general mask mandates. The effectiveness of infection control measures may thus be affected by their timing. Future studies could vary the ordering of interventions in each class using a cross-over design. This would allow exploiting variation in the data between classes and reduce the influence of timing. Finally, our study design only allowed us to analyze variation within classes over time. We, therefore, did not analyze variation between classes, although we observed some differences such as lower CO2 levels and transmission in School 2. These differences may be explained by school-specific factors not measured in our study. In conclusion, using epidemiological, environmental, and molecular data, our study suggests that considerable transmission of SARS-CoV-2 occurred in the participating schools. General face mask wearing was associated with reduced SARS-CoV-2 transmission and prevented infections. The effectiveness of interventions was supported by significant decreases in the concentration of aerosols. Taken together, our results suggest that infection control measures can reduce the transmission of respiratory infections in school rooms. Future studies may use our multiple-measurement approach to assess the effectiveness of infection control measures in reducing the transmission of respiratory infections. Ideally, these data should be collected routinely in sentinel schools, thus continuously monitoring transmission risks and alerting health authorities when infection control measures should be taken. Supporting information S1 STROBE Checklist. STROBE Checklist. https://doi.org/10.1371/journal.pmed.1004226.s001 (PDF) S1 Appendix. Supplementary information. Text A. Details on laboratory and molecular analyses. Text B. Summary of case and molecular data. Text C. Longitudinal case data. Text D. Probabilistic simulation of case data. Text E. Estimating transmission and the effects of infection control measures. Text F. Model parameter estimation. Text G. Detailed results from transmission model. Text H. Estimating changes in particle concentrations. Text I. Computing rebreathed air volume and ventilation rate. Text J. Results for changes in environmental variables. Text K. Modeling transmission risk of SARS-CoV-2 using a modified Wells–Riley equation. Fig A. Proportion of suspected cases being actual cases of COVID-19. Fig B. Empirical and fitted distribution for the delay from symptom onset to absence. Fig C. Comparison of reported and estimates cases of COVID-19. Fig D. Prior for the probability of getting infected without interventions. Fig E. Choices of prior for incubation period. Fig F Estimated incidence over time. Fig G. Model- and simulation-based estimates of the number of COVID-19 cases. Fig H. Estimated number of avoided infections with interventions. Fig I. Boxplot of environmental variables by school and study condition. Fig J. Estimated transmission risk using a modified Wells–Riley equation. Table A. Reported cases of COVID-19, saliva, and airborne samples in School 1. Table B. Reported cases of COVID-19, saliva, and airborne samples in School 2. Table C. Overview of the study population, number of COVID-19 cases, and person-days of absences in each study class. Table D. List of confirmed and suspected cases over the study period in School 1. Table E. List of confirmed and suspected cases over the study period in School 2. Table F. Prior choices for model parameters. Table G. Estimation results from transmission model. Table H. Estimated reduction in aerosol and particle concentrations with interventions. https://doi.org/10.1371/journal.pmed.1004226.s002 (PDF) Acknowledgments We would like to thank the schools, teachers, and students participating in the study. We are also grateful to the Educational Department of the Canton of Solothurn for their support throughout the study. Finally, we are indebted to the student assistants who helped with the data collection at the schools.