Genome-wide CRISPR Screens Reveal Host Factors Critical for SARS-CoV-2 Infection

Abstract

Identification of host genes essential for severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) infection may reveal novel therapeutic targets and inform our understanding of coronavirus disease 2019 (COVID-19) pathogenesis. Here we performed genome-wide CRISPR screens in Vero-E6 cells with SARS-CoV-2, Middle East respiratory syndrome CoV (MERS-CoV), bat CoV HKU5 expressing the SARS-CoV-1 spike, and vesicular stomatitis virus (VSV) expressing the SARS-CoV-2 spike. We identified known SARS-CoV-2 host factors, including the receptor ACE2 and protease Cathepsin L. We additionally discovered pro-viral genes and pathways, including HMGB1 and the SWI/SNF chromatin remodeling complex, that are SARS lineage and pan-coronavirus specific, respectively. We show that HMGB1 regulates ACE2 expression and is critical for entry of SARS-CoV-2, SARS-CoV-1, and NL63. We also show that small-molecule antagonists of identified gene products inhibited SARS-CoV-2 infection in monkey and human cells, demonstrating the conserved role of these genetic hits across species. This identifies potential therapeutic targets for SARS-CoV-2 and reveals SARS lineage-specific and pan-CoV host factors that regulate susceptibility to highly pathogenic CoVs.

Keywords: COVID-19; CRISPR screen; Epigenetics; HMGB1; MERS-CoV; Middle East Respiratory Syndrome; SARS-CoV-2; SWI/SNF complex; Severe Acute Respiratory Syndrome.

Graphical abstract

Figure 1

Genome-wide CRISPR Screens Identify Genes Critical for CoV-Induced Cell Death (A) Schematic of the pooled screen. Vero-E6-Cas9 cells transduced with the genome-wide C. sabaeus library received mock treatment or were challenged with SARS-CoV-2, rcVSV-SARS-CoV-2-S, HKU5-SARS-CoV-1-S, MERS-CoV, or MERS-CoV T1015N. Surviving cells from each virus infection were isolated, and the sgRNA sequences were amplified by PCR and sequenced. (B) Volcano plot showing top genes conferring resistance and sensitivity to SARS-CoV-2. The gene-level Z score and −log10 (FDR) were calculated using the mean of the five Cas9-v2 conditions. Non-targeting control sgRNAs were grouped randomly into sets of 4 to serve as “dummy” genes and are shown in green. (C) Heatmaps of the top gene hits for SARS-CoV-2 resistance (20) and sensitization (20), ranked by mean Z score. The top 5 hits for MERS-CoV are also included and indicated by an asterisk. ARID1A was a top resistance gene for SARS-CoV-2 and MERS-CoV. (D–F) Correlation between gene enrichment in SARS-CoV-2 and rcVSV-SARS-CoV-2-S (D), HKU5-SARS-CoV-1-S (E), and MERS-CoV (F) screens. R, pearson correlation. (G) Venn diagram of the top 100 pro-viral genes from SARS-CoV-2, rcVSV-SARS-CoV-2-S, HKU5-SARS-CoV-1-S, and MERS-CoV screens.

Figure S1

Quality Control Metrics for CRISPR Screen, Related to STAR Methods (A) Correlation matrix depicting the Pearson correlation between the guide-level log-fold change values relative to the plasmid DNA. Cells were cultured in DMEM with 2% FBS (D2), 5% FBS (D5), 10% FBS (D10), plated at 2.5 × 10⁶ or 5.0 × 10⁶ cells per T150 flask and infected at a MOI 0.1 (hi) or MOI 0.01 (lo). (B) Receiver-operator characteristic (ROC) curve for the recovery of guides targeting essential genes in the mock-treated condition of the Cas9-v1 and Cas9-v2 screens. True positives are n = 1,528 essential genes (n = 6,178 guides); true negative genes are n = 622 non-essential genes (n = 2,504 guides). We mapped essential and non-essential genes, which were derived for human cell lines, to the African green monkey genome simply by matching gene symbols. AUC = area under curve. (C) Correlation between gene enrichment in Cas9-v1 and Cas9-v2 screens. Pearson correlation is reported. (D-E) GFP-based Cas9 activity assay in Vero-E6 cells stably expressing either Cas9-v1 (D) or Cas9-v2 (E). The pXPR_047 construct expresses GFP and an sgRNA targeting GFP; therefore, cells without Cas9 activity will express GFP, whereas cells with high Cas9 activity will knock out GFP and resemble parental cells. (F) Approach to calculate residuals from log-fold change data, using ACE2 and the 5% FBS, 5 × 10⁶ cells/flask, MOI 0.1 condition as an example. A natural cubic spline with four degrees of freedom is shown in blue, and a residual for each sgRNA is calculated to be the vertical distance from the fit spline.

Figure S2

A Genome-wide CRISPR Screen Identifies Genes Critical for SARS-CoV-2-Induced Cell Death, Related to Figure 1 (A) Performance of individual sgRNAs targeting ACE2, SMARCA4, CTSL, and TMPRSS2. The mean residual across the five Cas9-v2 conditions is plotted for the full library (top) and for the 4 guide RNAs targeting each gene. (B)Heatmaps of the top 25 gene hits for resistance and sensitivity, ranked by mean z-score in the Cas9-v2 conditions. Genes that are included in one of the gene sets labeled in (Figure 2A) are colored accordingly. Condition A: Cas9-v2 D5 (DMEM+5%FBS) 2.5 × 10⁶ cells/flask MOI 0.1; B: Cas9-v2 D5 5 × 10⁶ cells/flask MO 0.1; C: Cas9-v2 D2 (DMEM+2%FBS) 5 × 10⁶ cells/flask MOI 0.1; D: Cas9-v2 D10 (DMEM+10%FBS) 5 × 10⁶ cells/flask MOI 0.1; E: Cas9-v2 D5 2.5 × 10⁶ cells/flask MOI 0.01. (C) Nodes represent significantly enriched gene sets. The size of each gene set is proportional to its mean absolute z-score. Gene sets are colored by the direction in which they score. Edges represent significant overlap between gene sets. The transparency of each edge is proportional to the fraction of genes shared by two gene sets. Gene sets were clustered using the infomap algorithm and the most central set by PageRank is labeled for each cluster. The Fruchterman–Reingold algorithm was used to lay out the network.

Figure S3

Comparison of All Viruses from Genome-wide CRISPR Screens, Related to Figure 1 Comparison of gene enrichment of (A) SARS-CoV-2 relative to MERS-CoV T1015N, (B) rcVSV-SARS-CoV-2-S relative to HKU5-SARS-CoV-1-S, (C) rcVSV-SARS-CoV-2-S relative to MERS-CoV WT, (D) rcVSV-SARS-CoV-2-S relative to MERS-CoV T1015N. (E) HKU5-SARS-CoV-1-S relative to MERS-CoV WT, (F) HKU5-SARS-CoV-1-S relative to MERS-CoV T1015N, and (G) MERS-CoV WT relative to MERS-CoV T1015N. Pearson correlation is reported.

Figure 2

Performance of Genes in the Top Gene Sets (A) The top three gene sets that score in the positive direction (resistance) and top gene set that scores in the negative direction (sensitization) or both, filtered for gene sets with at least five genes, and that are most central to a given module (Figure S2C) and then ranked by mean absolute Z score. The number of genes in each set is indicated in parentheses. (B) For each gene in the “SWI/SNF complex” gene set from STRING, the Z score in each virus screen is shown. (C–F) Similarly, the genes in the gene sets (C) “RUNX3 regulates CDKN1A transcription” from Reactome, (D) “Cystatin, and endolysosome lumen” from STRING, (E) “Viral translation” from GO, and (F) “NURF complex” from GO.

Figure 3

CRISPR Subpool Screens Validate Primary Genome-wide Screens and Demonstrate the Specificity of Hits for CoVs A CRISPR subpool was generated with 10 sgRNAs per gene for each of the top 250 and bottom 250 genes from the SARS-CoV-2 genome-wide screen along with non-targeting controls and other genes of interest, including DPP4. (A) Correlation between gene enrichment in primary genome-wide and secondary subpool SARS-CoV-2 subscreens. Pearson correlation is reported. (B) Correlation matrix depicting the Pearson correlation between the guide-level log-fold change values relative to the plasmid DNA for the 13 subpool screens with the indicated viruses. All viruses were screened in duplicate (#1 and #2), except IAV-WSN. VSV was also screened, but no cells survived infection. (C) Principal-component analysis (PCA) plot of all viruses reveals clustering and overlap of gene hits among SARS-CoV-2, rcVSV-SARS-CoV-2-S, HKU5-SARS-CoV-1-S, and MERS-CoV WT and the T1015N cluster. IAV/WSN/1933 (IAV-WSN) and encephalomyocarditis virus (EMCV) are outliers among the CoV screens. (D–H) Comparison of gene enrichment in SARS-CoV-2 relative to rcVSV-SARS-CoV-2-S (D), HKU5-SARS-CoV-1-S (E), MERS-CoV (F), IAV-WSN (G), and EMCV (H). Pearson correlation is reported. (I) We generated a CRISPR subpool targeting 32 genes (inclusive of control genes) in the human lung cancer cell line Calu-3. Gene enrichment from the primary SARS-CoV-2 screen correlates with results from Calu-3 cells. Pearson correlation is reported.

Figure S4

Comparison of Secondary CRISPR Subpool Screens, Related to Figure 3 (A) Heatmap depicting genes with a z-score > 10 in any of the secondary subpool screens. (B) Heatmap showing genes involved in the SWI/SNF chromatin remodeling complex. (C) Heatmap showing genes in the “Runx3 regulates CDKN1A transcription” pathway. (D) Heatmap showing genes in the HUCA histone H3.3 chaperone complex.

Figure 4

Arrayed Validation of 18 Resistance and 7 Sensitization Hit Genes (A) Performance in the pooled screen of sgRNAs targeting the 25 genes selected for further validation. The mean residual across the five Cas9-v2 conditions is plotted for the full library (top) and for the 3–4 sgRNAs targeting each gene. Genes that scored as resistance hits are shown in red; genes that scored as sensitization hits are shown in blue. The dashed line indicates a residual of 0. (B) 42 unique sgRNAs targeting 25 genes were introduced into Vero-E6-Cas9-v2 cells. SARS-CoV-2 was added at MOI 0.2, and cell viability was measured at 3 dpi. (C) Western blot for ACE2, SMARCA4, KDM6A, and SMAD3 expression in control and the respective gene-disrupted Vero-E6 cells. (D) Z scores from the genome-wide CRISPR screen correlate with cell viability of individually disrupted genes. Genes with multiple sgRNAs from (B) are averaged to generate one point per gene Data were analyzed by one-way ANOVA with Tukey’s multiple comparison test. Shown are means ± SEM. ns, not statistically significant; ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001.

Figure 5

Small Molecules Protect Cells from SARS-CoV-2-Induced Cell Death (A–C) Vero-E6 cells were pretreated with the indicated concentrations of the Cathepsin L inhibitor Calpain inhibitor III (A), the SMARCA4 inhibitor PFI-3 (B), or the SMAD3 inhibitor SIS3 (C) for 48 h and then infected with SARS-CoV-2 at a MOI of 0.2. Cell viability was measured at 3 dpi and compared with mock-infected controls. Red, infected; blue, mock-infected. (D and E) Vero-E6 cells were pretreated with 10 μM Calpain inhibitor III, PFI-3, or SIS3 for 48 h and then infected with icSARS-CoV-2 mNG at a MOI of 1. Infected cell frequencies were measured by mNeonGreen expression at 2 dpi. Scale bars, 300 μm. (F–H) Vero-E6 (F), Huh7.5 (G), and Calu-3 (H) cells were pretreated with 10 μM SIS3 and 40 μM PFI-3 for 48 h and then infected with SARS-CoV-2 at a MOI of 0.1. Virus production, as measured by plaque-forming units (PFU) per milliliter, was determined by plaque assay. LOD, limit of detection. Data were analyzed by one-way ANOVA with Tukey’s multiple comparison test. Shown are means ± SEM. ^∗p < 0.05, ^∗∗p < 0.01, ^∗∗∗p < 0.001.

Figure 6

HMGB1 Is a Novel Regulator of ACE2 (A) Performance of individual guide RNAs targeting LOC103214541 (HMGB1-like). The mean residual across the five Cas9-v2 conditions is plotted for the full library (top) and for the 3 guide RNAs targeting that gene. (B) Western blot for HMGB1 expression in control and HMGB1-disrupted Vero-E6, Huh7.5, and Calu-3 cells. (C) Control and HMGB1-disrupted Vero-E6, Huh7.5, and Calu-3 cells were infected with SARS-CoV-2 at a MOI of 0.2. Cell viability relative to an uninfected control was measured 3 dpi (Vero-E6 and Calu-3) or 4 dpi (Huh7.5) with CellTiter Glo. (D) Vero-E6 cells were infected with SARS-CoV-2 at a MOI of 0.1. Virus production was measured by plaque assay. (E) Correlation between CRISPR screen Z score and gene expression in control and HMGB1-disrupted Vero-E6 cells reveals downregulation of ACE2 in HMGB1-disrupted cells. (F) qPCR and western blot were performed in HMGB1 knockout and complemented Vero-E6 cells. (G) Genome tracks of RNA sequencing (RNA-seq), ChIP-seq for H3K27ac, and ATAC-seq at the ACE2 locus in control and HMGB1-disrupted Vero-E6 cells. The p values for ChIP-seq and ATAC-seq are for the genomic region indicated by the black bar below the tracks. (H) HMGB1 knockout and complemented Vero-E6 cells were infected with VSV pseudoparticles (VSVpp): VSVpp-SARS-CoV-1-S, VSVpp-SARS-CoV-2-S, VSVpp-NL63-S, VSVpp-MERS-S, and VSVpp-VSV-G. Luciferase relative to a VSVpp-VSV-G control was measured 1 dpi. Data were analyzed by one-way ANOVA with Tukey’s multiple comparison test. Shown are means ± SEM. ^∗p < 0.05, ^∗∗p < 0.01, ^∗∗∗p < 0.001.

Figure S5

HMGB1 Acts Cell-Intrinsically to Regulate Susceptibility to SARS-CoV-2 Infection, Related to Figure 6 (A) Vero-E6 cells were mock-treated or infected with SARS-CoV-2 at a MOI of 1 for 24 hours before cell fractionation was performed. (B-C) Infection resulted in release of HMGB1 protein in the supernatant which was quantified by ELISA from Vero-E6 (B) and Huh7.5 (C) cells infected with SARS-CoV-2 for the indicated times. (D) Vero-E6 cells were pre-treated with the indicated concentration of recombinant HMGB1 (rHMGB1) for 24 hours and then infected with SARS-CoV-2 at a MOI of 0.2. Cell viability was measured at 3 dpi and compared to mock infected controls. (E) Vero-E6 cells were pre-treated with rHMGB1 for 24 hours and then infected with icSARS-CoV-2 mNG at a MOI of 1. Infected cell frequencies were measured by mNeonGreen expression at 1 dpi. (F) Vero-E6 cells were pre-treated with the indicated concentration of rHMGB1 for 24 hours and then infected with VSVpp-SARS-CoV-2-S and VSVpp-VSVG pseudovirus. Luciferase relative to VSVG control was measured at 1 dpi. Data were analyzed by one-way ANOVA with Tukey’s multiple comparison test. Shown are means ± SEM ns, not statistically significant; ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001.

Figure S6

Effects of HMGB1 Loss on Chromatin States across the Vero-E6 Genome, Related to Figure 6 (A) Volcano plot for RNA sequencing of control and HMGB1 disrupted cells. The x axis shows log2 fold-change and the y axis shows −log10 of the adjusted P value (adj. P) as calculated by DESeq2. (B) Top gene sets, which significantly enriched in the upregulated, downregulated or both of differentially expression genes (fold change > 1.5 and p < 0.05) from GO. (C-D)Volcano plots for ATAC-seq (C) and H3K27ac ChIP-seq (D) of control and HMGB1 disrupted cells. Fold change and adjusted p value for each called peak was calculated by DESeq2. (E) Correlation between changes in overlapping ATAC-seq and H3K27ac ChIP-seq peaks upon HMGB1 disruption. Dashed lines represent p = 0.01. Data were analyzed by one-way ANOVA with Tukey’s multiple comparison test. Shown are means ± SEM ns, not statistically significant; ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001.