Skeletal muscle eQTL meta-analysis implicates genes in the genetic architecture of muscular and cardiometabolic traits

Abstract

Identifying genetic variants that regulate gene expression can help uncover mechanisms underlying complex traits. We performed a meta-analysis of skeletal muscle expression quantitative trait locus (eQTL) using data from 1,002 individuals from two studies. A stepwise analysis identified 18,818 conditionally distinct signals for 12,283 genes, and 35% of these genes contained two or more signals. Colocalization of these eQTL signals with 26 muscular and cardiometabolic trait genome-wide association studies (GWASs) identified 2,252 GWAS-eQTL colocalizations that nominated 1,342 candidate genes. Notably, 22% of the GWAS-eQTL colocalizations involved non-primary eQTL signals. Additionally, 37% of the colocalized GWAS-eQTL signals corresponded to the closest protein-coding gene, while 44% were located >50 kb from the transcription start site of the nominated gene. To assess tissue specificity for a heterogeneous trait, we compared colocalizations with type 2 diabetes (T2D) signals across muscle, adipose, liver, and islet eQTLs; we identified 551 candidate genes for 309 T2D signals representing 36% of T2D signals tested and over 100 more than were detected with any one tissue alone. We then functionally validated the allelic regulatory effect of an eQTL variant for INHBB linked to T2D in both muscle and adipose tissue. Together, these results further demonstrate the value of skeletal muscle eQTLs in elucidating mechanisms underlying complex traits.

Keywords: GWAS; INHBB; allelic heterogeneity; colocalization; complex trait; eQTL; signal identification; skeletal muscle; transcriptomics; type 2 diabetes.

Figure 1.. Skeletal muscle eQTL gene and signal characteristics

(A) Study design for eQTL discovery. eQTL signals were identified based on a threshold of p ≤ 1 × 10⁻⁵. (B) Number of signals identified per gene. (C–E) Characteristics of the primary, secondary, and tertiary+ signals for 1,339 eGenes that contained three or more eQTL signals. (C) absolute effect size, (D) minor allele frequency (MAF), and (E) distance to the gene’s transcription start site (TSS).

Figure 2.. Three skeletal muscle eQTL signals for SH3RF2 identified in the meta-analysis

(A–C) SH3RF2 was identified as an eGene with (A) three conditionally distinct signals at rs340057, rs2097969, and rs4913059 in the meta-analysis, (B) one signal at rs340057 in GTEx, and (C) one signal at rs72818497 (r² = 0.48 with rs2097969) in FUSION. (D–F) Each of the three signals identified in the meta-analysis was isolated by conditioning on the other two signals. All plots are colored by linkage disequilibrium (LD) with the lead variants of the three signals as identified in the meta-analysis. Red lines indicate the significance threshold (p ≤ 1 × 10⁻⁵).

Figure 3.. GWAS signals did not always colocalize with the nearest protein-coding gene

(A) For the 1,713 colocalizations with an eQTL for a protein-coding gene, bars show the distance from the GWAS lead variants to the TSS of the nearest protein-coding gene and the colocalized eGene. Plot truncated at 5500 kb. (B) Example T2D GWAS signal colocalized with a muscle eQTL for the third-closest protein-coding gene, PCK2.

Figure 4.. Colocalizations between T2D and eQTLs across four tissues

The 862 T2D signals we identified in the GGI-European analysis (Suzuki et al.) were tested for colocalization with eQTLs from four tissue (A) Summary of the total number of colocalized (posterior probability PP_H4 ≥ 0.7) GWAS-eQTL signal pairs and the numbers of eGenes and T2D signals with at least one colocalized eQTL signal. The “Total across tissues” row shows the number of unique eGenes/T2D signals with at least one colocalization in at least one tissue. (B) UpSet plot showing eGenes linked to T2D signals across tissues. For each gene that had an eQTL signal colocalized with at least one T2D signal in at least one tissue, bars show which tissue(s) identified a colocalization with that gene. Many of the eGenes that only show evidence of colocalization in one tissue may not have been detected as eQTLs in the other tissues due to limited eQTL discovery power, low gene-expression levels, or differences in data sources or quality.

Figure 5.. Allelic effects of rs11688682 in putative regulatory elements in myoblasts and adipocytes

(A) Colocalized T2D GWAS (top) and skeletal muscle eQTL for INHBB (bottom), both with lead variant rs11688682. (B) rs11688682 is located 244 kb downstream of INHBB in a region of accessible chromatin (ATAC-seq) in at least myoblasts and adipocytes. (C and D) Transcriptional reporter assay results of a 610 bp fragment surrounding rs11688682, which showed enhancer activity in (C) LHCN-M2 myoblasts and (D) SGBS-derived differentiated adipocytes. EV, empty vector. The 7–10 points per allele or EV represent independent transfections and show the average of triplicate luciferase values. Bars show standard deviations; p values correspond to two-sided t tests. (E) Cartoon summarizing the directions of effect of the rs11688682^G allele associations with T2D risk and INHBB expression, and the observed effect on transcriptional activity.