Figures
Abstract
The human airway epithelium is the initial site of SARS-CoV-2 infection. We used flow cytometry and single cell RNA-sequencing to understand how the heterogeneity of this diverse cell population contributes to elements of viral tropism and pathogenesis, antiviral immunity, and treatment response to remdesivir. We found that, while a variety of epithelial cell types are susceptible to infection, ciliated cells are the predominant cell target of SARS-CoV-2. The host protease TMPRSS2 was required for infection of these cells. Importantly, remdesivir treatment effectively inhibited viral replication across cell types, and blunted hyperinflammatory responses. Induction of interferon responses within infected cells was rare and there was significant heterogeneity in the antiviral gene signatures, varying with the burden of infection in each cell. We also found that heavily infected secretory cells expressed abundant IL-6, a potential mediator of COVID-19 pathogenesis.
Author summary
SARS-CoV-2 infects the respiratory tract, targeting cells of the diverse airway epithelium. Infection outcomes depend on several factors that may vary in this heterogenous population including viral tropism, antiviral immunity, and response to antiviral therapies like remdesivir. We found that SARS-CoV-2 infects an array of airway epithelial cells, relying on the host protease TMPRSS2 for entry. Ciliated epithelial cells were the dominant target, and remdesivir blocked viral replication across multiple cell types. We uncovered cellular heterogeneity in early antiviral immunity to SARS-CoV-2 and identified cell type-specific ISGs associated with either high levels of viral replication or protection from infection.
Citation: Fiege JK, Thiede JM, Nanda HA, Matchett WE, Moore PJ, Montanari NR, et al. (2021) Single cell resolution of SARS-CoV-2 tropism, antiviral responses, and susceptibility to therapies in primary human airway epithelium. PLoS Pathog 17(1): e1009292. https://doi.org/10.1371/journal.ppat.1009292
Editor: Kanta Subbarao, The Peter Doherty Institute and Melbourne University, AUSTRALIA
Received: October 20, 2020; Accepted: January 7, 2021; Published: January 28, 2021
Copyright: © 2021 Fiege et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All sequencing data are available from NCI GEO accession number: GSE157526. Code is available from our github page: https://github.com/heznanda/scrnaseq-hybrid-cov2.
Funding: This work was supported by National Institutes of Health (https://www.nih.gov/) R01 AI148669 to RAL, R01 AI153602 to VDM and UL1TR002494 to SSS. WEM was supported by the National Institutes of Health (https://www.nih.gov/) T32 HL007741. JMT was supported by the National Institutes of Health (https://www.nih.gov/) T32 AI055433.The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist
Introduction
SARS-CoV-2, the virus responsible for COVID-19, primarily infects cells of the respiratory tract. The cellular tropism of SARS-CoV-2 may impact several aspects of the disease, including viral spread within and between hosts, mechanisms of immune control of infection or tissue pathology, and the therapeutic response to promising antivirals. Normal human tracheal bronchial epithelial (nHTBE) cells represent a diverse mix of ciliated epithelial cells, secretory cells, and basal cells that form a pseudostratified epithelium when cultured at the air-liquid interface, phenocopying the upper airway in humans [1,2]. Importantly, cells in this culture system also express endogenous levels of critical host factors including ACE2 and host proteases such as TMPRSS2 that are needed for SARS-CoV-2 viral entry [3–7]. This model also demonstrates key aspects of host antiviral epithelial immunity [8,9]. Recently, several studies using primary human lung cell cultures and respiratory cells isolated from SARS-CoV-2 infected patients have identified SARS-CoV-2 tropism for ciliated and secretory cells in the upper airway [10–14]. However, the heterogeneity of virus replication and induction of antiviral genes and proinflammatory cytokines within these cells is still unknown.
Remdesivir (GS-5734) has emerged as a promising direct antiviral therapy against SARS-CoV-2, with potent activity demonstrated against several coronaviruses [15,16]. A landmark clinical trial found that remdesivir treatment of hospitalized individuals with COVID-19 improved median recovery time [17], and this drug is now approved for COVID-19 under emergency use authorization by the U.S. Food and Drug Administration. Remdesivir is a prodrug that is metabolized in cells to the nucleotide analog remdesivir triphosphate, which interferes with coronavirus replication through delayed RNA chain termination [10,18–20]. Recent studies have identified differential efficacy of remdesivir against SARS-CoV-2 in a range of cell culture systems linked to metabolism of the prodrug to the active form [10]. In addition to differential metabolism, other factors that may impact the variable efficacy of this drug in different cell types include differential drug uptake and heterogeneous permissibility of each cell type to viral entry and replication. While remdesivir clearly exhibits antiviral activity against SARS-CoV-2 in nHTBE cultures, it is not known if there are cell type-dependent differences in drug efficacy.
Following infection, coronaviruses are recognized by MDA5 and RIGI leading to the production of type I and III interferons (IFNs), which induce transcriptional programs that mobilize cellular antiviral defenses. Coronaviruses use several mechanisms to successfully evade detection resulting in rare and heterogeneous IFN production [21,22], similar to influenza virus infected cells [23–26]. During influenza A virus infection, we have previously identified interferon stimulated genes (ISGs) specifically induced in cells supporting high levels of virus replication and we have defined cell type-specific ISGs [27,28]. Additionally, we and others have found significant heterogeneity in antiviral responses across different cell types [28,29]. Cell type-specific responses and the degree of heterogeneity in antiviral responses can dictate the outcome of immune responses and infection.
Here, we use nHTBE cells infected with SARS-CoV-2 to demonstrate that remdesivir reduces viral replication uniformly in all susceptible cell types within the upper respiratory tract. Additionally, we demonstrate that TMPRSS2 is the primary host protease used for SARS-CoV-2 entry across cell types in the upper airway. Using single cell RNA sequencing, we further define SARS-CoV-2 tropism and the induction of antiviral and proinflammatory immune responses. We also uncover cellular heterogeneity in SARS-CoV-2 replication within the human airway epithelium and determine the relationship between productive viral infection and induction of IFN. Additionally, we define respiratory epithelial ISG expression in response to SARS-CoV-2 infection and uncover epithelial cell type-specific ISGs, and ISGs associated with either high levels of viral replication or protection from infection. Finally, we demonstrate that the proinflammatory cytokine IL-6 is primarily upregulated within heavily infected secretory cells. Together, these data help to define the early SARS-CoV-2 replication and antiviral landscape within the respiratory tract.
Results
SARS-CoV-2 tropism and cell type-specific efficacy of remdesivir in differentiated tracheal bronchial epithelial cells
To determine how SARS-CoV-2 tropism varies within the upper airway we infected nHTBE cells differentiated at air-liquid interface with a SARS-CoV-2 infectious clone generated to express the fluorescent protein mNeonGreen (mNG) [30]. Using flow cytometry, we detected 12.5% of untreated nHTBE cells expressing mNG after icSARS-CoV-2-mNG infection (Fig 1A and 1B); treatment with remdesivir or IFNβ substantially reduced the percentage of mNG positive cells showing clear antiviral effect. We subsequently identified different cell types using the expression of several surface markers and SiR-tubulin (Fig 1C) [31] and found that the ratio of different cell types in each culture remained similar during infection, regardless of anti-viral treatments (Fig 1D). Among mNG+ SARS-CoV-2 infected cells, we discovered a preponderance of ciliated (SiR-tubulin+ CD271-) and few infected secretory cells (CD66c+ SiR-tubulin- CD271-) (Fig 1E and 1F). We also identified a subset of infected cells lacking cell type-specific markers (termed “triple negative”) and a population of cells expressing intermediate levels of CD271 and CD66c, but negative for SiR-tubulin, termed “other”. We did not identify mNG+ basal cells (CD271+ SiR-tubulin- CD66c-). We also observed that ciliated cells expressed the highest mNG MFI, suggesting that this cell type is the most permissive for SARS-CoV-2 replication (Fig 1G). These data are consistent with recent reports demonstrating ciliated cells as a dominant target cell of SARS-CoV-2 [10–14,32], and they suggest that basal cells are relatively resistant to infection. We next evaluated if the antiviral activity of remdesivir varies in different cell types depending on cellular uptake and metabolism [10], which could differ dramatically within the diverse population of cells in the upper respiratory tract. To determine how remdesivir impacts viral replication across all susceptible cells in the upper airway, we pretreated nHTBE cells with remdesivir at 0.1μM, a dose expected to inhibit 80% of SARS-CoV-2 replication [10] and infected with icSARS-CoV-2-mNG. As a positive control for antiviral activity, we also pretreated cells with IFNβ, which completely abrogated infection in nHTBE cells (Fig 1A and 1B), consistent with recent results demonstrating sensitivity of SARS-CoV-2 to interferon [33–35]. Remdesivir reduced both the number of SARS-CoV-2 infected cells by ~5-fold and the expression level of mNG within infected cells (Fig 1B, 1F and 1G). Importantly, among infected cells, we found that remdesivir similarly inhibited replication regardless of cell type, reducing both the frequency of infection and the total mNG MFI (Fig 1F and 1G). These data suggest that cell type heterogeneity within the respiratory tract does not impact remdesivir efficacy.
nHTBE cells differentiated at air-liquid interface were untreated or treated with IFN-β or remdesivir prior to infection with icSARS-CoV-2-mNG at a MOI of 2.5. Cells were analyzed 48 hpi by flow cytometry. (A) Representative plots of live cells (B) Frequency of infected cells (mNG+). (C) Representative plots of cell subset gating. Ciliated (SiR-Tubulin+ CD271-), secretory (SiR-Tubulin- CD271- CD66c+), basal (SiR-Tubulin- CD271+ CD66c-), Triple negative (SiR-Tubulin- CD271- CD66c-) and other cells (SiR-Tubulin- CD217int CD66cint) (D) Frequency of cell subsets from live cells. (E) Frequency of cell subsets of infected cells (mNG+). (F) Total frequency of mNG+ infected cell types. (G) mNG gMFI within each cell subset. The data (B, D-F) are 3 biological replicates per group +/-SEM, 1 of 2 independent experiments.
Single cell RNA sequencing reveals heterogeneity in SARS-CoV-2 infection across differentiated tracheal bronchial epithelial cells
Tracheal-bronchial epithelial cells represent the first targets of SARS-CoV-2 and are key sentinels in the initiation of early antiviral responses. To further determine viral tropism and heterogeneity of replication and antiviral responses in the presence and absence of remdesivir, we infected nHTBE cells with SARS-CoV-2 and analyzed single cell transcriptomes via the 10X Genomics platform. We captured an average of 1729 cells per experimental condition and 1836 variable genes per cell. We observed significant heterogeneity in the level of virus replication supported per cell (Fig 2A, left). The frequency of cells with abundant viral RNA (defined as >0.1% of total reads per cell (S1A Fig) closely matched what we observed using flow cytometry (Fig 2A, right). Using this threshold for infection enables us to discriminate between cells that are productively infected from those with low-level viral RNA. Additionally, we confirmed viral replication kinetics with and without remdesivir by qRT-PCR, using a SARS-CoV-2 standard to calculate total viral genome copies (S1B Fig). Together these results show a distinct impact of remdesivir treatment on SARS-CoV-2 viral reads per cell, consistent with our flow cytometry results.
(A) Total viral reads per cell in each condition (left) and frequency of cells expressing at least 0.01% virus reads (right). (B) Cells from all samples were combined and analyzed by UMAP. (C) Cells with greater than 0.1% of total reads mapping to SARS-CoV-2. (D) Cells colored according to experimental condition.
To determine SARS-CoV-2 tropism and the impact of remdesivir, we analyzed cells combined from all experimental conditions, including mock-infected and those at 24 or 48 hours post-infection (hpi), with or without remdesivir treatment, by uniform manifold approximation and projection (UMAP) [36,37]. While some clusters of cells in this analysis were defined by cell type (basal, ciliated, and ionocytes), the majority were undefined (Fig 2B and S1 Table). We further interrogated the transcriptional changes in these clusters by examining the gene ontology of the most significantly upregulated genes. This analysis confirmed the identify of ciliated and ionocytes cells and identified secretory cells as having a strong antiviral signaling response (S1C Fig). However, this analysis did not further illuminate the identity of the undefined cells. We observed several clusters that were characterized by a preponderance of infected cells (Fig 2C); each of these clusters contained cells from different experimental conditions (Fig 2D). For example heavily infected cells in clusters 2 and 7 are derived from cells infected at 24hpi, 24hpi with remdesivir and 48hpi. Consistent with our flow cytometry results, cells of a variety of phenotypes localized to each of these dominant infected cell clusters. Interestingly, the majority of cells treated with remdesivir (Fig 2A) clustered most closely with mock-infected cells. This clearly demonstrates the antiviral efficacy of remdesivir and suggests that the presence of active SARS-CoV-2 replication is the strongest driver of transcriptomic phenotype in this cell population. These data suggest that SARS-CoV-2 infected airway epithelial cells adopt diverse transcriptional profiles dominated by both viral RNA level and cell type.
Single cell RNA sequencing identifies SARS-CoV-2 tropism in differentiated tracheal bronchial epithelial cells
Because UMAP clustering of combined samples was dominated by experimental condition, we next used t-squared neighbor embedded (t-SNE) clustering on individual samples to identify infected cell types and uninfected bystander cells (Fig 3A–3E and S2 Table). Ciliated cells were a major infected group at 24 and 48 hpi (Fig 3B and 3D), consistent with our flow cytometry data. We also observed heavily infected cells that did not cluster by cell type at 48 hpi (Fig 3D), possibly due to loss of cell markers during SARS-CoV-2 infection, as has also been observed with influenza infection [38]. To further evaluate the transcriptional changes in these clusters by examining the gene ontology of the most significantly upregulated genes. This analysis confirmed cell lineage for some clusters but was unable to further resolve the undefined cell types (S2A Fig). Similar to what others have observed [11,12,39], we found low level expression of ACE2 in our primary human epithelial cells (Fig 3A–3E), likely reflecting the limited sensitivity of scRNA sequencing. We therefore validated these findings using qRT-PCR and detected low ACE2 (30.4 ct) which increased after IFNβ treatment (26.07 ct, ~25 fold) (S2B Fig), consistent with other studies [7,35]. Interestingly, we also observed augmented expression of the host protease TMPRSS2, implicated in SARS-CoV-2 entry [5], within ciliated cells, secretory cells, and ionocytes, but reduced or no expression in basal cells (Fig 3A–3E). To determine if TMPRSS2 is the primary protease utilized for SARS-CoV-2 entry in differentiated nHTBE cells we treated cells with the TMPRSS2 inhibitor camostat. Camostat has been shown to block pseudovirus infection in undifferentiated primary human lung [5] and SARS-CoV-2 infection in immortalized Calu3 cells [5,40]. Consistent with these results, blocking TMPRSS2 resulted in a near complete abrogation of infection in both ciliated and secretory cells (Fig 3F). These data demonstrate that TMPRSS2 is the dominant protease entry factor for cells of the upper airway and may help explain the relative resistance of basal cells to infection with SARS-CoV-2. nHTBE cells pre-treated with remdesivir had reduced numbers of infected cells (Fig 3C and 3E), similar to our flow cytometry data. Cell tropism of SARS-CoV-2 was not altered with remdesivir treatment, as infected cells were predominately ciliated cells and a few secretory cells at 24 hpi. Remdesivir pre-treatment had a stronger effect at 48 hpi as there were almost no infected cells detected.
tSNE plots demonstrating clustering (left) and infected cells (right, blue) and volume plots for ACE2 and TMPRSS2 in mock (A), 24 hpi (B), 24 hpi with remdesivir (C), 48 hpi (D), and 48 hpi with remdesivir (E). (F) nHTBE cells were pretreated with vehicle or camostat, then infected with icSARS-CoV-2-mNG at MOI = 2.5 and analyzed for mNG in ciliated (SiR-tubulin+ CD271-), secretory cells (CD66c+ SiR-tubulin- CD271-), and triple negative cells. Numbers indicate fold change between vehicle and camostat treated cells.
Early antiviral immune response to SARS-CoV-2 in human airway epithelium
To evaluate the relationship between SARS-CoV-2 infection and antiviral immune responses we measured the induction of type I and III IFNs among individual cells in the human airway epithelium. To determine if there is a relationship between the levels of virus and the induction of IFN, we plotted any cell expressing IFNB, IFNL1, IFNL2, IFNL3 or IFNL4 by the level of SARS-CoV-2 infection. There was a strong positive correlation between levels of SARS-CoV-2 replication and induction of IFN (Fig 4A and 4B). Interestingly, we found cells with high levels of virus and low levels of IFN expression, potentially representing the successful evasion by SARS-CoV-2 of virus sensing pathways triggered by MDA5 or RIGI. At 48 hpi we observed two clusters of undefined cells, both with high SARS-CoV2 replication but divergent IFN responses: cluster 3 with high IFN induction and cluster 5 with reduced IFN. This disparate induction of IFN could indicate opposite outcomes of SARS-CoV-2 infection that drive differential transcriptomic clustering. We identified very few (~0.5% at 24hpi) cells actively transcribing IFN, consistent with observations in other viral infections where a small number of IFN-producing cells is sufficient to establish the antiviral state [25,26]. Despite this, we detected robust induction of ISGs (Fig 4C). We evaluated ISGs with the greatest differences across sample condition and found several discrete patterns of expression. This analysis revealed ISGs that were only present in clusters of cells characterized with high levels of SARS-CoV-2 replication (Fig 4C, MICB and GEM). We also observed an ISG (MT1F) that was preferentially induced within ciliated cells (Fig 4C). At 48 hpi cells with the highest levels of viral RNA (clusters 3 and 5) had significantly lower levels of expression of several ISGs (DDIT, IFITM3, LY6E, TNFSF10), suggesting an impaired cellular antiviral immune response that enables unchecked SARS-CoV-2 replication (Fig 4C).
Cells from Fig 3 were analyzed at 24hpi (A) and 48hpi (B) cells expressing interferon (IFNB1, IFNL1, IFNL2, IFNL3, IFNL4) were assessed for levels of combined IFN transcripts plotted against combined SARS-CoV-2 reads. Colors indicate clusters from Fig 3B and 3D, respectively. (C) Volume plots of cells from Fig 3 for ISGs with statistically significant expression between samples. (D) Volume plot of IL-6 expression.
Notably with SARS-CoV-2 infection we observed induction of IL-6, a pro-inflammatory cytokine implicated in the pathogenesis of COVID-19 [41], only within infected cells at 24 and 48 hours post-infection. Heavily infected secretory cells expressed the highest levels of IL-6 at 24 hpi (Fig 4D). To determine the impact of early antiviral therapy on the hyperinflammatory response that drives COVID-19 pathology, we also analyzed infected cells pretreated with remdesivir. Remdesivir blocked the expression of IL-6 but not the induction of ISGs, suggesting that direct antiviral therapy can prevent induction of potentially pathogenic inflammation while preserving protective IFN responses.
Discussion
The SARS-CoV-2 pandemic continues to spread, causing significant worldwide morbidity and mortality. Understanding early events in virus-host interactions in relevant model systems will be critical for understanding the pathophysiology of the disease as well as for evaluating antiviral and immunomodulatory drugs. Here we demonstrate that SARS-CoV-2 has preferential tropism for ciliated cells and secretory cells in the upper airway. Within these primary cells SARS-CoV-2 is only sensed by a small number of infected cells leading to IFN activation. Despite the low number of cells expressing IFN, there is a robust antiviral ISG response established. Finally, we demonstrate that remdesivir is capable of restricting virus replication uniformly across all susceptible cells within the upper airway.
A recent study has demonstrated increased accumulation of remdesivir triphosphate in human airway epithelial cultures compared to Calu-3 and Vero cell cultures [18]. This study only addressed the accumulation and efficacy in readouts from total airway epithelial cell cultures. Here we are able to observe remdesivir efficacy in all of our individually defined cell populations using actual viral replication as the readout. In addition to blocking virus replication, remdesivir may also provide benefits from blocking induction of proinflammatory cytokines while maintaining antiviral signaling. We uncovered differential expression of the entry factor TMPRSS2 across multiple cell types with higher expression in ciliated cells and near absent expression in basal cells. The protease inhibitor camostat blocked SARS-CoV-2 infection across all susceptible cells, suggesting TMPRSS2 is the dominant protease that enables infection in the upper airway. Differential expression of this and other host entry factors may explain the resistance of basal cells to infection.
Coronaviruses employ diverse mechanisms to evade detection and induction of IFN. The exact mechanisms SARS-CoV-2 uses are unclear. However, SARS-CoV-2 clearly deploys successful strategies as induction of IFN within infected cells is rare (Fig 4A and 4B). While SARS-CoV-2 is proficient at evading IFN induction in most cells, it seems to be poor at countering effects of ISGs. Consistent with several recent reports, we demonstrate that SARS-CoV-2 is highly sensitive to IFN treatment (Fig 1) [33–35]. We have previously identified ISGs in epithelial cells associated with high levels of virus replication [28]. Congruently, cells with high levels of SARS-CoV-2 reads had significantly increased expression of some of these same ISGs including GEM (Fig 4D), RIPK2, PIM3, and ODC1 (S3 Table). We have previously demonstrated a number of cell type specific ISGs during influenza A virus infection in mice [28]. Interestingly, we also identified a cell type-specific ISG (MT1F) induced during SARS-CoV-2 infection, which was primarily expressed in ciliated epithelial cells. A better understanding of both the heterogeneity of ISG responses and cell-type specific induction will help to determine transcriptional combinations that are successful in blunting infection versus those that lead to establishment of infection, robust replication, and spread to new cells.
Together the results presented here demonstrate that ciliated airway epithelial cells are the predominant cell type initially infected by SARS-CoV-2 and reveal heterogeneity within infected cells and in antiviral responses. Importantly, remdesivir is capable of reducing viral replication in all infected cell types within this culture system. We also discovered that after viral spread at 48 hpi, IFN responses and TMPRSS2 expression are the major determinants of virus tropism in susceptible cells. Our results demonstrate that infected epithelial cells at the earliest stages can produce IL-6, which potentially contributes to the cascade of inflammatory disease that occurs later during infection. Defining virus-host interactions at the single cell level reveals broad tropism and successful evasion of virus sensing, key attributes that likely allowed the zoonotic emergence of SARS-CoV-2.
Materials and methods
Primary nHTBE cell culture
nHTBE cells (Lonza Bioscience, CC-2540S) harvested from healthy individuals were expanded using Pneumacult-Ex Plus (Stem Cell Technologies) and seeded on 24 well hanging inserts, 0.4 μm pore (Corning, USA) and were maintained at ALI for 21+ days in a Pneumacult ALI growth medium (Stem Cell Technologies) at 37°C/5% CO2 in a humidified incubator. Cells were used 21+ days post ALI for experiments.
Drug treatments and infection with SARS-CoV-2
Basal layer of medium was replaced with infection media (DMEM with 2% FBS, 100 U/mL penicillin and streptomycin, 10 mM HEPES, 0.1 mM nonessential amino acids, 4 mM L-glutamate, and 1mM sodium pyruvate). Wells were left untreated or treated in the basal compartment with 100 ng of IFNβ or 0.1 μM remdesivir for 1 hr prior to infection, or with 50 μM of camostat mesylate for 2 hrs prior to infection. For infection, SARS-CoV-2 strain 2019-nCoV/USA_WA1/2020 (provided by World Reference Center for Emerging Viruses and Arboviruses at the University of Texas Medical Branch) or icSARS-CoV-2-mNG, were added directly to the apical layer of cells for 1 hr at 37°C with rocking every 10 min. The virus was removed from the apical compartment and cells were incubated at 37°C until harvest.
Flow cytometry analysis
Prior to harvest, cells were incubated in culture medium containing SiR-Tubulin (Spirochrome/Cytoskeleton Inc.) for 1 hr at 37°C. Cells were washed and trypsinized using ReagentPack Subculture Reagents (Lonza Bioscience) per the manufacturer’s instructions. Single cell suspensions were washed with 1X PBS and stained with Ghost Dye Violet 450 (Tonbo) for 30 min on ice. Cells were washed once with FACS buffer (cold HBSS supplemented with 2% bovine serum), stained with surface Abs, then fixed with 4% PFA in PBS, before flow cytometric detection on a BD LSR Fortessa (Becton Dickinson). Abs used include the following: CD66c (B6.2/CD66) (BD Biosciences); CD271 (ME20.4) (Biolegend); and TSPAN8 (458811) (Biotechne).
Sample preparation for single-cell RNA sequencing
Cells were washed and trypsinized using ReagentPack Subculture Reagents (Lonza Bioscience). Single cell suspensions were washed with PBS + 5% FBS, resuspended in 200 μL PBS + 5% FBS, passed through a 70 μm filter, and placed on ice. Cells were counted using trypan blue exclusion on a Luna II (Logos Biosystems). The appropriate numbers of cells to achieve a targeted cell input of 5,000 cells per condition were used to generate the GEMs (Gel Bead-In Emulsion). The Chromium Next GEM Single Cell 3’ Gel beads v3.1 kit (10X Genomics, Pleasanton, CA, USA) was used to create GEMs following manufacturer’s instruction. All samples and reagents were prepared and loaded into the chip and ran in the Chromium Controller for GEM generation and barcoding. GEMs generated were used for cDNA synthesis and library preparation using the Chromium Single Cell 3’ Library Kit v3.1 (10X Genomics) following the manufacturer’s instruction.
Single-cell RNA sequencing analysis
The 10x Chromium Single Cell sequencing raw data of experimental samples were reference mapped using Cell Ranger Count (version 3.0.1) for alignment to a hybrid human+SARS-CoV-2 viral genome reference (Homo_sapiens.GRCh38; MN985325). The hybrid genome reference, annotation file analysis scripts are available from our github site. The resulting raw count matrix for each experimental data set was imported into an R pipeline using Seurat 2.3.4 [42], where the filter criteria for empty droplets are minimum 1000 genes per cell, for genes that are presented in a minimum of three cells and for mitochondria genes percentage is no more than 30%. The default parameters were used for subsequent data normalization and scaling. The Seurat default method was used to select variable genes for PCA and clustering analysis, which the possible positive viral and mitochondria variable genes were excluded from final variable gene lists for all samples. The Seurat functions and utilities were used for subsequent cluster marker gene selection and results visualization. The ISG genes in the all marker list (FindAllMarkers function) of integrated analysis was ranked by their adjusted p value and ISGs within the top 30 that demonstrated disparate expression patterns are plotted with Seurat ViolinPlot function. To assess the relationship between all cells across all experimental conditions we used the FindClusters function with default parameters (except dim = 1:16 and resolution = 0.15) and the RunUMAP function with default parameters (dims = 1:16) of Seurat (version 3.2.2) was used for the integrated analysis for all five samples using only host genes. Gene ontology analysis was performed using Panther on genes in each cluster with an adjusted p value <0.01 and LogFC >0.5.
Overall cellular infection status (infected/uninfected) was determined by examining the distributions of percentage of viral counts within each cell on the log2 scale to magnify the differences at the low end. Thresholds for calling cells “infected” were set by calculating the percentage of summed viral counts within each cell over total counts, which is 0.1%. Cell types were defined using genes identified by Plasschaert et al [43]. A cutoff of at least 15 cell-type-specific genes was used to identify cell types for individual clusters.
Quantitative RT-PCR
RNA was extracted from cells using RNeasy Plus Micro Kit (Qiagen) and cDNA was generated by reverse transcription using Superscript II (Invitrogen). Genes in each sample were amplified using gene-specific primers: Human_ACE2_F 5-CGAGTGGCTAATTTGAAACCAAGAA-3; Human_ACE2_R: 5-ATTGATACGGCTCCGGGACA-3. RNA was normalized to the endogenous alpha-tubulin primer probe set: 5’-GCCTGGACCACAAGTTTGAC-3’ and 5’-TGAAATTCTGGGAGCATGAC-3’. To determine the SARS-CoV-2 genome copy number in each sample, a standard curve was generated using serial dilutions of quantitative synthetic RNA (BEI #NR-52358) and virus specific nsp14 primers: 5’-TGGGGYTTTACRGGTAACCT-3’ and 5’-AACRCGCTTAACAAAGCACTC-3’. The genome copy number in each sample was derived from this curve.
Supporting information
S1 Fig. Evaluation of SARS-CoV-2 infection heterogeneity.
(A) Frequency of cells with the indicated percent of total reads mapping to SARS-CoV-2. (B) Total RNA from cells used in Fig 2 analyzed for SARS-CoV-2 (log10(RNA copies per sample)) detected in cells samples collected at 24 or 48 hours after infection with SARS-CoV-2 in the presence or absence of remdesivir calculated by qRT-PCR. (C) Gene ontology analysis of genes in each cluster with an adjusted p value <0.01 and LogFC >0.5 using Panther. Plotted GO terms had FDR <0.01 and FE >20.
https://doi.org/10.1371/journal.ppat.1009292.s001
(TIF)
S2 Fig. Gene ontology of cell clusters defined by tSNE.
(A) Gene ontology analysis of genes in each cluster with an adjusted p value <0.01 and LogFC >0.5 using Panther. Plotted GO terms had FDR <0.01 and FE >20. (B) Total RNA from cells in Fig 2 and cells treated with or without IFNβ analyzed by qRT-PCR for ACE2.
https://doi.org/10.1371/journal.ppat.1009292.s002
(TIF)
S1 Table. List of cell type-specific genes for each cluster from the UMAP in Fig 2B.
Only cell types with at least 5 cell type-specific genes are shown.
https://doi.org/10.1371/journal.ppat.1009292.s003
(XLSX)
S2 Table. List of cell type-specific genes for each cluster from the tSNE plots for all experimental conditions in Fig 3.
Only cell types with at least 5 cell type-specific genes are shown.
https://doi.org/10.1371/journal.ppat.1009292.s004
(XLSX)
S3 Table. List of all differentially expressed genes for each cluster from the tSNE for all experimental conditions in Fig 3.
https://doi.org/10.1371/journal.ppat.1009292.s005
(XLSX)
Acknowledgments
We thank the Minnesota Super Computing Institute, the World Reference Center for Emerging Viruses and Arboviruses and the University of Minnesota Medical School Deans office for support.
References
- 1. de Jong PM, van Sterkenburg MA, Kempenaar JA, Dijkman JH, Ponec M. Serial culturing of human bronchial epithelial cells derived from biopsies. In Vitro Cell Dev Biol Anim. 1993;29A(5):379–87. pmid:7686141.
- 2. de Jong PM, van Sterkenburg MA, Hesseling SC, Kempenaar JA, Mulder AA, Mommaas AM, et al. Ciliogenesis in human bronchial epithelial cells cultured at the air-liquid interface. Am J Respir Cell Mol Biol. 1994;10(3):271–7. pmid:8117445.
- 3. Sungnak W, Huang N, Becavin C, Berg M, Queen R, Litvinukova M, et al. SARS-CoV-2 entry factors are highly expressed in nasal epithelial cells together with innate immune genes. Nat Med. 2020;26(5):681–7. pmid:32327758.
- 4. Jia HP, Look DC, Shi L, Hickey M, Pewe L, Netland J, et al. ACE2 receptor expression and severe acute respiratory syndrome coronavirus infection depend on differentiation of human airway epithelia. Journal of virology. 2005;79(23):14614–21. pmid:16282461; PubMed Central PMCID: PMC1287568.
- 5. Hoffmann M, Kleine-Weber H, Schroeder S, Kruger N, Herrler T, Erichsen S, et al. SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor. Cell. 2020;181(2):271–80 e8. pmid:32142651; PubMed Central PMCID: PMC7102627.
- 6. Walls AC, Park YJ, Tortorici MA, Wall A, McGuire AT, Veesler D. Structure, Function, and Antigenicity of the SARS-CoV-2 Spike Glycoprotein. Cell. 2020;181(2):281–92 e6. pmid:32155444; PubMed Central PMCID: PMC7102599.
- 7. Ziegler CGK, Allon SJ, Nyquist SK, Mbano IM, Miao VN, Tzouanas CN, et al. SARS-CoV-2 Receptor ACE2 Is an Interferon-Stimulated Gene in Human Airway Epithelial Cells and Is Detected in Specific Cell Subsets across Tissues. Cell. 2020;181(5):1016–35.e19. Epub 2020/04/27. pmid:32413319; PubMed Central PMCID: PMC7252096.
- 8. Monteagudo PL, Muñoz-Moreno R, Fribourg M, Potla U, Mena I, Marjanovic N, et al. Differential Modulation of Innate Immune Responses in Human Primary Cells by Influenza A Viruses Carrying Human or Avian Nonstructural Protein 1. J Virol. 2019;94(1). Epub 2019/12/12. pmid:31597767; PubMed Central PMCID: PMC6912104.
- 9. Haye K, Burmakina S, Moran T, García-Sastre A, Fernandez-Sesma A. The NS1 protein of a human influenza virus inhibits type I interferon production and the induction of antiviral responses in primary human dendritic and respiratory epithelial cells. J Virol. 2009;83(13):6849–62. Epub 2009/04/29. pmid:19403682; PubMed Central PMCID: PMC2698524.
- 10. Pruijssers AJ, George AS, Schäfer A, Leist SR, Gralinksi LE, Dinnon KH, et al. Remdesivir Inhibits SARS-CoV-2 in Human Lung Cells and Chimeric SARS-CoV Expressing the SARS-CoV-2 RNA Polymerase in Mice. Cell Rep. 2020;32(3):107940. Epub 2020/07/07. pmid:32668216; PubMed Central PMCID: PMC7340027.
- 11. Hou YJ, Okuda K, Edwards CE, Martinez DR, Asakura T, Dinnon KH 3rd, et al. SARS-CoV-2 Reverse Genetics Reveals a Variable Infection Gradient in the Respiratory Tract. Cell. 2020. pmid:32526206; PubMed Central PMCID: PMC7250779.
- 12. Chua RL, Lukassen S, Trump S, Hennig BP, Wendisch D, Pott F, et al. COVID-19 severity correlates with airway epithelium-immune cell interactions identified by single-cell analysis. Nat Biotechnol. 2020;38(8):970–9. Epub 2020/06/26. pmid:32591762.
- 13. Zhu N, Wang W, Liu Z, Liang C, Ye F, Huang B, et al. Morphogenesis and cytopathic effect of SARS-CoV-2 infection in human airway epithelial cells. Nat Commun. 2020;11(1):3910. Epub 2020/08/06. pmid:32764693; PubMed Central PMCID: PMC7413383.
- 14. Ravindra NG, Alfajaro MM, Gasque V, Wei J, Filler RB, Huston NC, et al. Single-cell longitudinal analysis of SARS-CoV-2 infection in human bronchial epithelial cells. bioRxiv. 2020. Epub 2020/05/07. pmid:32511382; PubMed Central PMCID: PMC7263511.
- 15. Wang M, Cao R, Zhang L, Yang X, Liu J, Xu M, et al. Remdesivir and chloroquine effectively inhibit the recently emerged novel coronavirus (2019-nCoV) in vitro. Cell Res. 2020;30(3):269–71. Epub 2020/02/04. pmid:32020029; PubMed Central PMCID: PMC7054408.
- 16. Sheahan TP, Sims AC, Graham RL, Menachery VD, Gralinski LE, Case JB, et al. Broad-spectrum antiviral GS-5734 inhibits both epidemic and zoonotic coronaviruses. Sci Transl Med. 2017;9(396). pmid:28659436; PubMed Central PMCID: PMC5567817.
- 17. Beigel JH, Tomashek KM, Dodd LE, Mehta AK, Zingman BS, Kalil AC, et al. Remdesivir for the Treatment of Covid-19—Final Report. N Engl J Med. 2020. Epub 2020/10/08. pmid:32445440; PubMed Central PMCID: PMC7262788.
- 18. Sheahan TP, Sims AC, Zhou S, Graham RL, Pruijssers AJ, Agostini ML, et al. An orally bioavailable broad-spectrum antiviral inhibits SARS-CoV-2 in human airway epithelial cell cultures and multiple coronaviruses in mice. Sci Transl Med. 2020;12(541). Epub 2020/04/06. pmid:32253226; PubMed Central PMCID: PMC7164393.
- 19. Agostini ML, Andres EL, Sims AC, Graham RL, Sheahan TP, Lu X, et al. Coronavirus Susceptibility to the Antiviral Remdesivir (GS-5734) Is Mediated by the Viral Polymerase and the Proofreading Exoribonuclease. mBio. 2018;9(2). Epub 2018/03/06. pmid:29511076; PubMed Central PMCID: PMC5844999.
- 20. Brown AJ, Won JJ, Graham RL, Dinnon KH, Sims AC, Feng JY, et al. Broad spectrum antiviral remdesivir inhibits human endemic and zoonotic deltacoronaviruses with a highly divergent RNA dependent RNA polymerase. Antiviral Res. 2019;169:104541. Epub 2019/06/21. pmid:31233808; PubMed Central PMCID: PMC6699884.
- 21. Kindler E, Thiel V. To sense or not to sense viral RNA—essentials of coronavirus innate immune evasion. Curr Opin Microbiol. 2014;20:69–75. Epub 2014/06/05. pmid:24908561; PubMed Central PMCID: PMC7108419.
- 22. Menachery VD, Eisfeld AJ, Schäfer A, Josset L, Sims AC, Proll S, et al. Pathogenic influenza viruses and coronaviruses utilize similar and contrasting approaches to control interferon-stimulated gene responses. mBio. 2014;5(3):e01174–14. Epub 2014/05/20. pmid:24846384; PubMed Central PMCID: PMC4030454.
- 23. Stifter SA, Bhattacharyya N, Sawyer AJ, Cootes TA, Stambas J, Doyle SE, et al. Visualizing the Selectivity and Dynamics of Interferon Signaling In Vivo. Cell Rep. 2019;29(11):3539–50.e4. pmid:31825834.
- 24. Killip MJ, Jackson D, Pérez-Cidoncha M, Fodor E, Randall RE. Single-cell studies of IFN-β promoter activation by wild-type and NS1-defective influenza A viruses. J Gen Virol. 2017;98(3):357–63. Epub 2017/03/20. pmid:27983470; PubMed Central PMCID: PMC5721924.
- 25. Kallfass C, Lienenklaus S, Weiss S, Staeheli P. Visualizing the beta interferon response in mice during infection with influenza A viruses expressing or lacking nonstructural protein 1. J Virol. 2013;87(12):6925–30. Epub 2013/04/10. pmid:23576514; PubMed Central PMCID: PMC3676098.
- 26. Russell AB, Elshina E, Kowalsky JR, Te Velthuis AJW, Bloom JD. Single-Cell Virus Sequencing of Influenza Infections That Trigger Innate Immunity. Journal of virology. 2019;93(14). pmid:31068418; PubMed Central PMCID: PMC6600203.
- 27. Sjaastad LE, Fay EJ, Fiege JK, Macchietto MG, Stone IA, Markman MW, et al. Distinct antiviral signatures revealed by the magnitude and round of influenza virus replication in vivo. Proc Natl Acad Sci U S A. 2018;115(38):9610–5. Epub 2018/09/04. pmid:30181264; PubMed Central PMCID: PMC6156629.
- 28. Fay EJ, Aron SL, Macchietto MG, Markman MW, Esser-Nobis K, Gale M, et al. Cell type- and replication stage-specific influenza virus responses in vivo. PLoS Pathog. 2020;16(8):e1008760. Epub 2020/08/13. pmid:32790753; PubMed Central PMCID: PMC7447048.
- 29. Ma JZ, Ng WC, Zappia L, Gearing LJ, Olshansky M, Pham K, et al. Unique Transcriptional Architecture in Airway Epithelial Cells and Macrophages Shapes Distinct Responses following Influenza Virus Infection. J Virol. 2019;93(6). Epub 2019/03/05. pmid:30626665; PubMed Central PMCID: PMC6401432.
- 30. Xie X, Muruato A, Lokugamage KG, Narayanan K, Zhang X, Zou J, et al. An Infectious cDNA Clone of SARS-CoV-2. Cell host & microbe. 2020;27(5):841–8 e3. pmid:32289263; PubMed Central PMCID: PMC7153529.
- 31. Bonser LR, Koh KD, Johansson K, Choksi SP, Cheng D, Liu L, et al. Flow Cytometric Analysis and Purification of Airway Epithelial Cell Subsets. Am J Respir Cell Mol Biol. 2020. Epub 2020/11/16. pmid:33196316.
- 32. Milewska A, Kula-Pacurar A, Wadas J, Suder A, Szczepanski A, Dabrowska A, et al. Replication of Severe Acute Respiratory Syndrome Coronavirus 2 in Human Respiratory Epithelium. Journal of virology. 2020;94(15). pmid:32434888; PubMed Central PMCID: PMC7375387.
- 33. Vanderheiden A, Ralfs P, Chirkova T, Upadhyay AA, Zimmerman MG, Bedoya S, et al. Type I and Type III Interferons Restrict SARS-CoV-2 Infection of Human Airway Epithelial Cultures. J Virol. 2020;94(19). Epub 2020/09/15. pmid:32699094; PubMed Central PMCID: PMC7495371.
- 34. Lokugamage KG, Hage A, de Vries M, Valero-Jimenez AM, Schindewolf C, Dittmann M, et al. Type I interferon susceptibility distinguishes SARS-CoV-2 from SARS-CoV. J Virol. 2020. Epub 2020/09/16. pmid:32938761.
- 35. Busnadiego I, Fernbach S, Pohl MO, Karakus U, Huber M, Trkola A, et al. Antiviral Activity of Type I, II, and III Interferons Counterbalances ACE2 Inducibility and Restricts SARS-CoV-2. mBio. 2020;11(5). Epub 2020/09/10. pmid:32913009; PubMed Central PMCID: PMC7484541.
- 36. Trapnell C, Cacchiarelli D, Grimsby J, Pokharel P, Li S, Morse M, et al. The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells. Nat Biotechnol. 2014;32(4):381–6. Epub 2014/03/23. pmid:24658644; PubMed Central PMCID: PMC4122333.
- 37. Cao J, Spielmann M, Qiu X, Huang X, Ibrahim DM, Hill AJ, et al. The single-cell transcriptional landscape of mammalian organogenesis. Nature. 2019;566(7745):496–502. Epub 2019/02/20. pmid:30787437; PubMed Central PMCID: PMC6434952.
- 38. Dumm RE, Fiege JK, Waring BM, Kuo CT, Langlois RA, Heaton NS. Non-lytic clearance of influenza B virus from infected cells preserves epithelial barrier function. Nature communications. 2019;10(1):779. pmid:30770807; PubMed Central PMCID: PMC6377627.
- 39. Brann DH, Tsukahara T, Weinreb C, Lipovsek M, Van den Berge K, Gong B, et al. Non-neuronal expression of SARS-CoV-2 entry genes in the olfactory system suggests mechanisms underlying COVID-19-associated anosmia. Sci Adv. 2020;6(31). Epub 2020/07/24. pmid:32937591.
- 40. Hoffmann M, Schroeder S, Kleine-Weber H, Müller MA, Drosten C, Pöhlmann S. Nafamostat Mesylate Blocks Activation of SARS-CoV-2: New Treatment Option for COVID-19. Antimicrob Agents Chemother. 2020;64(6). Epub 2020/05/21. pmid:32312781; PubMed Central PMCID: PMC7269515.
- 41. Lucas C, Wong P, Klein J, Castro TBR, Silva J, Sundaram M, et al. Longitudinal analyses reveal immunological misfiring in severe COVID-19. Nature. 2020. pmid:32717743.
- 42. Butler A, Hoffman P, Smibert P, Papalexi E, Satija R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nature biotechnology. 2018;36(5):411–20. pmid:29608179; PubMed Central PMCID: PMC6700744.
- 43. Plasschaert LW, Zilionis R, Choo-Wing R, Savova V, Knehr J, Roma G, et al. A single-cell atlas of the airway epithelium reveals the CFTR-rich pulmonary ionocyte. Nature. 2018;560(7718):377–81. pmid:30069046; PubMed Central PMCID: PMC6108322.