Single-cell transcriptomic analysis reveals the immune landscape of lung in steroid-resistant asthma exacerbation
Edited by Ariel Munitz, Tel Aviv University, Tel Aviv, Israel, and accepted by Editorial Board Member Ruslan Medzhitov November 27, 2020 (received for review April 15, 2020)
Significance
Asthma exacerbation is not prevented by standard corticosteroid-based therapy and is the major burden in terms of morbidity, mortality, and health care costs associated with asthma. Etiology of disease exacerbation is highly heterogeneous. Using single-cell RNA deep sequencing, we undertook a large-scale, high-dimensional analysis of the immune landscape in lung following dexamethasone treatments in a mouse model of asthma exacerbation. We found multicellular signaling pathways closely associated with disease progression during the exacerbation. IL-13 produced by CD8+ memory T cells, ILC2, and basophils is key a cytokine in driving the pathogenesis of asthma exacerbation in this model. These data provide important insights into how the immune landscape in lung during asthma exacerbation can shape the progress of disease.
Abstract
Exaggerated airway hyperresponsiveness and inflammation are hallmarks of asthma, and lipopolysaccharide (LPS) exposure is linked to the severity of the disease and steroid resistance. To investigate the mechanisms underlying asthma exacerbation, we established a mouse model of LPS-induced steroid-resistant exacerbation on the background of house dust mite (HDM)-induced asthma to profile the immune cells in lung by using single-cell RNA deep sequencing. Twenty immune subsets were identified by their molecular and functional properties. Specific cell clusters of basophils, type 2 innate lymphoid cells (ILC2), and CD8+ memory T cells were the predominant sources of interleukin (IL)-4 and IL-13 transcripts whose expressions were dexamethasone resistant. Production of IL-13 by these cells was validated by IL-13-reporter mice. Neutralization of IL-13 abolished HDM/LPS-induced airway hyperresponsiveness, airway inflammation, and decreased mucus hypersecretion. Furthermore, using Ingenuity Pathway Analysis systems, we identified canonical pathways and upstream regulators that regulate the activation of basophils, ILC2, and CD8+ memory T cells. Our study provides mechanistic insights and an important reference resource for further understanding of the immune landscape during asthma exacerbation.
Sign up for PNAS alerts.
Get alerts for new articles, or get an alert when an article is cited.
Asthma exacerbations are the major cause of hospitalization of asthma patients, with a substantial economic burden on healthcare systems (1). Factors that trigger asthma exacerbation include house dust mite (HDM), pollen, respiratory infection, cigarette smoke, and pollutants (2). Accumulative evidence suggests that the interplay between the innate immune response and the underlying type 2 immunity microenvironment in lung alters inflammatory sequelae and triggers exacerbation (3). Type 2 cytokine-regulated airway inflammation and clinical hallmarks of asthma that are induced by allergen exposure can be effectively controlled by mainstay therapy with glucocorticoids and β2 agonists (4). However, it is increasingly clear that asthma exacerbation is strongly associated with the aberrant activation of immune responses in lung and is refractory to glucocorticoid treatment (5, 6). Such activation may be of particular relevance in the respiratory infection- and allergen-triggered acute exacerbations of asthma and severe asthma (7).
HDM allergens are among the leading causative factors for the induction of allergic lung inflammation in patients with asthma (8). Proteases derived from HDM contribute to the pathogenesis of asthma by prompting aberrant innate and adaptive immune responses (9). Lipopolysaccharide (LPS) is commonly derived from the membrane of colonizing Gram-negative bacteria and environmental contamination (10), and the level of LPS directly correlates with the severity of asthma and decline in lung function (11). Furthermore, clinical studies and experimental models of allergic asthma and asthma exacerbation indicate that LPS, together with allergens, may play important roles in disease exacerbation (12–15), skewing a glucocorticoid-responsive phenotype to a glucocorticoid-nonresponsive phenotype (16–18).
Episodes of asthma exacerbations are characterized by exaggerated inflammation with infiltration of eosinophils, macrophages/monocytes, neutrophils, and T lymphocytes to the airways (19). Of note, neutrophil recruitment is one of the prominent features of acute exacerbations in patients who poorly respond to glucocorticoid treatments (20). The infiltration of neutrophils into the airways indicates that the innate host defense pathway against infection has been activated (21). Although there are several lines of evidence that respiratory infection and allergen exposure contribute to the development of asthma exacerbations, the intricate molecular network between neutrophil infiltration, glucocorticoid resistance, and pathogenesis of airway obstruction during asthma exacerbations remains unclear. In recent years, immunotherapies targeting asthma have substantially improved the treatment strategies (22); however, their efficacy is not concordant among asthma patients or subtypes. In this regard, it is important to determine specific biomarkers and molecular networks not only to predict therapy outcomes but also to unravel the intricate reciprocation between the host immune cells.
Single-cell RNA deep sequencing (scRNA-seq) has been recently applied to understand the development and interplay of heterogeneous cells and examine the expression of vast genes and transcript isoforms at a genome-wide scale (23). This allows us to inspect the highly complex tissue microenvironment of chronic disease with unprecedented detailed information. To shed light on the complexity of infiltrating immune cells in lung, we established a mouse model of asthma exacerbation induced by HDM and LPS to understand the pathogenesis. We then undertook a large-scale, high-dimensional analysis of immune cells isolated from lung and determined the landscape of major clusters of immune cells in lung following dexamethasone treatments with scRNA-seq. We identified 20 major immune cell subsets with distinct patterns. Signature genes for these cells were examined in detail, and we found that expression of interleukin (IL)-4 and IL-13 by basophils, group 2 innate lymphoid cells (ILC2), and CD8+ memory cells is largely steroid resistant. This unprecedented detailed resource is essential for studying the characteristics of adaptive and innate immune cells in lung and for developing effective therapy strategies for patients with asthma exacerbations.
Materials and Methods
Mice.
BALB/c mice (6 to 8 wk) were obtained from the specific-pathogen–free (SPF) facilities of the University of Newcastle and Zhengzhou University. IL-13-tdTomato-reporter mice were provided by Andrew McKenzie, Medical Research Council, Cambridge, UK. Experiments were approved by the animal ethics committees of the University of Newcastle (#A-2017-721) and by Zhengzhou University (#ZZURIB20180120). Mice were housed in approved SPF containment facilities.
Induction of Allergic Airways Disease and Exacerbation.
Mice were sensitized and challenged by intranasal (i.n.) exposure to HDM extract (50 μg crude Dermatophagoides pteronyssinus extract in 50 μL sterile saline; Greer Laboratories) daily on days 1, 2, and 3 for sensitization, followed by i.n. rechallenges (5 μg in 50 μL sterile saline) daily on days 14, 15, 16, and 17. Nonsensitized mice received sterile saline only. Where indicated, mice received i.n. LPS (50 ng; Sigma-Aldrich) in phosphate-buffered saline (PBS) on days 20 and 22. Endpoints were assessed on day 25. To assess responsiveness to corticosteroid treatment, mice were treated with dexamethasone (1 mg/kg intraperitoneal [i.p.]; Sigma-Aldrich) or PBS (as vehicle) on days 20 and 22. To neutralize IL‐13, mice were injected with anti–IL‐13 antibody (150 μg, i.p, clone eBio1316H; eBioscience) or isotype control (rat IgG1K, clone eBRG1, eBioscience) on days 19 and 21.
Lung Function.
Airway resistance (Raw) in response to methacholine (Sigma-Aldrich) was measured using a Flexivent apparatus (FX1 system; Scireq) (24). Briefly, mice were anesthetized by injection of xylazine (2 mg/mL i.p.; Troy Laboratories), ketamine (40 mg/ml; Parnell), and PBS (1:4:5). A cannula was inserted into the trachea, and mice were ventilated with a tidal volume of 8 ml/kg at a frequency of 450 breaths/minute. Mice were challenged with saline aerosol, followed by increasing concentrations of methacholine (0.3, 1, 3, 10, and 30 mg/mL). Measurements were excluded if the coefficient of determination was lower than 95%. Airway resistance was recorded and presented as the percentage increase over baseline.
Bronchoalveolar Lavage Fluid Differential Cell Counts.
Bronchoalveolar lavage fluid (BALF) was collected immediately following lung function measurements, as previously described (25). The left lobe of the lung was tied off, and the right lung was flushed twice with 700 µl Hanks Balanced Salt Solution (Invitrogen). BALF samples were centrifuged to yield a pellet. Red blood cells were lysed using 0.86% (wt/vol) ammonium chloride. Total cell counts were then determined by hemocytometer, and the remaining cells were cytospun onto glass slides. Differential leukocyte counts were determined using morphological criteria and light microscopy (×100) on May-Grünwald and Giemsa-stained slides, counting 300 cells/slide/sample.
Lung Immune Cell Purification.
Lung immune cells were prepared and pooled for scRNA-seq from the lungs of six mice of each group. Mice were euthanized by i.p. injection of an overdose of pentobarbital. After thoracotomy, the pulmonary circulation was perfused with 37 °C PBS to remove intravascular cells. Lung tissues were mechanically minced using fine scissors into RPMI1640 containing digestion enzymes (3 mg/mL Collagenase IV; Worthington Biochemical) and 40 mg/mL DNase I (Sigma-Aldrich) and incubated for 30 min at 37 °C and 5% CO2 with gentle shaking every 5 to 10 min. Samples were depleted of erythrocytes using 0.86% (wt/vol) ammonium chloride. After washing, a single-cell solution with 1 × 107 cells/ml was prepared by filtering through a 100-μm filter. Cells were then washed with ice-cold 1% fetal calf serum (FCS)/PBS. Fc receptor blocker (anti-mouse CD16/32, 1:100) was added, followed by anti-mouse CD45 antibody (1:300) and washing with ice-cold 1% FCS/PBS twice. Cells were stained with Fixable Viability Dye 7-AAD (1:40) for exclusion of dead cells. Lung immune cells were purified with a BD AriaIII with a high degree of purity (routine over 95%).
In some experiments, single-cell suspension of lung cells of IL-13-tdTomato-reporter mice were then washed in PBS and resuspended in Zombie Yellow Fixable Viability Kit (BioLegend) for exclusion of dead cells. Fc receptor block (anti-CD16/32) was added followed by anti‐mouse CD3e, CD4, CD8a, CD45, CD90.2, ICOS, lineage marker, CD49b, FcεRIα, cKit, and ST2 antibodies (BD Pharmingen, BioLegend) and washed with fluorescence-activated cell sorting (FACS) buffer (PBS +1% FCS +2 mmol/L ethylenediaminetetraacetic acid) twice. Cells were then incubated on ice in 4% paraformaldehyde, permeabilized using BD Perm/Wash buffer (BD Pharmingen), and stained for intracellular marker (anti–IL‐13). Analysis was performed using Flowjo (BD Pharmingen). IL‐13–producing ILC2 cells were defined as CD45+Lin−ICOS+CD90.2+ST2+IL‐13+, IL‐13–producing CD8+ T cells as CD3+CD8a+IL‐13+, IL‐13–producing CD4+ T cells as CD3+CD4+IL13+, and IL‐13–producing basophils as CD49b+CD200R3+ FcεRIα+IgE+cKit−IL‐13+.
Droplet-Based scRNA-Seq.
Immediately post sorting, CD45+7-AAD− single cells were run on a 10× Chromium system (10 × Genomics) (26) and then through a library preparation by LC Sciences, following the recommended protocol for the Chromium Single Cell 30 Reagent Kit (v2 Chemistry). Libraries were run on the HiSeq4000 for Illumina sequencing. Postprocessing and quality control were performed using a 10× Cell Ranger package (v1.2.0; 10 × Genomics). Reads were aligned to the mm10 reference assembly (v1.2.0; 10 × Genomics). Primary assessment with the 10× Cell Ranger for the sample treated with saline reported 4,502 cell-barcodes with 1,105 median genes per cell sequenced to 93.4% sequencing saturation with 91,009 mean reads per cell. Primary assessment with this software for the sample treated with HDM/LPS + vehicle reported 4,034 cell-barcodes with 1,018 median genes per cell sequenced to 93.5% sequencing saturation with 129,637 mean reads per cell. Primary assessment with this software for the sample treated with HDM/LPS + dexamethasone reported 5,823 cell-barcodes with 1,098 median genes per cell sequenced to 93.7% sequencing saturation with 85,703 mean reads per cell.
Bioinformatic Analysis of scRNA-Seq Data.
Gene expression matrices generated by the 10× Cell Ranger aggregate option were further analyzed with R package Seurat (version 3.0) with default parameters. The data were filtered according to the following thresholds: less than 200 or greater than 2,500 as unique expressed genes (nFeature_RNA) and greater than 5% as the percentage of mitochondrial genome content. The data were then normalized by converting with a scale factor (default as 10,000) and log-transformed with the Seurat embedded function. A correlation analysis was performed by employing the RunPCA function of the Seurat package, followed by an integrated analysis of the three datasets. Clustering analysis was carried out with standard Seurat package procedures with a resolution at 1.2. The identified clusters were then visualized using t-distributed Stochastic Neighbor Embedding (tSNE) of the principal components in Seurat.
Average gene expression matrixes were then retrieved for each cluster, and differential expression among clusters was performed to identify the top markers at a high level by each cluster with the FindAllMarkers implemented function (parameters: only.pos = FALSE, min.pct = 0.2, thresh.use = 0.2).
Cell clusters were annotated with SingleR for unbiased cell-type recognition of scRNA-seq and compared with the Immunological Genome Project (reference mouse dataset) (27). Spearman correlation analysis was performed on variable genes after comparing with each sample in the reference dataset. The coefficient of multiple correlation (per cell type) was aggregated to provide a single value per cell type per single cell. All analyzed features were binned on the basis of averaged expression, and the control features were randomly selected from each bin. For some tSNE plots, feature plots were demonstrated for selected marker genes using the Seurat function FeaturePlot.
Pathway enrichment analysis examining enriched processes in clusters was performed using Ingenuity Pathway Analysis (IPA) (28). Gene ontology (GO) enrichment analyses were performed with topGO package in R (Bioconductor) (29). To uncover the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway that is potentially linked to the identified steroid-resistant genes in clusters 9, 11, 13, 15, and 16, pathway enrichment analysis was performed by ClueGO plugin of Cytoscape software (30). A right-side hypergeometric test was employed for calculation of the P value, followed by a Benjamini–Hochberg adjustment. A pathway with the adjusted P value of <0.05 was considered significant. The ShinyGO v0.41 tool was used to determine potential promoters for those identified steroid-resistant genes in each cluster (31).
Statistical Analysis.
Statistical analysis was performed using Prism version 8.3 (GraphPad Software). Two-way ANOVA was used to identify differences between two or more experimental groups, and Student’s unpaired t tests were used where comparisons were made between two treatment groups. All results are presented as mean ± SEM, and P values of <0.05 were considered statistically significant.
Results
LPS Induces Steroid-Resistant Airway Hyperresponsiveness and Airway Inflammation in a Mouse Model of Asthma Exacerbation.
To investigate the mechanisms of steroid resistance during asthma exacerbation, we first treated mice on days 20 and 22 with LPS i.n. challenge in the presence of HDM-induced allergic airways disease (AAD) (Fig. 1A). HDM-induced AAD features in lung were significantly suppressed by dexamethasone treatment (Fig. 1B). Three days after the final LPS challenge or dexamethasone treatment, lung function and inflammatory responses were examined. First, we demonstrated that the HDM/LPS group displayed significantly exacerbated and prolonged airway hyperresponsiveness (AHR), elevated infiltration of macrophages, neutrophils, mucus hypersecretion, and inflammation compared with the HDM or LPS alone group (SI Appendix, Fig. S1). In addition, LPS instillation did not affect the levels of HDM-induced eosinophil and lymphocyte infiltration in the HDM/LPS group. Importantly, dexamethasone administration did not suppress HDM/LPS-exaggerated AHR (Fig. 1C); elevated infiltration of eosinophils, lymphocytes, macrophages, and neutrophils in BALF (Fig. 1C); nor eosinophils and inflammatory scores, but rather only reduced a small amount of mucus-producing cells in the airways (Fig. 1D) and heightened pathological changes in lung (Fig. 1E). Saline challenge had no impact on the above-mentioned asthma features. Dataset S1 summarizes the samples that were analyzed.
Fig. 1.

scRNA-Seq in Lung Immune Cells Identifies Multiple Clusters of Immune Cells.
We performed scRNA-seq using a 10× Genomics platform (32) in pooled immune cells from six lungs from each individual group of three groups including saline (SAL)-, HDM/LPS+vehicle (VEH)-, or HDM/LPS+dexamethasone (DEX)-treated mice, as this tissue compartment contains abundant activated immune cells (33). Single CD45+ cells from the lungs of mice were purified by FACS with an AriaIII system. The majority of cells passed quality control. Approximately 1,000 genes were detected per cell, and over 23,000 different genes were detected overall, with 3,139 highly variable genes (HVGs) detected. Outliers expressing less than 200 or greater than 2,500 of nFeature_RNA were omitted, and doublets were removed with Doublet-Finder (34). Principal component analysis was performed on cells after normalization for the gene numbers detected in each group. Total RNA molecules and unique genes detected in each cell, and RNA molecules per gene within one cell, were analyzed, showing that some cell clusters may actively express a higher amount of RNA transcripts than observed with other clusters (SI Appendix, Fig. S2).
By analyzing the CD45+ immune cells in mouse lung tissue, we identified 102 distinct clusters with specific molecular markers by using SingleR (SI Appendix, Fig. S3). Heatmap gene expression of these cells was performed (SI Appendix, Fig. S4). By examining these cells in detail, 20 major clusters of cells were defined as known immune cell types, including four subpopulations of monocytes, ILC2, Tregs, and basophils (Fig. 2A). The expression of classical markers for lymphocytes (CD3, CD4, CD8a, and CD19), monocytes and macrophages (CD14, CD68, and Ly6C), and NK cells (NKG7) were shown in corresponding cell clusters (Fig. 2B). We also generated a gene expression heatmap and identified the top 20 marker genes for each of those clusters (SI Appendix, Fig. S5A and Dataset S2). Weighted gene coexpression network analysis was further used to evaluate those clusters (35). The resulting networks grouped differentially expressed transcripts into six coexpressed modules (SI Appendix, Fig. S6 and Dataset S3), each showing a specific phenotype for T and B lymphocyte, NK cell, CD11b− macrophage, monocyte, and neutrophil clusters by unsupervised hierarchical clustering.
Fig. 2.

By using topGO, enriched GO terms were those where P < 0.05 in the topGO Fisher’s exact test for the 3,136 HVGs, and the enrichment factors of the top 20 GO terms are shown (SI Appendix, Fig. S5B). The top GO term for the 3,136 HVGs was associated with carbohydrate binding, components of membrane, regulation of protein processing and maturation, regulation of IFNγ and complement function, and regulation of membrane functions were also identified as GO terms of these genes (SI Appendix, Fig. S5B and Dataset S4).
CD11b+ Macrophages, Multiple Dendritic Cells (DCs) Subpopulations, and ILC2 May Be Associated with Steroid-Resistant AHR and Airway Inflammation.
A set of genes, consisting of 1,769 transcripts, were up-regulated more than two-fold when using the adjusted false discovery rate P < 0.05 cutoff for significance after integration of the scRNA-seq data of three treatments (Fig. 3A). The molecular landscapes in the SAL, VEH, and DEX groups were different, with 848, 1,335, and 1,379 HVGs detected as marker genes in the SAL, VEH, and DEX groups, respectively (Fig. 3A and Dataset S5). LPS treatment induced a slight decrease in the numbers of lymphocyte clusters while inducing a slight increase in macrophages and CD103− DC clusters in the VEH group (Fig. 3 B and C). Interestingly, dexamethasone treatment augmented the levels of many clusters, including CD4+ T cells, CD19+ B cells, neutrophils, three monocyte clusters, and Tregs, but decreased CD11b− macrophages. The levels of other cells were not affected by dexamethasone treatment.
Fig. 3.

Levels of CD11b+ macrophages increased 72 h after LPS treatment in the VEH group and were steroid resistant. There were 382 genes expressed in both VEH and DEX groups (SI Appendix, Fig. S7A). ClueGO plugin of Cytoscape were used to study and visualize the functional enrichment of those 382 genes (SI Appendix, Fig. S7B). The complex networks were predominantly enriched in KEGG pathways including FcεRI signaling, chemokine signaling, AGE-RAGE signaling, C-type lectin receptor signaling, and NF-κB signaling pathways. By using ShinyGO that is a web-based bioinformatic tool, the top 20 promoters for the expression of the 382 genes were determined, including Klf factors, Patz1, Sp factors, two Zinc finger proteins, and Egr1 (SI Appendix, Fig. S7C). The gene enrichment KEGG pathway and promoter analysis were done in a similar manner with CD103− DCs (SI Appendix, Fig. S8), CD8− plasmacytoid DCs (pDCs) (SI Appendix, Fig. S9), CD8− DCs (SI Appendix, Fig. S10), and ILC2 (SI Appendix, Fig. S11), whose expressions also increased and were steroid resistant (Dataset S6). Importantly, only the apoptosis pathway was associated with the 121 steroid-resistant genes in CD8− DCs. Although complex, our data suggest that multiple pathways are activated to drive the inflammatory response that underpins the exacerbation. These interactions reflect the cellular heterogeneity of the response and the complex intercellular signaling arrangement. However, looking at the nature of cellular response and clinical data (36–38) suggests that ILC2s may be a critical link to the exacerbation response.
IL-4 and IL-13 Are Differentially Expressed by Basophils, ILC2, and CD8+ Memory T Cells.
It is known that multiple proinflammatory factors play a key role in the pathogenesis of asthma and asthma exacerbation (19, 39). We therefore evaluated the expression of classical factors associated with Th1, Th2, ILC2s, anti-infection responses, and airway remodeling by those cell clusters. Almost all clusters in the SAL, VEH, and DEX groups highly expressed CCL5, whose expression was not significantly influenced by dexamethasone treatment. Several genes were among the most highly expressed factors, including IL-4, IL-13, and IFNγ (Fig. 4A). Of note, IL-4 and IL-13 are key cytokines driving type 2 responses. Dexamethasone treatment did not suppress IL-4 expression by CD8− DCs and basophils; instead, it promoted Ly6C+ and Ly6C− monocytes, Tregs, and CD103− DCs to generate this cytokine. Interestingly, IL-13 was differentially expressed by CD8+ memory T cells, ILC2, and basophils, and this expression was glucocorticoid insensitive. Surprisingly, dexamethasone treatment slightly promoted expression of IL-13 by CD8+ memory T cells. IFNγ was expressed by NK1.1+ NK cells, CD11b+ macrophages, and CD8+ memory T cells, and this expression was steroid resistant. Several clusters, including Ly6C+ and Ly6C− monocytes and CD103− and CD8− DCs, expressed CCL24, which is one of the major eosinophil-attracting chemokines. IL-1β was highly expressed by CD11b− macrophages, CD8+ memory T cells, and CD103− and CD8− DCs, and this expression was steroid resistant. The major cellular sources of IL-18 appeared to be Ly6C− monocytes, CD11b+ macrophages, and CD8− DCs upon HDM/LPS challenges. Interestingly, IL-10 was detected mainly in Ly6C− monocytes, CD11b+ macrophages, CD103− DCs, and CD8− DCs. Tregs also expressed IL-10, but at significantly lower levels. Differential expression of key inflammatory factors by multiple key immune cells was observed for IL-1α, IL-2, IL-6, IL-12a, CCL2, TNF, TGF-β1, TGF-β2, and TGF-β3. We were not able to detect the expression of other cytokines including IL-3, IL-5, IL-17a, IL-25, IL-33, and TSLP in any cluster of the SAL, VEH, and DEX groups.
Fig. 4.

GO analysis was performed with the ClueGO plugin of Cytoscape to show overrepresented functional categories of aforementioned cytokines and chemokines (Fig. 5B and Dataset S7). The highest number of genes demonstrating significant enrichment was involved in positive regulation of leukocyte differentiation, positive regulation of tyrosine phosphorylation of STAT protein, and others. Of note, neutrophil migration and the lipopolysaccharide-mediated signaling pathway were also identified.
Fig. 5.

IL-13 Critically Regulated HDM/LPS Induced Steroid-Resistant AHR and Airway Inflammation.
IL-13 has been shown to play a central role in the pathogenesis of allergic airway diseases (40, 41). To validate our findings, IL-13-tdTomato-reporter mice were employed to determine its cellular source. CD8+ T cells, ILC2, and basophils were confirmed to produce IL-13, whose production by CD4+ T cells was almost negligible at this time point (Fig. 5 and SI Appendix, Fig. S12). Other immune cells did not produce detectable levels of IL-13.
To determine the functional role of IL-13 during HDM/LPS and/or DEX exposure, we neutralized IL-13 by the administration of anti-IL-13 mAb after AAD was established in mouse (Fig. 6A). Anti–IL-13 mAb abolished AHR in both VEH and DEX groups, as compared to the isotype control (Fig. 6B). Furthermore, IL-13 neutralization significantly suppressed the infiltration of eosinophils, lymphocytes, and neutrophils in the BALF of the aforementioned two groups. Goblet cell hyperplasia, mucus hypersecretion, increased numbers of eosinophils, and heightened inflammation scores in lung tissues were significantly reduced by anti–IL-13 mAb treatment in both VEH and DEX groups (Fig. 6C). Treatment with anti–IL-13 mAb significantly reduced the levels of pathological changes in lung (Fig. 6D). Collectively, these data suggest that IL-13 is essential for the induction of LPS-regulated steroid-resistant AHR and airway inflammation and that this cytokine is predominantly produced by CD8+ T cells, ILC2, and basophils during LPS-induced asthma exacerbation.
Fig. 6.

Canonical Pathways and Upstream Regulators of Basophils, ILC2, and CD8+ Memory T Cells Are Up-Regulated during Asthma Exacerbation.
To understand the activation of the intracellular molecular network of basophils, ILC2, and CD8+ memory T cells during asthma exacerbation, we performed gene set enrichment analysis using IPA with the “Canonical Pathway” and “Upstream Regulator” options. A total of 794 genes of basophils, 796 genes of ILC2, and 645 genes of CD8+ memory T cells that were differentially expressed by at least five-fold in the VEH group were selected. Meanwhile, a total of 729 genes of basophils, 801 genes of ILC2, and 665 genes of CD8+ memory T cells that were differentially expressed by at least five-fold in the DEX group were selected. Only canonical pathways and upstream regulators with significant predicted activation states (Z-score > 2) were presented (Fig. 7 and Datasets S8 and S9). The subset of genes of the canonical pathways with the highest activation Z-score was primarily involved with “EIF2 Signaling,” followed by “Oxidative Phosphorylation” and “Signaling by Rho Family GTPase” (Fig. 7A). Importantly, the activation of these pathways, which regulate cellular function including stress-related response, generation of reactive oxygen species, and cell motility, was steroid resistant and was shared by these three cell clusters.
Fig. 7.

Upstream regulator analysis predicted the heightened activation of a number of upstream regulators including STK11, NFE2L2, Hbb-b1, IFNγ, CD38, and IL-5 in response to HDM/LPS challenge (Fig. 7B). Interestingly, the activation of most upstream regulators was steroid resistant, except that of TLR4 in CD8+ memory T cells, TLR3/IL-1β in ILC2, and IL-4/Myc in basophils. MyD88, TLR3, TLR4, EGR1, Retnlb, and NF-κB were uniquely activated in basophil clusters. Likewise, EGFR and CHUK were uniquely activated in ILC2, TCF3, TCF4, and SMAD3 in CD8+ memory T cells. Surprisingly, dexamethasone treatment promoted the activation of STAT1/BCR/PPARGC1A in CD8+ memory T cells, IKBKG/MKNK1PITX2/CD28 in ILC2, and TICAM1/IFNAR in basophils. By using dot plots, we demonstrated that the transcripts of these upstream regulators were substantially up-regulated in these three cell clusters of SAL, VEH, and DEX groups (Fig. 7C).
Discussion
Due to the heterogeneous nature of asthma exacerbation, previous approaches to understanding the cellular framework of the immune response in lung have clear limitations as there are no suitable tools to investigate differences between histologically indistinguishable cells. Our single-cell profiling of lung immune cells adds to the understanding of the cellular makeup of the lung organ and provides a vital information resource of the immune cellular network in lung during asthma exacerbation and steroid resistance. Here, we have employed scRNA-seq to distinguish different transcriptional profiles present in the immune cells of lung at the single-cell level and have decoded the network in steroid resistance of asthma exacerbation.
We first established a mouse model of HDM/LPS-induced, steroid-resistant asthma exacerbation with exaggerated and prolonged AHR and airway inflammation. Single-cell transcriptome analysis was able to distinguish multiple clusters on the basis of their messenger RNA expression profiles. We were able to annotate the scRNA samples into 20 major cell clusters by using SingleR (42). Surprisingly, eosinophils were not annotated among these cell clusters. This is likely due to the fact that differentiated eosinophils in general have low levels of RNA content and high levels of cellular RNases and RNA inhibitory compounds. During the quality control of data processing, activated eosinophils cannot be distinguished from other cells with low quality and quantity of RNA transcripts so they are filtered out. This limitation could be overcome with highly improved droplet scRNA-seq techniques in the future. Treatments with vehicle and dexamethasone in the presence of HDM/LPS challenge did not alter the profile of the identified cell clusters, although dexamethasone treatment did suppress HDM/LPS-induced increased levels of CD11b− macrophages. It appears that the infiltration of the majority of immune cells upon HDM/LPS challenge is resistant to dexamethasone treatment, which is consistent with previous reports of glucocorticoid function (43–47). Our observation of increased levels of many immune cells are likely to reflect the timing of sampling where we observed a rebound effect and the types of cytokines and inflammatory mediators present during inflammation.
DCs in lung are a highly heterogeneous population, containing conventional DCs (cDC), plasmacytoid DCs, and monocyte-derived DCs (48). We were able to determine the increase of two cDC populations, CD103− DCs and CD8− DCs, which had increased expression of proinflammatory cytokines including IL-1α, IL-1β, IL-2, and IL-10. The expression of these cytokines was steroid resistant. Interestingly, CD8− DCs produced far more IL-4 than CD103− DCs following HDM/LPS administration, suggesting the crucial role of former cells in driving Th2 cell proliferation and activation. Our results have revealed that airway exposure to inflammatory stimuli may have long-lasting effects on the pulmonary immune landscape for subsequent induction of exacerbations.
Dexamethasone treatment promoted the IL-4 production by CD103− DCs, further indicating the pleiotropic effect of dexamethasone function. Although glucocorticoids are still widely used clinically to suppress inflammation, the long-term and unwanted side effect of this medicine should be carefully considered. By contrast, CD8− pDC was less active when compared to the two identified DC populations. These cells expressed only IL-6, IL-12a, CCL5, TNF, and TGFβ1, suggesting that their role in the disease pathogenesis remains to be further elucidated.
Residential macrophages are CD11b negative in healthy mouse lung. When activated by inflammatory stimulants, the expression of CD11b on these cells is highly increased (49). Indeed, an increase in CD11b+ alveolar macrophages has been linked to asthma, chronic obstructive pulmonary disease, and respiratory infection (50). Although cytokine profile analysis revealed that CD11b+ macrophages expressed high levels of IL-10, they also produced high levels of proinflammatory factors, including IL-1β, IL-6, IL-18, CCL2, IFNγ, and TNF. Interestingly, evidence has shown that CD11b+ macrophages not only generate IL-10 to suppress and control inflammation, but also produce proinflammatory cytokines including IL-6, IL-1β, and TNF (51). Although further studies are needed to characterize these cells in detail, our investigation indicates that CD11b+ macrophages are able to differentially express both anti- and proinflammatory factors, suggesting the complex nature of these cells in the pathogenesis of asthma exacerbations.
Inflammatory factors are essential in driving the development and exacerbation of asthma (19). We were able to identify the expression of multiple key factors that drive persistence of inflammation, antiviral responses, and tissue remodeling. These proinflammatory factors critically regulate the migration, differentiation, and cytokine/chemokine production of leukocytes, STAT intracellular signaling, and LPS-induced immune responses. Both ILC2-activiting cytokines including IL-25, IL-33, TSLP, and Th17 were not detected in the samples at these later time points after allergen exposure. Infection-associated cytokines such as IL-18, IFNγ, and TNF significantly increased. TGFβ1 modulates allergic airway inflammation and airway remodeling (52). Interestingly, innate immune cells including monocyte, macrophage, and DC cell clusters consistently express TGFβ1 but not TGFβ2 or TGFβ3, adaptive immune cells including CD4+ and CD8+ T cells and CD19+ B cell clusters also express TGFβ1 in response to HDM/LPS treatment. This suggests that both innate and adaptive immune cells play an important role in airway remodeling in asthma. We were not able to detect IL-5 expression in any cell cluster, which is likely to reflect the time since the last allergen challenge (7 d). However, CCL24 that was produced by three monocyte clusters and CD103− and CD8− DC clusters is able to attract eosinophils in BALF in the absence of IL-5. Of note, we were able to identify two type 2 cytokines, IL-4 and IL-13, that are key players in the pathogenesis of asthma (53). In particular, IL-13 plays a central role in orchestrating AHR, airway inflammation, mucus hypersecretion, and airway remodeling (41). We identified that IL-13 was predominantly expressed by CD8+ memory T cells, ILC2, and basophils upon HDM and LPS challenge, whose expression was largely steroid resistant in the former two cells. By contrast, CD4+ T cells expressed very low levels of IL-4 and IL-13, emphasizing the critical role of non-Th2 immune cells in driving the exaggerated inflammatory response inducing the exacerbation. We were able to validate our analysis by IL-13 reporter mice to demonstrate the predominant cellular sources of this cytokine: CD8+ T cells, ILC2, and basophils. It appears that the percentage of IL-13+ILC2 correct is greater than those of IL-13 producing CD8+ T cells and basophils. Classic IL-13–producing cells, CD4+ Th2 cells, ceased to be the producers during LPS-induced asthma exacerbation. The key role of IL-13 in the pathogenesis was also proven by anti–IL-13 mAb treatment during exacerbation phase, as the neutralization of IL-13 resulted in diminished AHR, infiltration of eosinophils and neutrophils, and suppressed mucus production.
Furthermore, upstream regulator and canonical pathway analysis of these three cell clusters by IPA unraveled multiple top regulators and pathways. By combining expression analysis, STK11, NFE2L2, IFNγ, CD38, and INSR were among the top upstream regulators whose activation (Z-scores) was not affected by the effect of dexamethasone on ILC2 and on CD8+ memory T cells. STK11, a serine/threonine kinase, is a tumor-suppressive gene and critically regulates cell motility, differentiation, and metabolism (54). NFE2L2 is a transcription regulator that drives the host defense to oxidative stress by promoting synthesis of antioxidant and detoxifying enzymes including glutathione S-transferase and sulfotransferase (55). Molecules underpinning proinflammatory NF-κB pathway activation, including NF-KB1, NF-KB2, Rela, Relb, and Rel, were also differentially identified in the three cell clusters. Interestingly, Rela and Relb in ILC2 and basophils were steroid resistant, and NF-KB1 and NK-KB2 in basophils differentially responded to dexamethasone treatment. Although the contribution of these pathways to the production of IL-4 and IL-13 and asthma exacerbation remains unclear, this reflects the complex nature of asthma exacerbation, to which both anti- and proinflammatory programs contribute.
Our analysis of the top 20 pathways of basophil, ILC2, and CD8+ memory T cell clusters also revealed that multiple conservative pathways are presumed to play a role in disease pathogenesis. These data strongly indicate that common or unique host responses at the single-cell level are modulated by a large number of low-risk intracellular events, which are not likely to be easily identified through the application of conventional analysis approaches. For example, activation of EIF2 signaling is crucial for stress-induced regulation of translation in mammalian cells and often occurs upon infection (56). Overreaction of oxidative phosphorylation leads to excessive generation of reactive oxygen species and thus cell damage and diseases (57). Other pathways, including PKC-, Ephrin receptor-, phospholipase C-, Rho-, and Rac-signaling pathways, that could not be suppressed by dexamethasone treatment, are likely conservatively employed by these three cell clusters for cell mobility, inflammatory responses, and cytokine production. By contrast, several unique pathways may critically regulate cell-specific functions, such as leukocyte extravasation signaling and Tec kinase signaling of CD8+ memory T cells, dendritic cell maturation signaling of basophils, and fMLP signaling in neutrophils of ILC2. Interestingly, dexamethasone treatment had different impacts on these cells. This reflects highly complicated cellular and molecular interactions in the lung microenvironment to address different aspects of environmental insults. Further investigation is needed to elucidate how these networks contribute to the course of airway chronic inflammation and disease pathogenesis.
In conclusion, our single-cell transcriptional analysis of mouse lung immune cells provides an insightful framework for understanding the shared and distinct expression patterns of inflammatory genes at the individual cell level. It provides a cellular resource for uncovering multicellular signaling pathways critical for disease progression and for assessing the importance of distinct cell populations that emerge during the course of disease. IL-13 produced by CD8+ memory T cells, ILC2, and basophils is a key cytokine in driving the pathogenesis of asthma exacerbation in this model. This detailed information serves as a benchmark for the development of cell-targeted therapies to treat and prevent asthma exacerbations.
Data Availability
All study data are included in the article and SI Appendix.
Acknowledgments
We thank the staff members of the Hunter Medical Research Institute Bioresources facility at the University of Newcastle and LC Sciences (Hangzhou, China) for their support and expertise.
Supporting Information
Appendix (PDF)
- Download
- 52.66 MB
Dataset_S01 (XLSX)
- Download
- 8.96 KB
Dataset_S02 (XLSX)
- Download
- 4.07 MB
Dataset_S03 (XLSX)
- Download
- 10.73 KB
Dataset_S04 (XLSX)
- Download
- 84.49 KB
Dataset_S05 (XLSX)
- Download
- 665.97 KB
Dataset_S06 (XLSX)
- Download
- 30.46 KB
Dataset_S07 (CSV)
- Download
- 47.99 KB
Dataset_S08 (XLSX)
- Download
- 202.07 KB
Dataset_S09 (XLSX)
- Download
- 132.96 KB
References
1
A. J. Woolcock, Steroid resistant asthma: What is the clinical definition? Eur. Respir. J. 6, 743–747 (1993).
2
C. Gautier, D. Charpin, Environmental triggers and avoidance in the management of asthma. J. Asthma Allergy 10, 47–56 (2017).
3
P. G. Holt, D. H. Strickland, Interactions between innate and adaptive immunity in asthma pathogenesis: New perspectives from studies on acute exacerbations. J. Allergy Clin. Immunol. 125, 963–972, quiz 973–974 (2010).
4
M. Becerra-Díaz, M. Wills-Karp, N. M. Heller, New perspectives on the regulation of type II inflammation in asthma. F1000 Res. 6, 1014 (2017).
5
L. Cohn, J. A. Elias, G. L. Chupp, Asthma: Mechanisms of disease persistence and progression. Annu. Rev. Immunol. 22, 789–815 (2004).
6
S. C. Eisenbarth, S. Cassel, K. Bottomly, Understanding asthma pathogenesis: Linking innate and adaptive immunity. Curr. Opin. Pediatr. 16, 659–666 (2004).
7
A. M. Singh, W. W. Busse, Asthma exacerbations. 2: Aetiology. Thorax 61, 809–816 (2006).
8
J.-Y. Wang, The innate immune response in house dust mite-induced allergic inflammation. Allergy Asthma Immunol. Res. 5, 68–74 (2013).
9
Y. Matsumura, Role of allergen source-derived proteases in sensitization via airway epithelial cells. J. Allergy (Cairo) 2012, 903659 (2012).
10
A. Trompette et al., Allergenicity resulting from functional mimicry of a Toll-like receptor complex protein. Nature 457, 585–588 (2009).
11
O. Michel et al., Severity of asthma is related to endotoxin in house dust. Am. J. Respir. Crit. Care Med. 154, 1641–1646 (1996).
12
N. H. ten Hacken et al., Elevated serum interferon-gamma in atopic asthma correlates with increased airways responsiveness and circadian peak expiratory flow variation. Eur. Respir. J. 11, 312–316 (1998).
13
M. Yoshida, R. M. Watson, T. Rerecich, P. M. O’Byrne, Different profiles of T-cell IFN-gamma and IL-12 in allergen-induced early and dual responders with asthma. J. Allergy Clin. Immunol. 115, 1004–1009 (2005).
14
M. E. Dahl, K. Dabbagh, D. Liggitt, S. Kim, D. B. Lewis, Viral-induced T helper type 1 responses enhance allergic disease by effects on lung dendritic cells. Nat. Immunol. 5, 337–343 (2004).
15
H. D. Held, S. Uhlig, Mechanisms of endotoxin-induced airway and pulmonary vascular hyperreactivity in mice. Am. J. Respir. Crit. Care Med. 162, 1547–1552 (2000).
16
Z. Zhu, S. Y. Oh, T. Zheng, Y. K. Kim, Immunomodulating effects of endotoxin in mouse models of allergic asthma. Clin. Exp. Allergy 40, 536–546 (2010).
17
S. C. Eisenbarth et al., Lipopolysaccharide-enhanced, toll-like receptor 4-dependent T helper cell type 2 responses to inhaled antigen. J. Exp. Med. 196, 1645–1651 (2002).
18
J. Daan de Boer et al., Lipopolysaccharide inhibits Th2 lung inflammation induced by house dust mite allergens in mice. Am. J. Respir. Cell Mol. Biol. 48, 382–389 (2013).
19
P. A. B. Wark, P. G. Gibson, Asthma exacerbations. 3: Pathogenesis. Thorax 61, 909–915 (2006).
20
H. Gao, S. Ying, Y. Dai, Pathological roles of neutrophil-mediated inflammation in asthma and its potential for therapy as a target. J. Immunol. Res. 2017, 3743048 (2017).
21
A. Ray, J. K. Kolls, Neutrophilic inflammation in asthma and association with disease severity. Trends Immunol. 38, 942–954 (2017).
22
J. W. Upham, L. P. Chung, Optimising treatment for severe asthma. Med. J. Aust. 209, S22–S27 (2018).
23
Y. Y. Zhu, E. M. Machleder, A. Chenchik, R. Li, P. D. Siebert, Reverse transcriptase template switching: A SMART approach for full-length cDNA library construction. Biotechniques 30, 892–897 (2001).
24
M. Yang, R. K. Kumar, P. S. Foster, Interferon-gamma and pulmonary macrophages contribute to the mechanisms underlying prolonged airway hyperresponsiveness. Clin. Exp. Allergy 40, 163–173 (2010).
25
M. Yang, R. K. Kumar, P. S. Foster, Pathogenesis of steroid-resistant airway hyperresponsiveness: Interaction between IFN-γ and TLR4/MyD88 pathways. J. Immunol. 182, 5107–5115 (2009).
26
G. X. Zheng et al., Massively parallel digital transcriptional profiling of single cells. Nat. Commun. 8, 14049 (2017).
27
Z. Xiaoshu, L. Hong-Dong, G. Lilu, W. Fang-Xiang, W. Jianxin, Analysis of single-cell RNA-seq data by clustering approaches. Curr. Bioinform. 14, 314–322 (2019).
28
A. Krämer, J. Green, J. Pollard, Jr, S. Tugendreich, Causal analysis approaches in ingenuity pathway analysis. Bioinformatics 30, 523–530 (2014).
29
A. Alexa, J. Rahnenfuhrer, topGO: Enrichment Analysis for Gene Ontology (R package version 2.40.0, Bioconductor, 2020). https://bioconductor.org/packages/release/bioc/html/topGO.html. Accessed 22 June 2020.
30
G. Bindea et al., ClueGO: A Cytoscape plug-in to decipher functionally grouped gene ontology and pathway annotation networks. Bioinformatics 25, 1091–1093 (2009).
31
S. X. Ge, D. Jung, R. Yao, ShinyGO: A graphical gene-set enrichment tool for animals and plants. Bioinformatics 36, 2628–2629 (2020).
32
X.-T. Huang et al., Technical advances in single-cell RNA sequencing and applications in normal and malignant hematopoiesis. Front. Oncol. 8, 582 (2018).
33
P. A. Szabo et al., Single-cell transcriptomics of human T cells reveals tissue and activation signatures in health and disease. Nat. Commun. 10, 4706 (2019).
34
E. A. K. DePasquale et al., DoubletDecon: Deconvoluting doublets from single-cell RNA-sequencing data. Cell Rep. 29, 1718–1727.e8 (2019).
35
P. Langfelder, S. Horvath, WGCNA: An R package for weighted correlation network analysis. BMC Bioinformatics 9, 559 (2008).
36
C. A. Christianson et al., Persistence of asthma requires multiple feedback circuits involving type 2 innate lymphoid cells and IL-33. J. Allergy Clin. Immunol. 136, 59–68.e14 (2015).
37
P. Nagakumar et al., Type 2 innate lymphoid cells in induced sputum from children with severe asthma. J. Allergy Clin. Immunol. 137, 624–626.e6 (2016).
38
S. G. Smith et al., Increased numbers of activated group 2 innate lymphoid cells in the airways of patients with severe asthma and persistent airway eosinophilia. J. Allergy Clin. Immunol. 137, 75–86.e8 (2016).
39
F. T. Ishmael, The inflammatory response in the pathogenesis of asthma. J. Am. Osteopath. Assoc. 111(11, suppl. 7)S11–S17 (2011).
40
M. Wills-Karp et al., Interleukin-13: Central mediator of allergic asthma. Science 282, 2258–2261 (1998).
41
G. Grünig et al., Requirement for IL-13 independently of IL-4 in experimental asthma. Science 282, 2261–2263 (1998).
42
D. Aran et al., Reference-based analysis of lung single-cell sequencing reveals a transitional profibrotic macrophage. Nat. Immunol. 20, 163–172 (2019).
43
D. C. Dale, A. S. Fauci, I. V. Guerry D, S. M. Wolff, Comparison of agents producing a neutrophilic leukocytosis in man. Hydrocortisone, prednisone, endotoxin, and etiocholanolone. J. Clin. Invest. 56, 808–813 (1975).
44
E. H. Dunsky, B. Zweiman, E. Fischler, D. A. Levy, Early effects of corticosteroids on basophils, leukocyte histamine, and tissue histamine. J. Allergy Clin. Immunol. 63, 426–432 (1979).
45
D. Franchimont et al., Tumor necrosis factor alpha decreases, and interleukin-10 increases, the sensitivity of human monocytes to dexamethasone: Potential regulation of the glucocorticoid receptor. J. Clin. Endocrinol. Metab. 84, 2834–2839 (1999).
46
T. J. Creed et al., The effects of cytokines on suppression of lymphocyte proliferation by dexamethasone. J. Immunol. 183, 164–171 (2009).
47
S. Liu et al., Steroid resistance of airway type 2 innate lymphoid cells from patients with severe asthma: The role of thymic stromal lymphopoietin. J. Allergy Clin. Immunol. 141, 257–268.e6 (2018).
48
M. Guilliams et al., Unsupervised high-dimensional analysis aligns dendritic cells across tissues and species. Immunity 45, 669–684 (2016).
49
P. J. Murray et al., Macrophage activation and polarization: Nomenclature and experimental guidelines. Immunity 41, 14–20 (2014).
50
M. Duan et al., CD11b immunophenotyping identifies inflammatory profiles in the mouse and human lungs. Mucosal Immunol. 9, 550–563 (2016).
51
G. Arango Duque, A. Descoteaux, Macrophage cytokines: Involvement in immunity and infectious diseases. Front. Immunol. 5, 491 (2014).
52
H. Li et al., Genetic polymorphisms in transforming growth factor beta-1 (TGFB1) and childhood asthma and atopy. Hum. Genet. 121, 529–538 (2007).
53
P. S. Foster et al., Modeling TH 2 responses and airway inflammation to understand fundamental mechanisms regulating the pathogenesis of asthma. Immunol. Rev. 278, 20–40 (2017).
54
D. B. Shackelford, R. J. Shaw, The LKB1-AMPK pathway: Metabolism and growth control in tumour suppression. Nat. Rev. Cancer 9, 563–575 (2009).
55
J. Lorsch, Methods in enzymology. Laboratory methods in enzymology: RNA. Preface. Methods Enzymol. 530, xxi (2013).
56
N. Shrestha et al., Eukaryotic initiation factor 2 (eIF2) signaling regulates proinflammatory cytokine expression and bacterial invasion. J. Biol. Chem. 287, 28738–28744 (2012).
57
M. Nita, A. Grzybowski, The role of the reactive oxygen species and oxidative stress in the pathomechanism of the age-related ocular diseases and other pathologies of the anterior and posterior eye segments in adults. Oxid. Med. Cell. Longev. 2016, 3164734 (2016).
Information & Authors
Information
Published in
Classifications
Copyright
© 2021. Published under the PNAS license.
Data Availability
All study data are included in the article and SI Appendix.
Submission history
Published online: January 4, 2021
Published in issue: January 12, 2021
Keywords
Acknowledgments
We thank the staff members of the Hunter Medical Research Institute Bioresources facility at the University of Newcastle and LC Sciences (Hangzhou, China) for their support and expertise.
Notes
This article is a PNAS Direct Submission. A.M. is a guest editor invited by the Editorial Board.
Authors
Competing Interests
The authors declare no competing interest.
Metrics & Citations
Metrics
Citation statements
Altmetrics
Citations
Cite this article
Single-cell transcriptomic analysis reveals the immune landscape of lung in steroid-resistant asthma exacerbation, Proc. Natl. Acad. Sci. U.S.A.
118 (2) e2005590118,
https://doi.org/10.1073/pnas.2005590118
(2021).
Copied!
Copying failed.
Export the article citation data by selecting a format from the list below and clicking Export.
Cited by
Loading...
View Options
View options
PDF format
Download this article as a PDF file
DOWNLOAD PDFLogin options
Check if you have access through your login credentials or your institution to get full access on this article.
Personal login Institutional LoginRecommend to a librarian
Recommend PNAS to a LibrarianPurchase options
Purchase this article to access the full text.