Multimodal CRISPR perturbations of GWAS loci associated with coronary artery disease in vascular endothelial cells

Abstract

Genome-wide association studies have identified >250 genetic variants associated with coronary artery disease (CAD), but the causal variants, genes and molecular mechanisms remain unknown at most loci. We performed pooled CRISPR screens to test the impact of sequences at or near CAD-associated genetic variants on vascular endothelial cell functions. Using CRISPR knockout, inhibition and activation, we targeted 1998 variants at 83 CAD loci to assess their effect on three adhesion proteins (E-selectin, ICAM1, VCAM1) and three key endothelial functions (nitric oxide and reactive oxygen species production, calcium signalling). At a false discovery rate ≤10%, we identified significant CRISPR perturbations near 42 variants located within 26 CAD loci. We used base editing to validate a putative causal variant in the promoter of the FES gene. Although a few of the loci include genes previously characterized in endothelial cells (e.g. AIDA, ARHGEF26, ADAMTS7), most are implicated in endothelial dysfunction for the first time. Detailed characterization of one of these new loci implicated the RNA helicase DHX38 in vascular endothelial cell senescence. While promising, our results also highlighted several limitations in using CRISPR perturbations to functionally dissect GWAS loci, including an unknown false negative rate and potential off-target effects.

Fig 1. Pooled CRISPR screens to identify CAD variants and genes that modulate vascular endothelial functions.

(A) From 92 loci associated with coronary artery disease (CAD) risk by genome-wide association studies (GWAS), we identified 2893 sentinel and linkage disequilibrium proxy variants for testing. For each of these variants, we attempted to design a maximum of five high-quality guide RNAs (sgRNAs) within a 100-bp window. In the design of the library, we also included sgRNAs that target genes essential for cell viability, as well as sgRNAs that target the coding sequence and promoter of genes that control endothelial cell functions (known genes, positive controls). (B) Number of sgRNAs per targeted variant that passed stringent quality-control filters. In total, we designed 7393 sgRNAs against 1998 CAD-associated variants (mean and median number of sgRNA per variant are 3.7 and 5, respectively). (C) Distribution of the absolute distance of the sgRNA cut-site relative to the targeted variant in base pairs (the vertical dashed line indicates mean sgRNA distance). (D) Fraction of variants at each locus that are successfully targeted by our pooled CRISPR screens. Each row represents one of the CAD loci that we tested. In green is the fraction of variants—including sentinel and LD proxies—for which we designed high-quality sgRNAs and obtained results for the endothelial function phenotypes. On average, 76% of variants at any given CAD locus are captured in the screens (vertical dashed line). (E) Most severe annotation for the 1998 CAD variants targeted by the lentiviral sgRNA libraries using ENSEMBL’s Variant Effect Predictor (VEP) module. (F) As a control step, we sequenced the plasmid library to ensure even representation of sgRNAs in the pool. Then, we produced four independent batches of lentiviruses which we used to infect teloHAEC cells that stably express Cas9, dCas9-KRAB (CRISPRi) or dCas9-VP64 (CRISPRa). Following antibiotic selection and TNFα treatment (for Cas9 and CRISPRi), we stained teloHAEC for cell surface markers (E-selectin, ICAM-1, VCAM-1) or intracellular signaling molecules (reactive oxygen species (ROS), nitric oxide (NO), calcium (Ca²⁺)). By flow cytometry, we sorted cells from the bottom and top 10 percentiles of the marker distributions, and sequenced sgRNAs found in each fraction.

Fig 2. Quality-controls of the pooled CRISPR screens for vascular endothelial cell phenotypes.

(A) Two-dimensional uniform manifold approximation and projection (UMAP) representation of 148 fluorescence-activated cell sorting (FACS) samples based on the normalized read counts of the top 10% most variable sgRNAs across all samples. (B) Density distributions of effect sizes (Beta, x-axis) across all Cas9 variants for essential genes and the rest of the sgRNA library. Positive betas indicate that sgRNA are enriched in the cell fractions when compared to the input library, while negative betas indicate a depletion of sgRNA across all samples. We observed a depletion of sgRNA targeting essential genes with all three Cas9 variants. (C-E) Rank of all control sgRNAs and targeted CAD variants in the (C) Cas9, (D) CRISPRi and (E) CRISPRa screens for three adhesion proteins: E-selectin (left), VCAM1 (middle) and ICAM1 (right). For each panel, the y-axis corresponds to the effect sizes (Beta, comparing top vs bottom FACS 10% fractions). For the Cas9 and CRISPRi experiments, we found an enrichment of sgRNAs targeting the coding and promoter sequences of genes encoding adhesion proteins in the bottom 10% cell fractions (negative Betas). In contrast, sgRNAs targeting the promoter of these genes were enriched in the top 10% cell fractions in the CRISPRa experiments. In green and blue, we highlight sgRNAs targeting coding exons and promoters, respectively. The number in front of the name of each control sgRNA indicates its rank in the corresponding analysis.

Fig 3. Discovery and validation of CRISPR perturbations that induce atheroprone vascular endothelial cell phenotypes.

(A) Heatmap of CAD-associated variants that are significant (false discovery rate (FDR) ≤10%) for at least one of six endothelial phenotypes tested in the teloHAEC pooled CRISPR screens. Each row corresponds to a combination of Cas9 variant and cellular readout, and each column corresponds to a CAD variant. For each variant, we added the name of a nearby gene to simplify locus identification, although we do not imply that these genes are causal. Dendrograms of rows and columns represent hierarchical clustering based on euclidean distance. The FDR is capped at 0.1%. (B) Validation by flow cytometry of six hits from the pooled CRISPR screens. For each validation, we used the top sgRNA from the pooled CRISPR screens to target the variant/locus with the corresponding Cas9 variant. We compared the distribution of the fluorescence intensity of the cellular markers (x-axis) between the sgRNA identified in the screens and a safe harbor negative control sgRNA. We assessed statistical significance using the Kolmogorov-Smirnov (KS) test, all validations shown are significant (KS P-value <2.2x10^-16). Validations were performed in at least three independent experiments for each sgRNA (S6 Table). For E-selectin and ICAM1, the fluorochrome is PE; for ROS, the fluorochrome is FITC.

Fig 4. Characterization of a CAD-associated regulatory variant located within an enhancer at the FURIN/FES locus.

(A) CRISPRa perturbations highlighted rs12906125 as a potential regulatory variant for FURIN and FES. The variant overlaps an ATAC-seq peak in the promoter of FES and a H3K27ac-defined enhancer that physically interacts with the FURIN promoter through chromosomal loops predicted by the ABC model applied to teloHAEC Hi-C data [7,33]. (B) Within a 2.5-Mb window, FES and FURIN are the top two differentially expressed genes when targeting rs12906125 by CRISPRa in teloHAEC. The inset plot shows the induction of both FES and FURIN expression with sgRNA_06939 when compared to the control safe harbor sgRNA. (C) teloHAEC are heterozygous (A/G) at rs12906125. We used base editing to change the genotype at rs12906125 to G/G. There was no significant difference in expression for FES and FURIN when comparing unstimulated teloHAEC with the A/G and G/G genotypes. However, upon activation with TNFα, we found that the reduction in FURIN levels was independent from the rs12906125 genotype whereas for FES, the reduction was genotype-dependent. Numbers above the bars are Student’s t-test P-values. We tested at least six clones of each genotype.

Fig 5. Disruption of DHX38 induces vascular endothelial cell senescence.

(A) Perturbations with the Cas9 nuclease highlighted two synonymous variants (rs2074626, rs2240243) in the DHX38 gene for several endothelial phenotypes. DHX38 is located downstream of the HP and HPR genes, which have previously been associated with LDL-C levels. However, the CAD and LDL-C GWAS signals are distinct based on co-localization analyses (posterior probability for two independent association signals (H3) = 80.9%). (B) Gene-set enrichment analysis results for differentially expressed genes identified by RNA-seq in teloHAEC between a sgRNA targeting a DHX38 coding exon and a safe harbor negative control sgRNA. Only pathways with a Benjamini-Hochberg-corrected P-value <0.05 and normalized enrichment scores (NES) <-1 or >1 are shown. (C) Experimental design for the characterization of DHX38 using the fluorescent marker CRIMSON in place of an antibiotic resistance gene. We did all experiments in teloHAEC that stably express Cas9. We monitored the impact of a DHX38 sgRNA on cell proliferation, indel induction, gene expression and senescence-associated β-galactosidase (SA-βGal) activity. (D) Comparison of endothelial cell proliferation between teloHAEC with a DHX38 sgRNA or a safe harbor negative control sgRNA. The differences in the number of CRIMSON⁺ cells were not significant two or four days post-infection. However, there were 27% less CRIMSON⁺ cells with DHX38 sgRNA relative to the safe harbor control at seven days post-infection (Student’s t-test P-value = 7.3x10^-8). Results are mean ± standard deviation for 6 replicates for safe harbor and three replicates for two DHX38 targeting sgRNA. (E) Quantification of DHX38 indels by tracking of indel by decomposition (TIDE) analysis. As expected, we found no indels in the CRIMSON^- cells (S8 Table). However, in CRIMSON⁺ cells that received a DHX38 sgRNA, we found indels with an average frequency of 15%, 42% and 40% at day 2, 4 and 7, respectively. Results are mean ± standard deviation for 6 replicates for safe harbor and three replicates for two DHX38 targeting sgRNA. (F) Expression levels of DHX38 and CDKN1A in CRIMSON^- and CRIMSON⁺ teloHAEC that have received a sgRNA that targets DHX38 or a safe harbor region (negative control). There were no significant differences in DHX38 expression levels at day 2. However, at day 4 and 7, DHX38 was significantly down-regulated and CDKN1A was significantly up-regulated in CRIMSON⁺ cells that received the DHX38 sgRNA. N.S., not significant. We provide Student’s t-test P-values when P<0.05. Bars are mean normalized expression and error bars represent one standard deviation. (G) Quantification of senescent teloHAEC by flow cytometry using senescence-associated β-galactosidase (SA-βGal) staining. At day 4 and 7 post-infection, there were significantly more senescent cells in the CRIMSON⁺ DHX38 sgRNA experiment than in the CRIMSON^- cells or in the CRIMSON⁺ cells that received the safe harbor sgRNA. We used the DNA damaging agent etoposide as a positive control to induce senescence. N.S., not significant. We provide Student’s t-test P-values when P<0.05. Results are mean percentage SA-βGal⁺ teloHAEC and error bars represent one standard deviation.

Fig 6. Validated CRISPRa effect at the CNNM2 and CCDC92/ZNF664 loci do not nominate candidate causal CAD genes.

(A) Locus view for the CAD locus with nearby gene CNNM2. We provide the position of the sentinel CAD variant (rs11191416) and the putative functional variant identified in the pooled CRISPR screen (rs78260931). The LD proxies and sgRNAs tested are also shown. ATAC-seq and RNA-seq data in resting teloHAEC are from ref. [7]. (B) Locus view for the CAD locus with nearby genes ZNF664 and CCDC92. We provide the position of the sentinel CAD variant (rs11057401) and the functional variant identified in the pooled CRISPR screen (rs12311848). The LD proxies and sgRNAs tested are also shown. ATAC-seq and RNA-seq data in resting teloHAEC are from ref. [7]. (C) Uniform manifold approximation projection (UMAP) for 11,756 cells from human right coronary arteries analyzed by single-cell RNA-sequencing [36]. We color-coded cells based on the level of expression of candidate causal CAD genes identified and characterized in this study. We used the expression of the endothelial cell marker gene CDH5 (encoding VE-Cadherin) to identify endothelial cells (circle in top left panel). All five candidate genes are expressed in human vascular endothelial cells from coronary arteries.