Whole genome sequencing based analysis of inflammation biomarkers in the Trans-Omics for Precision Medicine (TOPMed) consortium Free

Abstract

Introduction

Chronic inflammation is a risk factor for many diseases including cardiovascular disease, asthma, cancer and diabetes [1–3]. Chronic inflammation has been assessed in human cohorts using a variety of immunoassay measured biomarker traits, particularly markers of innate immune system activation such as C-reactive protein (CRP) and interleukin 6 (IL-6) [2]. Though there is a strong influence of social and environmental factors, previous analyses, including genome-wide association studies (GWAS), have demonstrated an underlying genetic component to variance in these traits [4, 5]. Heritability of biomarkers of inflammation have been estimated, for instance, to be 25%–60% [6, 7] for IL-6 and 30%–45% [8–12] for CRP. However, most studies have only analyzed relatively small and ancestrally homogenous (mostly European ancestry) populations and as such have not fully elucidated the genetic influence on these traits [4, 13–16].

The National Heart Lung and Blood Institute’s Trans-Omics for Precision Medicine (TOPMed) initiative has now generated whole genome sequencing data on > 150 000 individuals from diverse population-based cohorts enriched for heart, lung, and blood relevant disease traits. Novel ancestry-differentiated variant associations for CRP [17] (including confirmation of regulatory impacts in vitro) and E-selectin [18] reported in earlier TOPMed publications demonstrated the potential for genetic discovery for inflammation traits in these diverse cohorts. Thus, analysis of more biomarkers across a larger, more diverse set of samples with the addition of rare variant aggregate tests may identify additional associated individual variants and genomic regions. Here, we perform single variant and aggregate rare variant analyses across 21 inflammation-related biomarkers, some of which are in moderate to low correlation (Fig. S1), assessed in TOPMed cohort studies (Table 1), including performing detailed conditional analyses to identify distinct genetic association signals. Our results both inform our understanding of inflammation trait biology and of the expected findings for sequencing-based analyses of complex traits, particularly protein quantitative biomarkers.

Results

Our analyses of 21 inflammation biomarkers, generally measured by ELISA, included 12 cohorts from the TOPMed Program (Table S1, Table S2); phenotype availability and sample characteristics varied by trait (Table S3). For example, for CRP as an example trait, average age was 57.6, and the participants analyzed were 64.8% female, 1.5% Asian, 22.9% non-Hispanic Black, 22.7% Hispanic/Latino, 52.9% non-Hispanic White. By contrast, due to different cohorts contributing, analysis of interleukin-6 included 70.9% non-Hispanic White participants. In brief, we performed single variant analysis to identify trait-associated loci, followed by stepwise conditional analysis to identify the total number of statistically distinct signals. We also conditioned on previously associated variants to identify distinct signals not identified in prior papers. We performed genetic region and gene centric rare variant set-based analyses for each trait and likewise conditioned on previously identified signals.

Of the 21 traits tested, CRP, E-selectin, intercellular adhesion molecule 1 (ICAM-1), interleukin 18 (IL-18), IL-6, lipoprotein-associated phospholipase A2 (Lp-PLA2) activity and mass, monocyte chemoattractant protein-1 (MCP-1), matrix metalloproteinase-9 (MMP-9), P-selectin, and tumor necrosis factor α receptor 2 (TNFR2) had at least 1 genome-wide significant locus in single variant analyses. Across these 11 traits there were a total of 30 genome-wide significant loci (P < 1.0 × 10⁻⁹ [21]) (Table S4, Figs S2–S32), for which stepwise conditional analysis revealed a total of 67 distinct signals (Table S5). After conditioning on previously identified associations (Table S6), 22 conditionally distinct variants across 8 loci remained locus-wide significant for 6 traits (Table S7 and Table 2, Fig. 1, significance thresholds listed in Table S7), and 1 trait (MMP-9) had a locus not reported in the GWAS catalog (Table 2, Fig. 1). We focus on these 7 traits with findings distinct from those already reported in the GWAS catalog below.

Figure 1

Single variant findings conditionally distinct from GWAS catalog variants. We report P-values for association for marginal and conditional results, reference population effect allele frequencies (EAF) by continental ancestry group as defined by 1000 genomes project (1000G) [19]—African (AFR), admixed American (AMR), East Asian (EAS), European (EUR), South Asian (SAS), as well as all participants in 1000G (ALL)—and the overall effect allele frequency for all participants included in our TOPMed analyses, and annotation principal components (aPCs) from FAVOR [20]. NA means the variant is not reported in the reference panel. We note that this information is available for all variants in close linkage disequilibrium with these lead variants in Table S8.

Open in new tabDownload slide

In aggregate rare variant analyses, we detected 51 significant gene-centric sets associated with 6 traits (Table S10A) and 214 significant 2-kb sliding windows associated with 7 traits (Table S11A). We observed 19 significant rare variant aggregate test associations (some in overlapping or adjoining regions) after conditioning on known variants from the GWAS catalog and single-variant signals in the present analysis (Table 3), with traits P-selectin, ICAM-1, CRP, Lp-PLA2 activity and mass, all of which also had conditionally distinct single variant results (Tables S10C and S11C). If possible, we attempted to replicate distinct single variant findings using semiquantitative inflammation biomarker measures (i.e. measures which do not give an exact protein concentration in their results, for example in mg/dl) from the SomaScan or Olink platforms in independent samples (Table 4), as well as, for CRP, with rare variant testing in UK Biobank. While not all rare variant signals could be replicated due to lack of data availability, we note that all distinct rare variant aggregation signals were in known regions, increasing the plausibility of their association with inflammation traits.

C-reactive protein

We identified genetic variants associated with CRP consistent with and expanding upon our previous analysis of CRP in 23 279 TOPMed participants [17]. All 8 distinct single variant signals at the CRP locus previously known in TOPMed [17] (in partially overlapping samples) were also found here. We identified 1 additional distinct signal, rare variant rs370370301 (TOPMed Effect Allele Frequency (EAF): 0.2%, 1000G EUR EAF: 0.1%, 1000G SAS EAF: 0.1%, and not available in all other populations in 1000G), which was still significant after all conditional analyses. This non-coding variant did not reach genome-wide significance in the previous TOPMed analysis (P = 5.0 × 10⁻⁶) but was significant (P = 1.4 × 10⁻¹¹) in the present analysis. This increase in significance is likely due to the larger number of participants, which is now 38 465. This variant was not previously identified in an analysis of CRP in UK Biobank (UKB) [22], likely due to the rare frequency. Rare variant analysis yielded 1 significantly associated gene-centric set of 54 missense rare variants (P = 3.6 × 10⁻²²) on CRP locus driven in part by rs77832441 (P = 7.8 × 10⁻¹⁶ for analysis of individual variant in TOPMed) (Table S10B). We also tested a similar gene-centric missense rare variant set for association in UKB (P = 6.4 × 10⁻³⁴ based on 116 variants, details in Table S15). rs77832441 (MAC = 153, EAF = 0.2%) was previously identified in Schick et al. [23] We note that rs77832441 was pruned from the conditional analysis list based on linkage disequilibrium (LD) (see Materials and Methods) but a variant in close LD, rs553202904 (r² = 0.97), was included (Tables S6, S10C and S11C), and the significance of the gene centric test was attenuated but still significant (missense set, P = 1.3 × 10⁻⁸, Tables 3 and S10A) when this signal was adjusted for, suggesting additional subthreshold CRP missense variants in particular remain to be identified as individually significant in larger analyses.

In addition to signals at the CRP locus, we also identified multiple loci in the single variant association analyses not previously detected in prior TOPMed analysis, including 3 with multiple distinct signals (LEPR, SALL1, APOE) (Table S5). Each of these signals were attenuated below the genome-wide significance threshold after adjusting for known associations from the GWAS catalog and other prior publications [17, 22].

E-selectin

There are 9 distinct signals at the SELL/SELE, FUT6, and ABO loci associated with E-selectin, and 2 distinct signals remaining at the ABO locus after conditioning on previously identified signals, including single variant signals from previous TOPMed analysis. This pair of signals, rs8176719 and rs374594061, were the second and third distinct signals in our marginal analysis. Variant rs8176719 is a frameshift insertion exonic variant common across all populations that tags blood group O [28]. We do note that in our prior work from TOPMed [18], while this was not captured as an independent genome-wide signal, associations of differential E-selectin levels across blood groups (with O treated as reference) were also observed. This variant’s association with E-selectin further illustrates the extensive pleiotropy of the ABO locus, which has been previously associated with many diseases and traits. E-selectin associated distinct variant rs374594061 is rare across all populations (TOPMed EAF: 0.9%, and not available in 1000G) and, likely as a consequence, has no previously reported associations in the GWAS catalog and was also not tested in available replication cohorts.

Intercellular adhesion molecule 1

For ICAM-1, we identified 9 distinct single variant signals at the ICAM1 and ABO loci; 5 distinct signals at ICAM1 remained (Table S7) after conditioning on known associations (Table S6). The GWAS conditionally significant association at rs5491, the fourth distinct signal in unconditional results at the ICAM1 locus, is an exonic variant (TOPMed EAF: 4.5%, 1000G AFR EAF: 25.0%, 1000G AMR EAF: 1.7%, 1000G EAS EAF: 5.3%, 1000G EUR EAF: 0.7%, 1000G SAS EAF: 2.0%) that is low frequency in most populations but common among African ancestry populations. There are 4 other conditionally distinct noncoding variants—rs11575071, rs139053442, rs28382777, rs5030400—at the ICAM1 locus (Table 2); most have low or rare frequency across all populations. As displayed in Fig. S17, there is some long-range LD for variants identified in the ICAM1 locus, notably for rs5491 (displayed in turquoise) in Fig. S17.

We also identified multiple conditionally significant rare variant set-based associations with ICAM-1 including 2 gene-centric sets (Table S10A) and 6 2-kb sliding windows (Table S11A, individual variants included in tests included in Table S11B), and 2 of them overlap the ICAM1 locus. We identify a set of missense rare variants at ICAM1, whose most significant variant was the identified rs139053442 association but which remains significant after conditioning on rs139053442 and other single variant findings from TOPMed and other studies (Table S6).

Matrix Metalloproteinase-9

We identified the MMP9 encoding gene for association with MMP-9 levels in single variant analysis. This cis pQTL locus included 1 distinct signal at intronic variant rs3918249 that was common in all populations, and it has repressed regulatory function with high H3K27me3 score 48 according to FAVOR [20]. Our identified variant rs3918249 (Fig. 1, TOPMed EAF 35.5%) is highly linked (r² = 0.938) with coding variant rs17576 (Table S8).

P-selectin

For P-selectin, we identified 5 distinct single variant signals at the SELP locus (Table S5), and 3 of them remain significant (Table S7) after conditioning on known associations (Table S6), and 1 distinct single variant signal at the ABO locus that is significant conditional on known associations (Tables S5 and S7). At the SELP locus, 2 of 3 conditionally significant signals are intronic (rs3917677, rs3917825). rs3917825 is relatively conserved (top 9.1% genome-wide aPC-conservation score) [20]. Both of these variants have low frequency in AFR ancestry participants (1.7% for rs3917677 in 1000G, 2.8% for rs3917825 in 1000G) and are not observed in EUR ancestry participants (from reference panels). The remaining significant signal in the SELP locus is the synonymous variant rs6128, which is more common in AFR ancestry (53.3%) than in EUR ancestry (16.6%) participants from 1000G. Variant rs6128 is a platelet splice QTL that alters SELP exon 14 skipping and soluble versus transmembrane P-selectin protein production [33].

For aggregate tests of rare variants, lead signals were detected at 2 consecutive 2-kb sliding windows in the SELP locus located at chr1:169615464-169617463 and chr1:169616464-169618463 (Table S11A), which are driven in part by rs7529463. This coding variant is highly conserved (top 1.6% genome-wide aPC-conservation score), very rare (TOPMed AF 0.1%), and has high aPC protein function scores (top 0.2% genome wide) [20].

At the ABO locus, the distinct signal (rs635634, which tags blood group A) remained significant after conditioning on known variants (Table S6); however, the P-value is significantly attenuated (from P = 1.0 × 10⁻⁵⁵ to P = 2.0 × 10⁻¹⁵, Table S7) when adjusting for known GWAS catalog variants.

Interleukin 6

We identified the IL6R locus in the marginal single variant analysis (lead variant rs61812598, P = 1.1 × 10⁻⁴⁹, Table S5). After conditioning on previous GWAS-identified variants (Table S6), the lead variant was rs568587329 (Table S7) (P = 1.2 × 10⁻⁶, when adjusted for known variants at the IL6R locus (GWAS catalog, [35]), which met our locus-wide significance threshold. This variant is rare in all populations (TOPMed EAF: 0.03%, 1000G AFR EAF: 0.4%, and not available in all other populations in 1000G) and has a high aPC-Transcription-Factor score 17.29 (top 1.87% genome wide) [20].