Skip main navigation
Close Drawer MenuOpen Drawer Menu

Common and Rare Coding Genetic Variation Underlying the Electrocardiographic PR Interval

Originally published18 Jun 2018https://doi.org/10.1161/CIRCGEN.117.002037Circulation: Genomic and Precision Medicine. 2018;11:e002037

    Abstract

    Background:

    Electrical conduction from the cardiac sinoatrial node to the ventricles is critical for normal heart function. Genome-wide association studies have identified more than a dozen common genetic loci that are associated with PR interval. However, it is unclear whether rare and low-frequency variants also contribute to PR interval heritability.

    Methods:

    We performed large-scale meta-analyses of the PR interval that included 83 367 participants of European ancestry and 9436 of African ancestry. We examined both common and rare variants associated with the PR interval.

    Results:

    We identified 31 genetic loci that were significantly associated with PR interval after Bonferroni correction (P<1.2×10−6), including 11 novel loci that have not been reported previously. Many of these loci are involved in heart morphogenesis. In gene-based analysis, we found that multiple rare variants at MYH6 (P=5.9×10−11) and SCN5A (P=1.1×10−7) were associated with PR interval. SCN5A locus also was implicated in the common variant analysis, whereas MYH6 was a novel locus.

    Conclusions:

    We identified common variants at 11 novel loci and rare variants within 2 gene regions that were significantly associated with PR interval. Our findings provide novel insights to the current understanding of atrioventricular conduction, which is critical for cardiac activity and an important determinant of health.

    Introduction

    See Editorial by Hedley et al

    Clinical Perspective

    The duration of PR interval is an important biomarker of the cardiac conduction system. Increasing evidence suggests that cardiac conduction measurements including PR interval are heritable. It is thus interesting to understand the biological and potential clinical implications of genetic variation underlying cardiac conduction. We performed a large-scale meta-analysis of PR interval that included 83 367 participants of European ancestry and 9436 of African ancestry using the Illumina exome chip. Thirty-one genetic loci were significantly associated with PR interval after Bonferroni correction, including 11 loci that have not been previously reported. Our findings provide new insights to the current understanding of atrioventricular conduction, which is critical for cardiac activity and an important determinant of health.

    Electrical conduction from the cardiac sinoatrial node to the ventricles is critical for normal heart function. Abnormalities of atrioventricular conduction can cause significant morbidity, and have been associated with atrial fibrillation (AF),1,2 cardiac malformations, and sudden death.3,4 Conduction from the sinus node through the atria, atrioventricular node, and His-Purkinje fibers is readily evaluated from surface ECG, by measurement of the duration of PR interval. Despite the critical role that the cardiac conduction system plays in cardiac physiology and disease, the formation and regulation of the conduction system remains incompletely understood.

    Recent data indicate that cardiac conduction measurements are heritable5–7 and have a genetic basis.8–11 To date, genetic studies of PR interval have been relatively modest-sized largely European-ancestry samples and have implicated cardiac-expressed ion channels, cardiac developmental transcription factors, signaling molecules, and novel pathways not previously known to be involved in cardiac conduction processes. Nevertheless, existing studies have focused on the role of common and predominantly noncoding genetic variants, which account for only a modest proportion of trait heritability.6

    To better understand the biological and potential clinical implications of genetic variation underlying cardiac conduction, there is a need to examine both common and rare variations underlying atrioventricular conduction in large, well-powered, multiethnic studies. Moreover, assessment of genetic variation that alters protein coding has the potential to more directly implicate genes involved in processes critical to cardiac conduction. We therefore sought to examine PR interval duration in relation to predominantly coding genetic variants, in large, multi-ethnic analyses using the exome chip.

    Methods

    The data, analytic methods, and study materials will be made available to other researchers for purposes of reproducing the results, subject to Data Use/Sharing Agreements adopted by individual participating cohorts. The summary results from the current article are available at the Broad Cardiovascular Disease Knowledge Portal (www.broadcvdi.org).

    Study Participants

    The current project included participants of European ancestry from 22 studies: AGES (Age, Gene/Environment Susceptibility Study); ARIC study (Atherosclerosis Risk in Communities); BRIGHT (British Genetics of Hypertension); Massachusetts General Hospital CAMP (Cardiology and Metabolic Patient) cohort; CHS (Cardiovascular Health Study); ERF (Erasmus Rucphen Family Study); FHS (Framingham Heart Study); GOCHA (Genes for Cerebral Hemorrhage on Anticoagulation); GRAPHIC (Genetic Regulation of Arterial Pressure in Humans in the Community); INTER99 (The Inter99 Study); KORA (Cooperative Health Research in the Region Augsburg); KORCULA (CROATIA-Korcula); LifeLines (LifeLines Cohort Study); MESA (Multi-Ethnic Study of Atherosclerosis); NEO (The Netherlands Epidemiology of Obesity); RS (Rotterdam Study); GS:SFHS (Generation Scotland: Scottish Family Health Study); SHIP (Study of Health in Pomerania); TwinsUK; UHP (Utrecht Health Project); WHI (Women’s Health Initiative); and YFS (Young Finns Study).

    In addition, we included participants of African ancestry from 5 studies. These studies included ARIC, CHS, JHS (Jackson Heart Study), MESA, and WHI.

    Institutional Review Boards or Ethics Committees approved study procedures at each contributing site. All participants provided written informed consent to participate in genetic research.

    Measurement of PR Interval

    PR interval duration, in milliseconds, was measured from the onset of the P wave to the onset of the QRS interval for each cohort. The following exclusions were applied: extreme PR values (≤80 or ≥320 ms); second- or third-degree heart block; AF on baseline ECG; history of myocardial infarction, heart failure, or Wolff–Parkinson–White syndrome; pacemaker placement; use of class I or III blocking medications (ATC code prefix C01B); digoxin use (ATC code C01AA05); or pregnancy.

    Genotyping

    Genotyping was performed independently in each study using the Illumina Human Exome BeadChip (v1.0, 1.1, or 1.2). Data were called and cleaned according to Cohorts for Heart and Aging Research in Genomic Epidemiology Exome Chip best practices.12 Detailed information for each study on genotyping platforms, variant calling, and quality control metrics is shown in Table I in the Data Supplement. All studies used the same set of reference alleles to recode variants to ensure consistency.

    Statistical Analyses

    Before association analysis, PR interval was first adjusted for covariates by taking residuals from a linear regression of PR on age, sex, height, body mass index, and RR interval. Each cohort additionally adjusted as necessary for cohort-specific variables, such as clinic sites, family structure, and population structure. To reduce sensitivity to extreme PR values, the residuals were inverse normal transformed and used as the outcome for association testing.

    Because single-marker–based analyses typically have low power to identify associations between rare variants and traits, we separated the analysis for common and rare variants based on minor allele frequency (MAF). Common variants were defined as those with MAF ≥1%, and the remaining variants were defined as rare variants (MAF<1%). For each of the common variants, we evaluated its association with the transformed PR interval and accounted for multiple testing by Bonferroni correction (P<0.05/42 075=1.2×10−6). For the rare variants, we restricted analyses to nonsynonymous or splicing variants with MAF <1% because such variants are more likely to be functional than synonymous or more common variants. As we expect some rare variants may act in the same or opposite directions even in the same gene region,13 we used a modified version of the Sequence Kernel Association Test,14 which avoids problems of signals canceling out each other in burden test results. Many gene regions had few or no rare nonsynonymous or splicing variants. Monomorphic variants from each study also were reported in the cohort-level results as they were used for the cumulative MAF computations in gene-based tests. Gene regions with a cumulative MAF of rare variants <1% were excluded, resulting in 5761 gene regions that were tested (Results). Therefore, Bonferroni-corrected significance threshold for our gene-based tests was P<0.05/5761=8.7×10−6. In secondary analyses, we limited the analysis to damaging variants, defined as nonsense variants or variants predicted to be damaging by PolyPhen-215 or SIFT.16

    Analyses were performed using the prepScores function of the seqMeta R package. Family-based studies implemented the kins option in prepScores to specify kinship matrices. Each study provided single-variant Z statistics from score tests, as well as genotype covariance matrices, which were then combined by fixed-effects meta-analysis. The heterogeneity across studies was assessed by the Cochran Q, which is a nonparametric statistical test defined as the weighted sum of squared differences between individual study effects and the pooled effect. We performed both race-stratified and race-combined meta-analyses, and the race-combined results were used for the remaining sections unless stated otherwise.

    Comparison With Genetic Loci Associated With AF and P-Wave Indices

    We also compared genetic loci associated with PR interval with those associated with AF and P-wave indices (PWI) to see whether there are any shared genetic mechanisms. AF loci were identified by a recent exome chip analysis that included 22 806 AF cases and 132 612 referents.17 PWI loci were identified from a meta-analysis of P-wave duration and P-wave terminal force that included 44 456 participants.18 In addition, for each of the top variants associated with PR, we also examined its association with AF and PWI.

    Examine Potential Function of PR-Related Variants for Gene Expression, Regulation, and Biological Pathways

    Pathway analysis was performed by MAGENTA19 with default settings. The summary result for the common variants was used as the input, and significant pathways were defined as those with a false discovery rate20 <0.05. The implication of genetic variants on cardiac gene expression (expression quantitative trait loci analysis) was performed by querying the GTEx database.21 At each PR-related locus, we identified the top variant and its neighboring variants that were within 500 kb and in linkage disequilibrium with the top variant (r2≥0.5). Four heart and vascular tissues were queried, including artery aorta, artery coronary, atrial appendage, and heart left ventricle. Significant expression quantitative trait locis were defined as those with false discovery rate <0.05. Regulatory regions were downloaded from the ENCODE Project22 and the NIH Roadmap Epigenomics Program.23 Four tracks were created: (1) included all 98 cell types from Roadmap epigenomics H3K27ac sites; (2) included only 4 heart tissues (aorta, right atrium, left ventricle, and right ventricle) from Roadmap epigenomics H3K27ac sites; (3) included all 125 cell lines from ENCODE Dnase HS sites; and (4) included only 3 heart-derived cell lines (cardiac fibroblasts, atrial fibroblasts, and cardiac myocytes). The enrichment of PR-related loci in regulatory regions was examined by the VSE R package.24 For comparison, we randomly created 1000 variant sets with MAF values and linkage disequilibrium structures similar to those seen for PR-related loci.

    Results

    The current analyses included a total of 92 803 individuals from 27 cohorts, with 83 367 individuals from 22 studies of European ancestry and 9436 individuals from 5 studies of African ancestry. Clinical characteristics of the study participants are in Table 1.

    Table 1. Clinical Characteristics of the Participating Studies

    AncestryStudyTotal NMen, N (%)Age, y, MeanPR Interval, ms, Mean±SDRR Interval, ms, Mean±SDBMI, kg/m2, Mean±SDHeight, cm, Mean±SDSBP, mm Hg, mean±SDβ-Blockers (%)Diuretics (%)Calcium Antagonists* (%)
    European ancestryAGES2052742 (36.2)75.9±5.4170.5±26.8895±12927.0±4.4166±9143±20635 (31.0)ND108 (5.3)
    ARIC98284528 (46.1)54.1±5.7160.3±23.3928±13626.9±4.7169±9118±17789 (8.0)1085 (11.0)176 (1.8)
    BRIGHT841324 (38.9)57.6±10.7161.1±19.9960±16927.5±3.8166±9153±24248 (29.5)260 (30.9)18 (2.1)
    CAMP24931394 (55.9)60.7±11.6163.0±26.8926±16628.5±5.8171±10NDExcludedNDExcluded
    CHS32471313 (40.4)72.4±5.4167.8±28.2956±15126.4±4.4165±9136±21366 (11.3)750 (23.1)206 (6.3)
    ERF982447 (45.5)48.2±14.3152.9±23.2982±15926.9±4.6168±10140±203 (0.3)133 (13.5)25 (2.5)
    FHS75803428 (45.2)39.3±9.8152.0±22.1910±17526.0±5.0169±10119±15ExcludedNDExcluded
    GOCHA355161 (45.4)73.2±8.2167.6±27.7913±17326.1±4.6169±10N/AExcludedNDND
    GRAPHIC1755893 (50.9)39.1±14.5153.0±24.0934±14526.1±4.6171±9128±1939 (2.2)NDND
    INTER9958362843 (48.7)46.1±7.9158.2±22.4921±15026.3±4.6172±9130±18NDNDND
    KORA26171247 (47.6)48.3±13.0162.1±22.2944±14926.9±4.4168±9127±19199 (7.6)152 (5.8)13 (0.5)
    KORCULA293106 (36.2)55.0±13.4159.8±24.0929±12728.0±4.3168±9139±148 (2.7)3 (1.0)6 (2.0)
    LifeLines1934781 (40.3)45.2±13.1156.7±24.7896±14525.9±4.5175±9122±1664 (3.3)39 (2.0)23 (1.2)
    MESA24551171 (47.7)62.8±10.2164.7±25.21047±15827.8±5.1169±10123±21ExcludedNDExcluded
    NEO57822717 (47.0)55.9±5.9164.5±23.4940±15130.0±4.8174±10133±17ExcludedNDND
    RS23581086 (46.1)68.6±8.1168.2±24.7871±14426.3±3.6168±9ND293 (12.4)NDND
    GS:SFHS91683786 (41.3)52.0±13.6164.1±24.9886±14626.9±5.1168±10134±18192 (2.1)NDND
    SHIP64932608 (40.2)49.2±15.3158.5±23.3897±14627.5±5.0170±9131±20NDNDND
    TwinsUK46532 (6.9)52.3±11.7159.6±22.6923±14826.8±5.4163±7119±16NDNDND
    UHP1735779 (44.9)39.1±13.0155.9±22.5950±15124.9±3.9175±10125±1769 (4.6)32 (1.8)18 (1.0)
    WHI13 2520 (0)66.0±6.5161.4±24.0921±13828.7±5.6162±6130±18735 (5.5)1715 (13.3)1230 (9.3)
    YFS1846824 (44.6)41.9±5.0156.2±22.61028±16526.4±4.9172±9119±1438 (2.1)24 (1.3)1 (0.1)
    African ancestryARIC33661291 (38.4)53.3±5.8171.2±26.8929±15129.4±6.1168±9128±22315 (9.4)717 (21.3)222 (6.6)
    CHS627232 (37.0)72.4±5.5170.2±28.1918±16128.4±5.5165±9142±2254 (8.6)217 (34.6)60 (9.6)
    JHS2220833 (37.5)52.7±12.5172.7±27.3956±15031.4±6.4169±9126±18ExcludedNDExcluded
    MESA1565718 (45.9)62.3±10.0170.9±26.31050±17230.2±5.9168±10132±21ExcludedNDExcluded
    WHI16580 (0)64.6±6.4167.1±24.8921±14831.1±5.8162±7134±1787 (5.2)393 (23.7)341 (20.6)

    Exclusion criteria are given in Table I in the Data Supplement. AGES indicates Age, Gene/Environment Susceptibility Study; ARIC study, Atherosclerosis Risk in Communities; BMI, body mass index; BRIGHT, British Genetics of Hypertension; CAMP, Cardiology and Metabolic Patient; CHS, Cardiovascular Health Study; ERF, Erasmus Rucphen Family Study; FHS, Framingham Heart Study; GOCHA, Genes for Cerebral Hemorrhage on Anticoagulation; GRAPHIC, Genetic Regulation of Arterial Pressure in Humans in the Community; GS:SFHS, Generation Scotland: Scottish Family Health Study; JHS, Jackson Heart Study; KORA, Cooperative Health Research in the Region Augsburg; KORCULA, CROATIA-Korcula; LifeLines, LifeLines Cohort Study; MESA, Multi-Ethnic Study of Atherosclerosis; N/A, not available; ND, not determined; NEO, The Netherlands Epidemiology of Obesity; RS, Rotterdam Study; SBP, systolic blood pressure; SHIP, Study of Health in Pomerania; UHP, Utrecht Health Project; WHI, Women’s Health Initiative; and YFS, Young Finns Study.

    *Nondihydropyridine calcium antagonists.

    Identification of 31 Loci Associated With PR Interval

    A total of 42 075 common variants were analyzed (MAF ≥1%). As shown in Figure 1 and Table 2, 31 loci were significantly associated with PR interval after Bonferroni correction (P<1.2×10−6), including 22 loci that reached the conventional genome-wide significance threshold (P<5×10−8). The results of the random effects meta-analysis were similar to those of the fixed-effects analysis (Table II in the Data Supplement). The most significant locus was tagged by rs6795970 (P=4.0×10−240), a missense variant in SCN10A, which encodes a sodium channel that has been associated previously with the PR interval (r2=0.97 with the top SNP rs6599250 reported previously).8 Highly associated variants clustered in the linker region between the second and third domains of SCN10A (Figure 2). The top variants at 12 loci are missense variants. In addition, the top variants at 4 loci (including 3 novel loci) are low-frequency variants (1%<MAF<5%), illustrating the power of exome chip analyses to identify low-frequency coding associations. Detailed information of the nearest gene to each genome-wide significant locus is given in Table III in the Data Supplement.

    Table 2. Common Variants Significantly Associated With PR Interval From Meta-Analysis of All Studies

    SNPLocusClosest GeneFunctionCoding AlleleCAF*βSEP ValueNo. of Studies†Prolong or Shorten PR IntervalNovel Locus
    rs67959703p22.2SCN10AMissenseA0.370.17050.00524.0×10−24027Prolong
    rs39228443p22.2SCN5AIntronicA0.34−0.10690.00539.3×10−9026Shorten
    rs38079897q31.2CAV1IntronicA0.430.09080.00503.0×10−7427Prolong
    rs76607024q21.23ARHGAP24IntronicC0.33−0.09210.00531.2×10−6827Shorten
    rs1728729312p12.1LINC00477IntergenicG0.14−0.10840.00711.9×10−5227Shorten
    rs118971192p14MEIS1IntronicC0.390.05660.00554.2×10−2525Prolong
    rs189631212q24.21TBX3IntergenicG0.280.05640.00558.7×10−2526Prolong
    rs88307912q24.21TBX53′UTRG0.290.05500.00544.5×10−2426Prolong
    rs1162023563p22.2DLEC1MissenseA0.02−0.19530.01991.0×10−2227Shorten
    rs2512535q35.1CREBRFIntergenicG0.42−0.04390.00514.7×10−1826Shorten
    rs111537306q22.31SLC35F1IntergenicC0.47−0.04200.00499.5×10−1827ShortenNovel
    rs356586965q21.1PAMMissenseG0.040.09560.01198.5×10−1627Prolong
    rs20704923p22.2SLC22A14MissenseT0.100.06240.00834.0×10−1427ProlongNovel
    rs258589713q12.11XPO4IntronicA0.170.04710.00642.8×10−1327Prolong
    rs20429952q31.2TTNMissenseC0.260.03750.00574.3×10−1127Prolong
    rs43996932p25.1ID2IntergenicA0.340.03740.00589.1×10−1125Prolong
    rs4130668813q34ADPRHL1MissenseC0.030.10020.01737.4×10−922ProlongNovel
    rs47451q22EFNA1MissenseT0.490.02990.00531.2×10−826Prolong
    rs1107807817p12LINC00670IntronicA0.400.02810.00502.2×10−827Prolong
    rs6063261010q22.2SYNPO2LMissenseT0.15−0.03710.00684.5×10−827ShortenNovel
    rs1184878514q24.2SIPA1L1IntronicG0.240.03170.00584.6×10−827Prolong
    rs37334144q35.2FAT1MissenseA0.380.02800.00514.8×10−827Prolong
    rs173625882q31.2CCDC141MissenseA0.08−0.04910.00905.5×10−827ShortenNovel
    rs22961721p34.3MACF1MissenseG0.200.03260.00611.1×10−727ProlongNovel
    rs93986526q22.31GJA1IntergenicA0.140.03900.00741.3×10−726ProlongNovel
    rs4421774q22.1AFF1IntronicC0.42−0.02620.00501.8×10−726ShortenNovel
    rs70020028q24.3PLECMissenseA0.38−0.02720.00522.1×10−725ShortenNovel
    rs17682083p22.1MOBPIntronT0.250.02880.00573.6×10−727ProlongNovel
    rs21197884q34.1HAND2IntergenicC0.52−0.02460.00495.6×10−727ShortenNovel
    rs17391905‡1p32.3C1orf185IntergenicG0.03−0.06940.01429.6×10−727Shorten
    rs52429510q24.1ALDH18A1IntergenicA0.40−0.02610.00539.7×10−726Shorten

    CAF indicates coding allele frequency; and UTR, Untranslated region.

    †Some variants did not reach pass the quality filtering in respective studies and thus were excluded.

    ‡SNP was not significant if African participants were excluded.

    Figure 1.

    Figure 1. Manhattan plot showing the association between common variants and PR interval from combined ancestry analysis. The x axis represents the chromosomal position for each SNP, and the y axis represents the −log10(P value) of the association with PR interval. The dashed line represents the genome-wide significance cutoff of 5×10−8, and the blue line represents the Bonferroni P value cutoff of 1.3×10−6. Black color represents known loci, whereas red color represents novel loci.

    Figure 2.

    Figure 2. Diagram of sodium voltage-gated channel α subunit 10 (SCN10A). Each yellow circle represents a genetic variant with a P value less than the significance cutoff (1.2×10−6). Each red circle represents a genetic variant with a P value greater than the significance cutoff, but <0.05.

    We then examined the associations between these top PR variants with AF and electrocardiographic PWI. Eight out of 31 PR loci identified in our analysis were associated with AF after Bonferroni correction (P<0.05/31=1.6×10−3), consistent with some shared mechanisms between the regulation of PR interval and AF. Variants in SCN10A most significantly associated with PR interval were also significantly associated with AF (Table IV in the Data Supplement). Among PR-related SNPs, rs60632610 at the SYNPO2L locus was most significantly associated with AF (odds ratio, 1.90 [0.87–0.93]; P=1.5×10−10). Figure I in the Data Supplement shows the overlap among loci associated with PR interval, AF, and PWI.

    We also performed a sensitivity analysis that separated samples of European and African ancestry. As shown in Table V in the Data Supplement and Figure II in the Data Supplement, all of the 31 loci except rs17391905 at the 1p32.3 locus (P=2.6×10−6) were also significant in the analysis of European-only samples. Table VI in the Data Supplement and Figure III in the Data Supplement show the result for the analysis of African ancestry-only samples. Three loci were significant: SCN5A (rs3922844), SCN10A (rs6795970), and TBX5 (rs883079) after Bonferroni correction; P<1.3×10−6. All 3 loci were also significant in the analysis of European-only samples. The result from each individual study is shown in Table VII in the Data Supplement.

    Rare Variations in MYH6 and SCN5A Are Associated With PR Interval

    We next examined the association between PR interval and rare variants (MAF<1%) in gene regions. Variation in 2 gene regions, MYH6 (P=5.9×10−11) and SCN5A (P=1.1×10−7), was associated with PR interval (Table 3). Tables VIII and IX in the Data Supplement show the association of each rare variant within MYH6 and SCN5A with PR interval, respectively. MYH6 encodes a cardiac myosin heavy chain subunit, and SCN5A encodes the major cardiac sodium channel and was previously found to be associated with PR interval.8MYH6 was also recently found to associate with PWI.18 We also performed an ancestry-stratified analysis in the same way as the combined analysis. The same 2 gene regions were significant using data from European samples alone (P=4.1×10−12 and 8.3×10−7 for MYH6 and SCN5A, respectively). These 2 genes did not reach the significance cutoff in African samples (P=0.03 and 0.01 for MYH6 and SCN5A, respectively). Two other genes, HEATR2 (P=2.2×10−6) and THRAP3 (P=4.2×10−6), were significantly associated in African samples alone. However, in the combined analysis, these 2 genes were not significant (P=0.02 and 0.06 for HEATR2 and THRAP3, respectively), probably because of a low cumulative allele frequency.

    Table 3. Top 10 Gene Regions Associated With PR Interval by the SKAT Test*

    GeneP ValueQmeta†CMAFNo. of VariantsNo. of Studies With At Least 1 Rare VariantAverage Number of Variants in Each Study
    MYH6‡5.9×10−11‡23 537 3400.0215322712
    SCN5A1.1×10−7‡16 604 8430.0289352713
    GORASP11.3×10−514 361 2520.030816276
    NEBL1.9×10−511 787 6990.0309362711
    TRIML21.2×10−410 173 9780.0223232710
    SLC22A111.5×10−46 539 6560.013611276
    MTRF12.8×10−49 073 0980.023510263
    CD363.5×10−48 001 7770.015628279
    CAPRIN23.7×10−46 886 3750.016915277
    PIK3R66.0×10−49 763 3360.031623268

    CMAF indicates cumulative minor allele frequency; and SKAT, Sequence Kernel Association Test.

    *The analysis included only nonsynonymous and splice site rare variants (minor allele frequency<1%) within the gene regions.

    †Qmeta: the SKAT Q statistic, defined as where Wj is the weight, and Sj is the squared score.

    ‡The significance level for gene-based tests after Bonferroni correction was P<0.05/5759=8.7×10−6; the 2 genes reached this significant cutoff.

    In our secondary analysis of pooled samples, we analyzed only damaging variants, defined as nonsense mutations or alternations predicted to be damaging by PolyPhen-215 or SIFT.16 Three genes reached the significance cutoff (P<0.05/2030=2.5×10−5), including GORASP1 (P=1.1×10−5), NEBL (P=1.9×10−5), and SCN5A (P=2.2×10−5; Table X in the Data Supplement).

    Expression Quantitative Trait Loci Analysis

    We also performed expression quantitative trait loci analysis to determine whether any of the novel PR-related variants were associated with cardiac gene expression using data from GTEx.21 Eight loci were associated with expression of at least 1 gene in the atrial appendage, left ventricle, coronary artery, or aorta, suggesting the importance of these loci in the regulation of gene expression in heart or vascular tissues (Table XI in the Data Supplement).

    Enrichment of PR-Related Variants in Regulatory Regions

    We examined involvement of PR-related variants in regulatory function. As shown in Figure IV in the Data Supplement, PR-related variants were significantly enriched in regulatory regions in both primary heart tissues (Padj=3.7×10−9) and heart-derived cell lines (Padj=0.002), but not in all tissues (Padj>0.05). The observed enrichment suggested involvement of these loci in tissue-specific regulatory functions. In addition, the variants also tended to locate within evolutionarily conserved regions (Padj=2.8×10−5 for primates and 6.4×10−5 for mammals).

    Enrichment of PR-Related Variants in Biological Pathways

    We examined the enrichment of PR-related variants in biological pathways by MAGENTA.19 Table XII in the Data Supplement shows the top pathways identified. The most significant pathway was heart morphogenesis (P=3.6×10−5; false discovery rate=0.049), suggesting that many PR-related genes might be involved in cardiac development. The pathway was the only significant pathway after correction for multiple testing (false discovery rate<0.05).

    Discussion

    We conducted large-scale analyses of the genetic determinants of atrioventricular conduction in 92 803 individuals by studying the electrocardiographic PR interval. In total, we observed 31 genetic loci that were associated with atrioventricular conduction, 11 of which are novel. In aggregate, the results implicate loci containing genes encoding ion channels in the heart, sarcomeric proteins, cardiac transcription factors, and other proteins with unknown cardiac function. Our findings provide new insights to the current understanding of atrioventricular conduction, which is critical for cardiac function.

    Interestingly, rare variants in SCN5A and MYH6 were associated with PR interval. A missense mutation (D1275N) in SCN5A has previously been reported in a large family with multiple members affected by dilated cardiomyopathy, conduction disorder, and arrhythmia.25 The mutation, together with several other mutations within the same gene, has also been associated with dilated cardiomyopathy,26 AF,27 and long-QT syndrome.28–31 Rare mutations within MYH6 were associated with sick sinus syndrome,28 congenital heart defects,32 and atrial septal defects.33

    Our observations support and extend prior analyses of cardiac conduction. Most previous genome-wide association studies involved the study of common genetic variation in smaller samples of up to 28 517 individuals.8,10,11 In keeping with those prior studies, we again observed that SCN10A is the most prominent gene involved in atrioventricular conduction. A recent genome-wide association study was based on 105K samples corroborates many of our current findings.34 However, our current study had greater power than earlier analyses for assessment of rare coding variation.

    Our study has 2 major implications. First, our results underscore the utility of assessing coding variation as an efficient way to identify functional molecular domains. In particular, our findings provide insights into the functional topology of SCN10A. The SCN10A sodium channel gene is widely expressed in the nervous system and heart,21 but it has only recently been implicated in cardiac conduction8,34–36 and arrhythmias such as AF35 and Brugada syndrome.37SCN10A encodes an α-subunit (with 6 transmembrane-spanning regions), which forms tetrameric, voltage-gated sodium channels responsible for the Nav 1.8 late sodium channel current.38,39 We found a collection of amino acid substitutions in the linker region between the second and third domains of SCN10A that were associated with PR duration (Figure 2). Variants in this linker region that were associated with the PR interval also were associated with AF, suggesting that function of this domain may have important clinical implications.

    Prior work on the homologous SCN5A cardiac sodium channel gene—which is also a cardiac conduction locus—indicates that this linker region is critical for sodium channel inactivation. Sodium influx is predominantly responsible for cardiomyocyte depolarization. Moreover, channel inactivation is essential for restoration of the hyperpolarized state needed for cyclic cardiomyocyte depolarization and contraction. Therefore, variations in this linker region might be involved in Nav 1.8 inactivation. Other data are necessary to identify relationships among variation in the linker region, the late sodium channel current, and channel inactivation in both healthy and diseased states.

    Together with previously discovered susceptibility genes, our findings implicate genes in different functional classes that regulate atrioventricular conduction such as ion channels and cardiac transcription factors. In many cases, anomalies in these genes have been found to cause human cardiac diseases, such as congenital heart defects, primary cardiac conduction abnormalities, and syndromes predisposing to sudden cardiac death (Table III in the Data Supplement). Interestingly, some of the genes are not expressed (in high abundance) in the right atrial appendage or the left ventricle, according to existing data sets—although most are active in the heart (Table XIII in the Data Supplement). Atrioventricular nodal conduction also can be influenced by external tone from the autonomic nervous system. Therefore, further work is necessary to determine the mechanisms by which identified genes that are not expressed in the heart influence the PR interval.

    We acknowledge several limitations of our study. Since the PR interval was measured across many cohorts, it is possible that there is some heterogeneity that would diminish our power to detect modest associations. We excluded individuals with extreme values of PR interval, which might have been gleaned from large variations in cardiac conduction. We also performed inverse normal transformation on the raw PR interval to reduce the heterogeneity, which on the other hand might reduce the interpretability. Although we performed single-variant and gene-based tests, we did not examine the association of haplotype patterns with PR interval, so it is unclear whether there are any haplotypes that might be associated with PR interval. Most of the genetic variants analyzed were in exons. Therefore, the effects of variants within regulatory regions were not investigated. We note that the variants identified may not be causally related to the studied phenotypes (PR interval, AF, and PWI) but may be in linkage disequilibrium with causal variants. We anticipate that future increases in sample size with additional replications and more comprehensive genotyping platforms, such as denser SNP arrays or genome sequencing, will help address these limitations.

    In conclusion, we studied genetic variants associated with PR-interval duration and identified 31 common loci—including 11 that were novel—and 2 rare variant regions. Our findings greatly expand our knowledge of the genes that underlie atrioventricular conduction in the heart.

    Sources of Funding

    This work was partly supported by grants from the National Institutes of Health (NIH) to Drs Ellinor and Benjamin (2RO1HL092577 and R01HL128914), Dr Ellinor (K24HL105780), and Arking and Sotoodehnia (R01HL116747). Dr Ellinor is also supported by an Established Investigator Award from the American Heart Association (13EIA14220013) and by the Fondation Leducq (14CVD01). Dr Lin was partly supported by Boston University Digital Health Initiative, and the National Center for Advancing Translational Sciences, National Institutes of Health, through BU-CTSI Grant Number 1UL1TR001430. N. Verweij is supported by NWO-VENI (016.186.125) and Marie Sklodowska-Curie GF (call: H2020-MSCA-IF-2014, Project ID: 661395). Dr Lubitz is supported by NIH grants K23HL114724 and a Doris Duke Charitable Foundation Clinical Scientist Development Award 2014105. Dr Sotoodehnia is supported by NIH grants HL116747 and HL111089. Dr Asselbergs is supported by a Dekker scholarship-Junior Staff Member 2014T001—Netherlands Heart Foundation and UCL Hospitals NIHR Biomedical Research Centre. Drs Hayward, Marten and Boutin were funded by MRC (UK) core funding to the MRCHGU QTL in Health and Disease programme. The views expressed in this article are those of the authors and do not necessarily represent the views of the National Heart, Lung, and Blood Institute; the National Institutes of Health; or the US Department of Health and Human Services.

    Disclosures

    Dr Ellinor is PI of a grant from Bayer AG to the Broad Institute focused on the genetics and therapeutics of atrial fibrillation. Dr Psaty serves on the Data and Safety Monitoring Board of a clinical trial funded by the manufacturer (Zoll LifeCor) and on the Steering Committee of the Yale Open Data Access Project funded by Johnson & Johnson. Dr Lubitz receives sponsored research support from Bayer HealthCare, Biotronik, and Boehringer Ingelheim and has consulted for St. Jude Medical and Quest Diagnostics. The other authors report no conflicts.

    Footnotes

    *Drs Ellinor, Lubitz, and Isaacs contributed equally to this work.

    The Data Supplement is available at http://circgenetics.ahajournals.org/lookup/suppl/doi:10.1161/CIRCGEN.117.002037/-/DC1.

    http://circgenetics.ahajournals.org

    Honghuang Lin, PhD, Section of Computational Biomedicine, Department of Medicine, Boston University School of Medicine, 72 E Concord St, B-616, Boston, MA 02118, E-mail or Steven A. Lubitz, MD, MPH, Cardiac Arrhythmia Service & Cardiovascular Research Center, Massachusetts General Hospital, 55 Fruit St, GRB 109, Boston, MA 02114, E-mail

    References

    • 1. Soliman EZ, Prineas RJ, Case LD, Zhang ZM, Goff DC. Ethnic distribution of ECG predictors of atrial fibrillation and its impact on understanding the ethnic distribution of ischemic stroke in the atherosclerosis risk in communities (ARIC) study.Stroke. 2009; 40:1204–1211. doi: 10.1161/STROKEAHA.108.534735.LinkGoogle Scholar
    • 2. Cheng S, Keyes MJ, Larson MG, McCabe EL, Newton-Cheh C, Levy D, et al. Long-term outcomes in individuals with prolonged PR interval or first-degree atrioventricular block.JAMA. 2009; 301:2571–2577. doi: 10.1001/jama.2009.888.CrossrefMedlineGoogle Scholar
    • 3. Xiao HB, Roy C, Fujimoto S, Gibson DG. Natural history of abnormal conduction and its relation to prognosis in patients with dilated cardiomyopathy.Int J Cardiol. 1996; 53:163–170.CrossrefMedlineGoogle Scholar
    • 4. Thiene G, Pennelli N, Rossi L. Cardiac conduction system abnormalities as a possible cause of sudden death in young athletes.Hum Pathol. 1983; 14:704–709.CrossrefMedlineGoogle Scholar
    • 5. Pilia G, Chen WM, Scuteri A, Orrú M, Albai G, Dei M, et al. Heritability of cardiovascular and personality traits in 6,148 Sardinians.PLoS Genet. 2006; 2:e132. doi: 10.1371/journal.pgen.0020132.CrossrefMedlineGoogle Scholar
    • 6. Silva CT, Kors JA, Amin N, Dehghan A, Witteman JC, Willemsen R, et al. Heritabilities, proportions of heritabilities explained by GWAS findings, and implications of cross-phenotype effects on PR interval.Hum Genet. 2015; 134:1211–1219. doi: 10.1007/s00439-015-1595-9.CrossrefMedlineGoogle Scholar
    • 7. Newton-Cheh C, Guo CY, Wang TJ, O’donnell CJ, Levy D, Larson MG. Genome-wide association study of electrocardiographic and heart rate variability traits: the Framingham Heart Study.BMC Med Genet. 2007; 8(suppl 1):S7. doi: 10.1186/1471-2350-8-S1-S7.CrossrefMedlineGoogle Scholar
    • 8. Pfeufer A, van Noord C, Marciante KD, Arking DE, Larson MG, Smith AV, et al. Genome-wide association study of PR interval.Nat Genet. 2010; 42:153–159. doi: 10.1038/ng.517.CrossrefMedlineGoogle Scholar
    • 9. Holm H, Gudbjartsson DF, Arnar DO, Thorleifsson G, Thorgeirsson G, Stefansdottir H, et al. Several common variants modulate heart rate, PR interval and QRS duration.Nat Genet. 2010; 42:117–122. doi: 10.1038/ng.511.CrossrefMedlineGoogle Scholar
    • 10. Butler AM, Yin X, Evans DS, Nalls MA, Smith EN, Tanaka T, et al. Novel loci associated with PR interval in a genome-wide association study of 10 African American cohorts.Circ Cardiovasc Genet. 2012; 5:639–646. doi: 10.1161/CIRCGENETICS.112.963991.LinkGoogle Scholar
    • 11. Smith JG, Magnani JW, Palmer C, Meng YA, Soliman EZ, Musani SK, et al; Candidate-gene Association Resource (CARe) Consortium. Genome-wide association studies of the PR interval in African Americans.PLoS Genet. 2011; 7:e1001304. doi: 10.1371/journal.pgen.1001304.CrossrefMedlineGoogle Scholar
    • 12. Grove ML, Yu B, Cochran BJ, Haritunians T, Bis JC, Taylor KD, et al. Best practices and joint calling of the HumanExome BeadChip: the CHARGE Consortium.PLoS One. 2013; 8:e68095. doi: 10.1371/journal.pone.0068095.CrossrefMedlineGoogle Scholar
    • 13. Lee S, Emond MJ, Bamshad MJ, Barnes KC, Rieder MJ, Nickerson DA, et al; NHLBI GO Exome Sequencing Project—ESP Lung Project Team. Optimal unified approach for rare-variant association testing with application to small-sample case-control whole-exome sequencing studies.Am J Hum Genet. 2012; 91:224–237. doi: 10.1016/j.ajhg.2012.06.007.CrossrefMedlineGoogle Scholar
    • 14. Wu MC, Lee S, Cai T, Li Y, Boehnke M, Lin X. Rare-variant association testing for sequencing data with the sequence kernel association test.Am J Hum Genet. 2011; 89:82–93. doi: 10.1016/j.ajhg.2011.05.029.CrossrefMedlineGoogle Scholar
    • 15. Adzhubei I, Jordan DM, Sunyaev SR. Predicting functional effect of human missense mutations using PolyPhen-2.Curr Protoc Hum Genet. 2013; Chapter 7:Unit7.20. doi: 10.1002/0471142905.hg0720s76.MedlineGoogle Scholar
    • 16. Kumar P, Henikoff S, Ng PC. Predicting the effects of coding non-synonymous variants on protein function using the SIFT algorithm.Nat Protoc. 2009; 4:1073–1081. doi: 10.1038/nprot.2009.86.