The impact of common and rare genetic variants on bradyarrhythmia development

Figures

• 

  Fig. 1. Study design. a, Typical anatomical regions in which conduction tissue affected by SND and DCD are localized in the heart. A dual-chamber PM is also demonstrated. Sample sizes are shown by common variant GWAS and rare variant association tests (RVAT) for all three outcomes. b, Overview of common and rare variant analyses for SND, DCD and PM implantation. Common-variant GWAS were performed in ten collaborating studies with genotyped and imputed data, and rare variant burden testing was performed for the same phenotypes in two studies with whole-exome sequencing. The results were combined in meta-analyses, including up to 1.3 million individuals for GWAS and 471k individuals for rare variant association testing. For common variant loci reaching genome-wide significance, follow-up evaluations included analyses of cardiac gene expression profiles, predicted transcriptomes, pleiotropic associations, genetic correlations, Mendelian randomization (MR), polygenic risk scores (PRS) derivation and relations with risk of PM implantation and phenome-wide association study (PheWAS). The interaction of rare variants and polygenic risk was further evaluated.

• 

  Fig. 2. Manhattan plot for bradyarrhythmias. a–c, GWAS results are shown separately for SND (a; n = 1,258,554), DCD (b; n = 1,314,957) and PM implantation (c; n = 1,304,231). Two-sided P values (on −log10 scale) for each association test between variants and bradyarrhythmias from fixed-effect meta-analyses of multi-ancestry individuals are shown on the y axis. Genome-wide significant association loci (P < 5 × 10−8 after Bonferroni correction; dashed line) are annotated with the name of the gene closest to the index variant.

• 

  Fig. 3. Associations of polygenic risk scores for bradyarrhythmias with outcomes in unrelated UKBB participants. a, Cumulative incidence of PM implantation in UKBB participants stratified by PRSs for SND, DCD and PM implantation. The shaded area surrounding each line refers to the two-sided 95% CI. PRSs were constructed using the clumping and thresholding method separately for each phenotype (nvariants for SND PRS = 28, nvariants for DCD = 57 and nvariants for PM = 51). A total of 327,702 unrelated participants without a history of PM implantation at study enrollment were included in the analyses. Participants were stratified into three groups based on the tertiles of residuals of each PRS after adjustment by the first ten PCs. The statistical significance of the associations between tertiles of each PRS and PM incidence was evaluated using analysis of variance comparing a Cox proportional hazards model with only sex, age and genotyping array as predictors and Cox proportional hazards models with each PRS tertile as an additional predictor. The median follow-up time was 11.2 (Q1–Q3, 10.5–11.7) years. b, Associations of PRS residuals (after adjustment by the first ten PCs) for…

• 

  Fig. 4. Causal links of bradyarrhythmias with potential causes. a–c, In each panel, a bubble plot shows results from the MR screen modeling the given bradyarrhythmia as the outcome, using the weighted median method. a, Results for SND. b, Results for DCD. c, Results for PM. The y axis represents the signed −log10 of the P value from the MR analysis, where each bubble represents a different disease/trait (exposure or outcome), and −log10(P) is signed by the direction of the MR effect estimate; the diseases/traits are ordered by their signed −log10(P) from high (left) to low (right). All plotted P values are two-sided. The full red and blue lines represent the Bonferroni-corrected significance level (P < 0.05/(74 × 2) tests) for a given bradyarrhythmia, while the dotted lines represent P < 0.01. Traits/diseases reaching two-sided Bonferroni significance in the weighted median screen are annotated with their names; traits/diseases annotated in bold with a red or blue color also passed all subsequent filters and reached significance using CAUSE (Methods), while traits/diseases annotated in black passed MR–Egger intercept filtering (two-sided P > 0.025) but showed only a nominal association…

• 

  Fig. 5. Rare-variant association tests. a–c, Gene-level results from rare variant burden tests are shown separately for SND (a; n = 460,813), DCD (b; n = 471,469) and PM implantation (c; n = 464,692). Two-sided Cauchy P values (on −log10 scale) from fixed-effect meta-analysis for each gene are shown on the y axis. Dashed lines indicate exome-wide significance thresholds (P < 2.7 × 10−6 after Bonferroni correction). Significant associations are presented in red, and gene names are annotated.

• 

  Fig. 6. PM implantations in carriers of protein-disrupting variants among 305,633 unrelated UKBB participants. a, Proportion of unrelated UKBB participants with PM implantations at study enrollment or during follow-up in participants who were carriers of a protein-disrupting variant (a LOF variant or a missense variant predicted to be pathogenic by at least 80% of bioinformatics tools) in any of the five genes (LMNA, SCN5A, MYBPC3, SMAD6 and TTN) that were significantly associated with at least one bradyarrhythmia phenotype in rare variant association tests. b, Proportion of participants with PMs in participants who were carriers of protein-disrupting variants in each gene mentioned above. c, Proportion of participants with PMs among carriers of protein-disrupting variants (right) and noncarriers (left), stratified by the tertiles of a PRS for PM implantation adjusted by the first ten PCs. The error bars in a–c correspond to binomial 95% CIs calculated using the Agresti–Coull method. Two-sided P values for the PRS were derived with logistic regression using PM implantation as the outcome and the PRS tertiles, sex and age at study enrollment as predictors. Numbers in parentheses are the…

• 

  Extended Data Fig. 1. Heritability enrichment for bradyarrhythmias in nine major cell types from human heart. Results of stratified LD score regression (s-LDSC) on the combined major cell types in heart. One-sided P values were derived by testing for GWAS heritability enrichment near cell-type-specific genes controlling for an annotation based on SNPs near any gene and 52 additional baseline annotations from s-LDSC. Dashed lines show statistical (red, with Bonferroni correction) and nominal (blue, P = 0.05) significant P-value thresholds.

• 

  Extended Data Fig. 2. The layered Cauchy combination pipeline in rare variant studies of bradyarrhythmia, as applied for each gene. a, Birds-eye view of the pipeline as applied to each transcript of a given gene. As shown on the far left, we defined 9 ‘transcripts’ that were defined as (1) all exons; (2) canonical gene transcripts as determined by Ensemble; (3) aortic expressed transcripts, formed by variants with aortic-specific pext values ≥ 0.8; (4) atrial appendage-expressed transcripts; (5) left ventricle-expressed transcripts; (6) coronary artery-expressed transcripts; (7) tibial artery-expressed transcripts; (8) whole blood-expressed transcripts; and (9) mean transcript expression across GTEx tissues. For each of the defined gene transcripts, burden-testing P values based on LOF variant masks were combined into a single P value using the Cauchy distribution test; burden-testing P values based on various missense masks were combined into a single P value using the Cauchy distribution test; and all LOF + missense burden P values were combined into a single P value using the Cauchy distribution test. Then the LOF, missense and LOF + missense P values were combined into a single transcript P value…