Genome-wide association studies (GWASs) have identified thousands of single nucleotide polymorphisms (SNPs) associated with a wide range of complex traits (MacArthur et al., 2017; Visscher et al., 2017). However, the path from GWAS to biology is not straightforward as most SNPs implicated by GWASs reside in non-coding regions of the genome (MacArthur et al., 2017) and do not directly inform on the functional mechanism through which variants exert their effect on phenotypes.

GWASs have been performed on gene expression (Võsa et al., 2021), DNA methylation (Min et al., 2021), protein (Sun et al., 2018), and metabolite (Shin et al., 2014; Lotta et al., 2021) levels, identifying genetic variants influencing molecular traits, commonly referred to as molecular quantitative trait loci (molQTLs). The large overlap between complex and molecular trait-associated variants suggests that integrating these data can help interpreting GWAS loci (Vandiedonck, 2018; Taylor et al., 2019; Ongen et al., 2017). Advances in the field of transcriptomics make gene expression the best studied molecular phenotype, thanks to the presence of large expression QTL (eQTL) studies (e.g., eQTLGen Consortium [Võsa et al., 2021] N>30,000). Availability of these datasets fostered the development of summary statistic-based statistical approaches aiming at identifying associations between transcripts and complex traits (Zhu et al., 2016; Porcu et al., 2019; Hormozdiari et al., 2014; Gusev et al., 2016), prioritizing genes from known GWAS loci for functional follow-up, and inferring the directionality of these relations (Porcu et al., 2019; Porcu et al., 2021b). However, the cascade of events that mediates the effect of genetic variants on complex traits involves more than one molecular trait. Although approaches used for gene expression can be extended to other molecular data, investigating whether these molecular traits reside along the same causal pathway remains under-explored and only recently have studies applied colocalization and Mendelian randomization (MR) to methylation, gene expression, and protein levels data (Giambartolomei et al., 2018; Wu et al., 2018; Gleason et al., 2020; Sadler et al., 2022) and to a lesser extent to metabolic QTLs (mQTL) (Yin et al., 2022).

Metabolites are often the final products of cellular regulatory processes and the most proximal omic layer to complex phenotypes. Their levels could thus represent the ultimate response of biological systems to genetic and environmental changes. For instance, the metabolic status of organisms reflects disease progression, as metabolic disturbances can often be observed several years prior to the symptomatic phase (Shah et al., 2010; Wang et al., 2011; Sabatine et al., 2005). Therefore, using metabolomics to identify early-stage biomarkers of complex phenotypes, such as prediabetes and COVID-19 susceptibility, has gained increased interest (Wang-Sattler et al., 2012; Julkunen et al., 2021). While two-sample MR approaches using metabolites as single exposure have revealed biomarkers for several diseases (Qian et al., 2021; Lord et al., 2021; Porcu et al., 2021a), these analyses focused on the prediction of disease risk rather than on deciphering the mechanisms of discovered associations.

Integrating transcriptomics with metabolomics data can provide insights into how metabolites are regulated, elucidating targetable functional mechanisms. Here, we develop a framework based on established MR methodology that hypothesizes a mediating role of the metabolome in the transcript-to-phenotype axis, with the primary exposure being defined as an upstream omic layer, namely gene expression. Specifically, our integrative MR analysis combines summary-level multi-omics (i.e., GWAS, eQTL, and mQTL) data to compute the indirect effect of gene expression on complex traits mediated by metabolites in three steps (Figure 1). First, we map the transcriptome to the metabolome by identifying causal associations between transcripts and metabolites. Next, we screen metabolites for downstream causal effects on 28 complex phenotypes, resulting in the identification of gene expression → metabolite → phenotype cascades (Figure 1A). In parallel, we prioritize trait-associated genes by testing the association of transcripts with phenotypes (Figure 1B). Third, for transcripts found to causally influence either a metabolite (A) or a complex phenotype (B), we test whether the identified target genes exert their effect on the phenotype through the metabolite using multivariable MR (MVMR; Figure 1C). Finally, we carried out extensive power analyses to determine under which conditions the mediation analysis (Figure 1C) outperforms the conventional exposure-outcome MR framework (Figure 1B).

Figure 1 with 1 supplement see all

Download asset Open asset

In this study, we combined MR approaches integrating eQTL, mQTL, and GWAS summary statistics to explore the role of the metabolome in mediating the effect of the transcriptome on complex phenotypes. Applied to 28 medically relevant traits, our approach revealed 216 causal transcript-metabolite-phenotype triplets. Our automated literature review indicates that while some detected associations were previously reported, a large fraction, especially among the triplets, appears to be novel. It should be noted that the number of previously reported associations is likely underestimated as our approach does not account for all synonyms of a given feature and requires the terms to appear in the title or abstract of the publication. This makes it more likely for hypothesis-driven studies, inherently biased toward well-studied genes and metabolites, to be identified. Conversely, high-throughput, hypothesis-free studies that report the given association in a supplemental table are likely to be missed. Furthermore, due to its automated nature, our search is context-blind, so that some of the identified studies might report negative results, associations only under specific conditions (e.g., different organisms, experimental settings), or usage of the search term with a different meaning. To attenuate the latter, we also performed manual review of the retained abstracts for transcript-to-phenotype searches. While flawed, this rough estimate of the amount of existing evidence supporting our findings can be interpreted in combination with other lines of evidence. For instance, among the 90 signals that were also identified through TWMR, 93% showed a directionally concordant effect between the transcript-to-phenotype, transcript-to-metabolite, and metabolite-to-phenotype estimates (i.e., sign of product of the transcript-to-metabolite and metabolite-to-phenotype effects agrees with the sign of the transcript-to-phenotype effect). In these situations, dissection of causal effects provides clues as to the molecular mechanism through which involved genes modify complex phenotypes. This information is particularly valuable to identify key molecular mediators of highly pleiotropic genetic regions, such as the TMEM258/FADS1/FADS2 locus (Figure 3). While transcript levels of these genes affected eleven metabolites, three complex lipids were highlighted as strong molecular mediators of the transcript-to-phenotype effects.

Strikingly, 58% of the 216 causal transcript-metabolite-phenotype triplets were missed by TWMR – an approach that only considers gene expression and GWAS data. We highlight two novel but biologically plausible mechanisms linking ANKH to calcium levels through modulation of citrate and SLC6A12 to serum creatinine levels through regulation of the renal osmolyte betaine. Simulation analyses showed that these signals were likely missed by TWMR due to lack of power, as mediation analysis is better suited to detect associations with a low direct to total effect proportion and stronger transcript-to-metabolite than metabolite-to-phenotype effect. Promisingly, our simulations showed that mediation analysis becomes increasingly powerful over a wider range of parameter settings as the sample size of the mediator QTL study increases, highlighting the importance of generating large and publicly available molQTL datasets that can help to unravel functional gene-to-phenotype mechanisms.

As illustrated through the selected examples, a large fraction of detected mediations involves genes encoding metabolic enzymes or transporters/channels, with an enrichment for ‘secondary active transmembrane transporter activity’, for example (GO:0015291; FDR=0.021; background: 7884 genes with ≥3 IVs assessed through TWMR; STRING database). Matching the finding that the most likely effector genes of mQTLs are enriched for pathway-relevant enzymes and transporters (Smith et al., 2022), these results are not surprising given that the proteins encoded by these genes directly interact with metabolites, making it more likely that the effect of changes in their expression is mediated by metabolites. While our method is well suited to detect such effects, interpretation of discovered mediations is limited by the lack of spatial resolution of the mQTL data. Access to metabolite concentrations in different cellular compartments (e.g., extracellular space, cytosol, mitochondrial matrix, etc.) would generate more fine-tuned mechanistic hypotheses that consider the directionality of metabolite fluxes.

The observation that 77% of the transcript’s effect on the phenotype is not mediated by metabolites suggests that either true direct effects are frequent or that other unassessed metabolites or molecular layers (e.g., proteins, post-translational modifications, etc.) play a crucial role in mediation. It is to note that in the presence of unmeasured mediators or measured mediators without genetic instruments, our mediation estimates are lower bounds of the total existing mediation. In addition, unmeasured mediators sharing genetic instruments with the measured ones can modify result interpretation as some of the observed mediators may simply be correlates of the true underlying mediators. While this is a limitation of all MR methods, metabolic networks may harbor particularly large number of genetically correlated metabolite species. Similarly, owing to linkage disequilibrium and regulatory variants affecting multiple genes, transcripts from adjacent genes might appear to be involved in the same signals, as exemplified with the TMEM258/FADS1/FADS2 locus (Figure 3). While literature supports the role of the FADS genes, one cannot exclude a role for TMEM258, nor disentangle the specific function of FADS1 and FADS2. Thanks to the flexibility of the proposed framework, we expect that in the future and upon availability of ever larger and more diverse datasets, our method could be applied to estimate the relative contribution of currently unassessed mediators in translating genotypic cascades.

Another consideration is that complex phenotypes can have a stronger impact on gene expression than the opposite (Porcu et al., 2021b). Due to the lack of genome-wide trans-eQTL association summary statistics, our method does not investigate reverse causality on metabolites and gene expression, nor the role of metabolites as regulators of gene expression. Metabolites might also integrate the effect of several transcripts (i.e., multiple transcripts causally impact the levels of the same metabolite) before affecting complex phenotypes (Supplementary file 1g) or multiple metabolites may jointly mediate the impact of a single transcript. Modelling the latter phenomenon, which is beyond the scope of our current work, requires the development of structural equation models accounting for such effects and will eventually lead to a more comprehensive modelling of causal relations in complex biological networks, nuancing the interpretation of the molecular mechanisms shaping complex traits.

In conclusion, we developed a modular MR framework that has increased power over classical MR approaches to detect causal transcript-to-phenotype relationships when these are mediated by alteration of metabolite levels and is likely to become increasingly powerful upon release of larger molQTL datasets.

All data used in this study are publicly available. GWAS summary statistics for outcome traits measured in the UK Biobank originate from the Neale Lab (http://www.nealelab.is/uk-biobank/). eQTL data originated from the eQTLGen Consortium (https://www.eqtlgen.org) and was published in Vosa et al., 2021 [3]. mQTL data originate from Shin et al. 2014 [6], and are available at the Metabolomics GWAS Server (http://metabolomics.helmholtz-muenchen.de/gwas/). The Human Metabolome Database (HMDB) was used to annotate metabolites and the v5.0 release from 2021-11-09 of the "All proteins" file was downloaded to extract transcript-metabolite interactions (https://hmdb.ca/downloads). PubMed was used for the automated literature review (https://pubmed.ncbi.nlm.nih.gov). The UCSC Genome Browser (https://genome.ucsc.edu/) was used to visualize the FADS locus, while the GWAS Catalog was used to assess the number of reported GWAS signals in the region (https://www.ebi.ac.uk/gwas/). The STRING database was used for the enrichment analysis (https://string-db.org/). Produced data is available as Supplementary File 1 and Source Data. Code used to perform analyses is freely available at https://github.com/eleporcu/Gene_Metab_Pheno; (copy archived at swh:1:rev:c6bff8d094e369ff0d399751fc85fcd5ea250134).

1. 1. Borden LA
   2. Smith KE
   3. Gustafson EL
   4. Branchek TA
   5. Weinshank RL

   (1995) Cloning and expression of a betaine/GABA transporter from human brain

   Journal of Neurochemistry 64:977–984.

   https://doi.org/10.1046/j.1471-4159.1995.64030977.x

   • PubMed
   • Google Scholar
2. 1. Bowden J
   2. Davey Smith G
   3. Burgess S

   (2015) Mendelian randomization with invalid instruments: effect estimation and bias detection through egger regression

   International Journal of Epidemiology 44:512–525.

   https://doi.org/10.1093/ije/dyv080

   • PubMed
   • Google Scholar
3. 1. Burgess S
   2. Butterworth A
   3. Thompson SG

   (2013) Mendelian randomization analysis with multiple genetic variants using summarized data

   Genetic Epidemiology 37:658–665.

   https://doi.org/10.1002/gepi.21758

   • PubMed
   • Google Scholar
4. 1. Burgess S
   2. Bowden J
   3. Fall T
   4. Ingelsson E
   5. Thompson SG

   (2017) Sensitivity analyses for robust causal inference from mendelian randomization analyses with multiple genetic variants

   Epidemiology 28:30–42.

   https://doi.org/10.1097/EDE.0000000000000559

   • PubMed
   • Google Scholar
5. 1. Bycroft C
   2. Freeman C
   3. Petkova D
   4. Band G
   5. Elliott LT
   6. Sharp K
   7. Motyer A
   8. Vukcevic D
   9. Delaneau O
   10. O’Connell J
   11. Cortes A
   12. Welsh S
   13. Young A
   14. Effingham M
   15. McVean G
   16. Leslie S
   17. Allen N
   18. Donnelly P
   19. Marchini J

   (2018) The UK Biobank resource with deep phenotyping and genomic data

   Nature 562:203–209.

   https://doi.org/10.1038/s41586-018-0579-z

   • PubMed
   • Google Scholar
6. 1. Chang CC
   2. Chow CC
   3. Tellier LC
   4. Vattikuti S
   5. Purcell SM
   6. Lee JJ

   (2015) Second-Generation PLINK: rising to the challenge of larger and richer datasets

   GigaScience 4:7.

   https://doi.org/10.1186/s13742-015-0047-8

   • PubMed
   • Google Scholar
7. 1. Giambartolomei C
   2. Zhenli Liu J
   3. Zhang W
   4. Hauberg M
   5. Shi H
   6. Boocock J
   7. Pickrell J
   8. Jaffe AE
   9. Pasaniuc B
   10. Roussos P
   11. CommonMind Consortium

   (2018) A bayesian framework for multiple trait colocalization from summary association statistics

   Bioinformatics 34:2538–2545.

   https://doi.org/10.1093/bioinformatics/bty147

   • PubMed
   • Google Scholar
8. 1. Gil RB
   2. Ortiz A
   3. Sanchez-Niño MD
   4. Markoska K
   5. Schepers E
   6. Vanholder R
   7. Glorieux G
   8. Schmitt-Kopplin P
   9. Heinzmann SS

   (2018) Increased urinary osmolyte excretion indicates chronic kidney disease severity and progression rate

   Nephrology, Dialysis, Transplantation 33:2156–2164.

   https://doi.org/10.1093/ndt/gfy020

   • PubMed
   • Google Scholar
9. 1. Gleason KJ
   2. Yang F
   3. Pierce BL
   4. He X
   5. Chen LS

   (2020) Primo: integration of multiple GWAS and omics QTL summary statistics for elucidation of molecular mechanisms of trait-associated snps and detection of pleiotropy in complex traits

   Genome Biology 21:236.

   https://doi.org/10.1186/s13059-020-02125-w

   • PubMed
   • Google Scholar
10. 1. Granchi D
    2. Baldini N
    3. Ulivieri FM
    4. Caudarella R

    (2019) Role of citrate in pathophysiology and medical management of bone diseases

    Nutrients 11:2576.

    https://doi.org/10.3390/nu11112576

    • PubMed
    • Google Scholar
11. 1. Gusev A
    2. Ko A
    3. Shi H
    4. Bhatia G
    5. Chung W
    6. Penninx BWJH
    7. Jansen R
    8. de Geus EJC
    9. Boomsma DI
    10. Wright FA
    11. Sullivan PF
    12. Nikkola E
    13. Alvarez M
    14. Civelek M
    15. Lusis AJ
    16. Lehtimäki T
    17. Raitoharju E
    18. Kähönen M
    19. Seppälä I
    20. Raitakari OT
    21. Kuusisto J
    22. Laakso M
    23. Price AL
    24. Pajukanta P
    25. Pasaniuc B

    (2016) Integrative approaches for large-scale transcriptome-wide association studies

    Nature Genetics 48:245–252.

    https://doi.org/10.1038/ng.3506

    • PubMed
    • Google Scholar
12. 1. Ho AM
    2. Johnson MD
    3. Kingsley DM

    (2000) Role of the mouse ANK gene in control of tissue calcification and arthritis

    Science 289:265–270.

    https://doi.org/10.1126/science.289.5477.265

    • PubMed
    • Google Scholar
13. 1. Hormozdiari F
    2. Kostem E
    3. Kang EY
    4. Pasaniuc B
    5. Eskin E

    (2014) Identifying causal variants at loci with multiple signals of association

    Genetics 198:497–508.

    https://doi.org/10.1534/genetics.114.167908

    • PubMed
    • Google Scholar
14. 1. Huang J
    2. Howie B
    3. McCarthy S
    4. Memari Y
    5. Walter K
    6. Min JL
    7. Danecek P
    8. Malerba G
    9. Trabetti E
    10. Zheng HF
    11. Gambaro G
    12. Richards JB
    13. Durbin R
    14. Timpson NJ
    15. Marchini J
    16. Soranzo N
    17. UK10K Consortium

    (2015) Improved imputation of low-frequency and rare variants using the UK10K haplotype reference panel

    Nature Communications 6:8111.

    https://doi.org/10.1038/ncomms9111

    • PubMed
    • Google Scholar
15. 1. Jouret F
    2. Leenders J
    3. Poma L
    4. Defraigne J-O
    5. Krzesinski J-M
    6. de Tullio P

    (2016) Nuclear magnetic resonance metabolomic profiling of mouse kidney, urine and serum following renal ischemia/reperfusion injury

    PLOS ONE 11:e0163021.

    https://doi.org/10.1371/journal.pone.0163021

    • PubMed
    • Google Scholar
16. 1. Julkunen H
    2. Cichońska A
    3. Slagboom PE
    4. Würtz P
    5. Nightingale Health UK Biobank Initiative

    (2021) Metabolic biomarker profiling for identification of susceptibility to severe pneumonia and COVID-19 in the general population

    eLife 10:e63033.

    https://doi.org/10.7554/eLife.63033

    • PubMed
    • Google Scholar
17. 1. Lord J
    2. Jermy B
    3. Green R
    4. Wong A
    5. Xu J
    6. Legido-Quigley C
    7. Dobson R
    8. Richards M
    9. Proitsi P

    (2021) Mendelian randomization identifies blood metabolites previously linked to midlife cognition as causal candidates in alzheimer’s disease

    PNAS 118:e2009808118.

    https://doi.org/10.1073/pnas.2009808118

    • PubMed
    • Google Scholar
18. 1. Lotta LA
    2. Pietzner M
    3. Stewart ID
    4. Wittemans LBL
    5. Li C
    6. Bonelli R
    7. Raffler J
    8. Biggs EK
    9. Oliver-Williams C
    10. Auyeung VPW
    11. Luan J
    12. Wheeler E
    13. Paige E
    14. Surendran P
    15. Michelotti GA
    16. Scott RA
    17. Burgess S
    18. Zuber V
    19. Sanderson E
    20. Koulman A
    21. Imamura F
    22. Forouhi NG
    23. Khaw K-T
    24. MacTel Consortium
    25. Griffin JL
    26. Wood AM
    27. Kastenmüller G
    28. Danesh J
    29. Butterworth AS
    30. Gribble FM
    31. Reimann F
    32. Bahlo M
    33. Fauman E
    34. Wareham NJ
    35. Langenberg C

    (2021) A cross-platform approach identifies genetic regulators of human metabolism and health

    Nature Genetics 53:54–64.

    https://doi.org/10.1038/s41588-020-00751-5

    • PubMed
    • Google Scholar
19. 1. MacArthur J
    2. Bowler E
    3. Cerezo M
    4. Gil L
    5. Hall P
    6. Hastings E
    7. Junkins H
    8. McMahon A
    9. Milano A
    10. Morales J
    11. Pendlington ZM
    12. Welter D
    13. Burdett T
    14. Hindorff L
    15. Flicek P
    16. Cunningham F
    17. Parkinson H

    (2017) The new NHGRI-EBI catalog of published genome-wide association studies (GWAS catalog)

    Nucleic Acids Research 45:D896–D901.

    https://doi.org/10.1093/nar/gkw1133

    • PubMed
    • Google Scholar
20. 1. Min JL
    2. Hemani G
    3. Hannon E
    4. Dekkers KF
    5. Castillo-Fernandez J
    6. Luijk R
    7. Carnero-Montoro E
    8. Lawson DJ
    9. Burrows K
    10. Suderman M
    11. Bretherick AD
    12. Richardson TG
    13. Klughammer J
    14. Iotchkova V
    15. Sharp G
    16. Al Khleifat A
    17. Shatunov A
    18. Iacoangeli A
    19. McArdle WL
    20. Ho KM
    21. Kumar A
    22. Söderhäll C
    23. Soriano-Tárraga C
    24. Giralt-Steinhauer E
    25. Kazmi N
    26. Mason D
    27. McRae AF
    28. Corcoran DL
    29. Sugden K
    30. Kasela S
    31. Cardona A
    32. Day FR
    33. Cugliari G
    34. Viberti C
    35. Guarrera S
    36. Lerro M
    37. Gupta R
    38. Bollepalli S
    39. Mandaviya P
    40. Zeng Y
    41. Clarke T-K
    42. Walker RM
    43. Schmoll V
    44. Czamara D
    45. Ruiz-Arenas C
    46. Rezwan FI
    47. Marioni RE
    48. Lin T
    49. Awaloff Y
    50. Germain M
    51. Aïssi D
    52. Zwamborn R
    53. van Eijk K
    54. Dekker A
    55. van Dongen J
    56. Hottenga J-J
    57. Willemsen G
    58. Xu C-J
    59. Barturen G
    60. Català-Moll F
    61. Kerick M
    62. Wang C
    63. Melton P
    64. Elliott HR
    65. Shin J
    66. Bernard M
    67. Yet I
    68. Smart M
    69. Gorrie-Stone T
    70. BIOS Consortium
    71. Shaw C
    72. Al Chalabi A
    73. Ring SM
    74. Pershagen G
    75. Melén E
    76. Jiménez-Conde J
    77. Roquer J
    78. Lawlor DA
    79. Wright J
    80. Martin NG
    81. Montgomery GW
    82. Moffitt TE
    83. Poulton R
    84. Esko T
    85. Milani L
    86. Metspalu A
    87. Perry JRB
    88. Ong KK
    89. Wareham NJ
    90. Matullo G
    91. Sacerdote C
    92. Panico S
    93. Caspi A
    94. Arseneault L
    95. Gagnon F
    96. Ollikainen M
    97. Kaprio J
    98. Felix JF
    99. Rivadeneira F
    100. Tiemeier H
    101. van IJzendoorn MH
    102. Uitterlinden AG
    103. Jaddoe VWV
    104. Haley C
    105. McIntosh AM
    106. Evans KL
    107. Murray A
    108. Räikkönen K
    109. Lahti J
    110. Nohr EA
    111. Sørensen TIA
    112. Hansen T
    113. Morgen CS
    114. Binder EB
    115. Lucae S
    116. Gonzalez JR
    117. Bustamante M
    118. Sunyer J
    119. Holloway JW
    120. Karmaus W
    121. Zhang H
    122. Deary IJ
    123. Wray NR
    124. Starr JM
    125. Beekman M
    126. van Heemst D
    127. Slagboom PE
    128. Morange P-E
    129. Trégouët D-A
    130. Veldink JH
    131. Davies GE
    132. de Geus EJC
    133. Boomsma DI
    134. Vonk JM
    135. Brunekreef B
    136. Koppelman GH
    137. Alarcón-Riquelme ME
    138. Huang R-C
    139. Pennell CE
    140. van Meurs J
    141. Ikram MA
    142. Hughes AD
    143. Tillin T
    144. Chaturvedi N
    145. Pausova Z
    146. Paus T
    147. Spector TD
    148. Kumari M
    149. Schalkwyk LC
    150. Visscher PM
    151. Davey Smith G
    152. Bock C
    153. Gaunt TR
    154. Bell JT
    155. Heijmans BT
    156. Mill J
    157. Relton CL

    (2021) Genomic and phenotypic insights from an atlas of genetic effects on DNA methylation

    Nature Genetics 53:1311–1321.

    https://doi.org/10.1038/s41588-021-00923-x

    • PubMed
    • Google Scholar
21. 1. Ongen H
    2. Brown AA
    3. Delaneau O
    4. Panousis NI
    5. Nica AC
    6. GTEx Consortium
    7. Dermitzakis ET

    (2017) Estimating the causal tissues for complex traits and diseases

    Nature Genetics 49:1676–1683.

    https://doi.org/10.1038/ng.3981

    • PubMed
    • Google Scholar
22. 1. Porcu E
    2. Rüeger S
    3. Lepik K
    4. Santoni FA
    5. Reymond A
    6. Kutalik Z
    7. eQTLGen Consortium
    8. BIOS Consortium

    (2019) Mendelian randomization integrating GWAS and eqtl data reveals genetic determinants of complex and clinical traits

    Nature Communications 10:3300.

    https://doi.org/10.1038/s41467-019-10936-0

    • PubMed
    • Google Scholar
23. 1. Porcu E
    2. Gilardi F
    3. Darrous L
    4. Yengo L
    5. Bararpour N
    6. Gasser M
    7. Marques-Vidal P
    8. Froguel P
    9. Waeber G
    10. Thomas A
    11. Kutalik Z

    (2021a) Triangulating evidence from longitudinal and mendelian randomization studies of metabolomic biomarkers for type 2 diabetes

    Scientific Reports 11:6197.

    https://doi.org/10.1038/s41598-021-85684-7

    • PubMed
    • Google Scholar
24. 1. Porcu E
    2. Sadler MC
    3. Lepik K
    4. Auwerx C
    5. Wood AR
    6. Weihs A
    7. Sleiman MSB
    8. Ribeiro DM
    9. Bandinelli S
    10. Tanaka T
    11. Nauck M
    12. Völker U
    13. Delaneau O
    14. Metspalu A
    15. Teumer A
    16. Frayling T
    17. Santoni FA
    18. Reymond A
    19. Kutalik Z

    (2021b) Differentially expressed genes reflect disease-induced rather than disease-causing changes in the transcriptome

    Nature Communications 12:5647.

    https://doi.org/10.1038/s41467-021-25805-y

    • PubMed
    • Google Scholar
25. Software

    1. Porcu E

    (2022) Gene_Metab_Pheno, version swh:1:rev:c6bff8d094e369ff0d399751fc85fcd5ea250134

    Software Heritage.

    https://archive.softwareheritage.org/swh:1:dir:d314533b7c8b4abe22b4fdf9be7499279fdc1e6e;origin=https://github.com/eleporcu/Gene_Metab_Pheno;visit=swh:1:snp:ab86f3b131ad58088c7db0744b97590acbde2121;anchor=swh:1:rev:c6bff8d094e369ff0d399751fc85fcd5ea250134
26. 1. Qian L
    2. Fan Y
    3. Gao F
    4. Zhao B
    5. Yan B
    6. Wang W
    7. Yang J
    8. Ma X

    (2021) Genetically determined levels of serum metabolites and risk of neuroticism: a Mendelian randomization study

    The International Journal of Neuropsychopharmacology 24:32–39.

    https://doi.org/10.1093/ijnp/pyaa062

    • PubMed
    • Google Scholar
27. 1. Reynolds LM
    2. Dutta R
    3. Seeds MC
    4. Lake KN
    5. Hallmark B
    6. Mathias RA
    7. Howard TD
    8. Chilton FH

    (2020) Fads genetic and metabolomic analyses identify the ∆5 desaturase (FADS1) step as a critical control point in the formation of biologically important lipids

    Scientific Reports 10:15873.

    https://doi.org/10.1038/s41598-020-71948-1

    • PubMed
    • Google Scholar
28. 1. Sabatine MS
    2. Liu E
    3. Morrow DA
    4. Heller E
    5. McCarroll R
    6. Wiegand R
    7. Berriz GF
    8. Roth FP
    9. Gerszten RE

    (2005) Metabolomic identification of novel biomarkers of myocardial ischemia

    Circulation 112:3868–3875.

    https://doi.org/10.1161/CIRCULATIONAHA.105.569137

    • PubMed
    • Google Scholar
29. 1. Sadler MC
    2. Auwerx C
    3. Porcu E
    4. Kutalik Z

    (2022) Quantifying the role of transcript levels in mediating DNA methylation effects on complex traits and diseases

    Nat Commun 13:7559.

    https://doi.org/10.1038/s41467-022-35196-3

    • Google Scholar
30. 1. Shah SH
    2. Bain JR
    3. Muehlbauer MJ
    4. Stevens RD
    5. Crosslin DR
    6. Haynes C
    7. Dungan J
    8. Newby LK
    9. Hauser ER
    10. Ginsburg GS
    11. Newgard CB
    12. Kraus WE

    (2010) Association of a peripheral blood metabolic profile with coronary artery disease and risk of subsequent cardiovascular events

    Circulation. Cardiovascular Genetics 3:207–214.

    https://doi.org/10.1161/CIRCGENETICS.109.852814

    • PubMed
    • Google Scholar
31. 1. Shin SY
    2. Fauman EB
    3. Petersen AK
    4. Krumsiek J
    5. Santos R
    6. Huang J
    7. Arnold M
    8. Erte I
    9. Forgetta V
    10. Yang TP
    11. Walter K
    12. Menni C
    13. Chen L
    14. Vasquez L
    15. Valdes AM
    16. Hyde CL
    17. Wang V
    18. Ziemek D
    19. Roberts P
    20. Xi L
    21. Grundberg E
    22. Waldenberger M
    23. Richards JB
    24. Mohney RP
    25. Milburn MV
    26. John SL
    27. Trimmer J
    28. Theis FJ
    29. Overington JP
    30. Suhre K
    31. Brosnan MJ
    32. Gieger C
    33. Kastenmüller G
    34. Spector TD
    35. Soranzo N
    36. Multiple Tissue Human Expression Resource Consortium

    (2014) An atlas of genetic influences on human blood metabolites

    Nature Genetics 46:543–550.

    https://doi.org/10.1038/ng.2982

    • PubMed
    • Google Scholar
32. 1. Smith CJ
    2. Sinnott-Armstrong N
    3. Cichońska A
    4. Julkunen H
    5. Fauman EB
    6. Würtz P
    7. Pritchard JK

    (2022) Integrative analysis of metabolite GWAS illuminates the molecular basis of pleiotropy and genetic correlation

    eLife 11:e79348.

    https://doi.org/10.7554/eLife.79348

    • PubMed
    • Google Scholar
33. 1. Sun BB
    2. Maranville JC
    3. Peters JE
    4. Stacey D
    5. Staley JR
    6. Blackshaw J
    7. Burgess S
    8. Jiang T
    9. Paige E
    10. Surendran P
    11. Oliver-Williams C
    12. Kamat MA
    13. Prins BP
    14. Wilcox SK
    15. Zimmerman ES
    16. Chi A
    17. Bansal N
    18. Spain SL
    19. Wood AM
    20. Morrell NW
    21. Bradley JR
    22. Janjic N
    23. Roberts DJ
    24. Ouwehand WH
    25. Todd JA
    26. Soranzo N
    27. Suhre K
    28. Paul DS
    29. Fox CS
    30. Plenge RM
    31. Danesh J
    32. Runz H
    33. Butterworth AS

    (2018) Genomic atlas of the human plasma proteome

    Nature 558:73–79.

    https://doi.org/10.1038/s41586-018-0175-2

    • PubMed
    • Google Scholar
34. 1. Szeri F
    2. Lundkvist S
    3. Donnelly S
    4. Engelke UFH
    5. Rhee K
    6. Williams CJ
    7. Sundberg JP
    8. Wevers RA
    9. Tomlinson RE
    10. Jansen RS
    11. van de Wetering K

    (2020) The membrane protein ANKH is crucial for bone mechanical performance by mediating cellular export of citrate and ATP

    PLOS Genetics 16:e1008884.

    https://doi.org/10.1371/journal.pgen.1008884

    • PubMed
    • Google Scholar
35. 1. Taylor DL
    2. Jackson AU
    3. Narisu N
    4. Hemani G
    5. Erdos MR
    6. Chines PS
    7. Swift A
    8. Idol J
    9. Didion JP
    10. Welch RP
    11. Kinnunen L
    12. Saramies J
    13. Lakka TA
    14. Laakso M
    15. Tuomilehto J
    16. Parker SCJ
    17. Koistinen HA
    18. Davey Smith G
    19. Boehnke M
    20. Scott LJ
    21. Birney E
    22. Collins FS

    (2019) Integrative analysis of gene expression, DNA methylation, physiological traits, and genetic variation in human skeletal muscle

    PNAS 116:10883–10888.

    https://doi.org/10.1073/pnas.1814263116

    • PubMed
    • Google Scholar
36. 1. Vandiedonck C

    (2018) Genetic association of molecular traits: a help to identify causative variants in complex diseases

    Clinical Genetics 93:520–532.

    https://doi.org/10.1111/cge.13187

    • PubMed
    • Google Scholar
37. 1. Verbanck M
    2. Chen C-Y
    3. Neale B
    4. Do R

    (2018) Detection of widespread horizontal pleiotropy in causal relationships inferred from Mendelian randomization between complex traits and diseases

    Nature Genetics 50:693–698.

    https://doi.org/10.1038/s41588-018-0099-7

    • PubMed
    • Google Scholar
38. 1. Visscher PM
    2. Wray NR
    3. Zhang Q
    4. Sklar P
    5. McCarthy MI
    6. Brown MA
    7. Yang J

    (2017) 10 years of GWAS discovery: biology, function, and translation

    American Journal of Human Genetics 101:5–22.

    https://doi.org/10.1016/j.ajhg.2017.06.005

    • PubMed
    • Google Scholar
39. 1. Võsa U
    2. Claringbould A
    3. Westra H-J
    4. Bonder MJ
    5. Deelen P
    6. Zeng B
    7. Kirsten H
    8. Saha A
    9. Kreuzhuber R
    10. Yazar S
    11. Brugge H
    12. Oelen R
    13. de Vries DH
    14. van der Wijst MGP
    15. Kasela S
    16. Pervjakova N
    17. Alves I
    18. Favé M-J
    19. Agbessi M
    20. Christiansen MW
    21. Jansen R
    22. Seppälä I
    23. Tong L
    24. Teumer A
    25. Schramm K
    26. Hemani G
    27. Verlouw J
    28. Yaghootkar H
    29. Sönmez Flitman R
    30. Brown A
    31. Kukushkina V
    32. Kalnapenkis A
    33. Rüeger S
    34. Porcu E
    35. Kronberg J
    36. Kettunen J
    37. Lee B
    38. Zhang F
    39. Qi T
    40. Hernandez JA
    41. Arindrarto W
    42. Beutner F
    43. BIOS Consortium
    44. i2QTL Consortium
    45. Dmitrieva J
    46. Elansary M
    47. Fairfax BP
    48. Georges M
    49. Heijmans BT
    50. Hewitt AW
    51. Kähönen M
    52. Kim Y
    53. Knight JC
    54. Kovacs P
    55. Krohn K
    56. Li S
    57. Loeffler M
    58. Marigorta UM
    59. Mei H
    60. Momozawa Y
    61. Müller-Nurasyid M
    62. Nauck M
    63. Nivard MG
    64. Penninx BWJH
    65. Pritchard JK
    66. Raitakari OT
    67. Rotzschke O
    68. Slagboom EP
    69. Stehouwer CDA
    70. Stumvoll M
    71. Sullivan P
    72. ’t Hoen PAC
    73. Thiery J
    74. Tönjes A
    75. van Dongen J
    76. van Iterson M
    77. Veldink JH
    78. Völker U
    79. Warmerdam R
    80. Wijmenga C
    81. Swertz M
    82. Andiappan A
    83. Montgomery GW
    84. Ripatti S
    85. Perola M
    86. Kutalik Z
    87. Dermitzakis E
    88. Bergmann S
    89. Frayling T
    90. van Meurs J
    91. Prokisch H
    92. Ahsan H
    93. Pierce BL
    94. Lehtimäki T
    95. Boomsma DI
    96. Psaty BM
    97. Gharib SA
    98. Awadalla P
    99. Milani L
    100. Ouwehand WH
    101. Downes K
    102. Stegle O
    103. Battle A
    104. Visscher PM
    105. Yang J
    106. Scholz M
    107. Powell J
    108. Gibson G
    109. Esko T
    110. Franke L

    (2021) Large-Scale cis- and trans-eqtl analyses identify thousands of genetic loci and polygenic scores that regulate blood gene expression

    Nature Genetics 53:1300–1310.

    https://doi.org/10.1038/s41588-021-00913-z

    • PubMed
    • Google Scholar
40. 1. Wang TJ
    2. Larson MG
    3. Vasan RS
    4. Cheng S
    5. Rhee EP
    6. McCabe E
    7. Lewis GD
    8. Fox CS
    9. Jacques PF
    10. Fernandez C
    11. O’Donnell CJ
    12. Carr SA
    13. Mootha VK
    14. Florez JC
    15. Souza A
    16. Melander O
    17. Clish CB
    18. Gerszten RE

    (2011) Metabolite profiles and the risk of developing diabetes

    Nature Medicine 17:448–453.

    https://doi.org/10.1038/nm.2307

    • PubMed
    • Google Scholar
41. 1. Wang-Sattler R
    2. Yu Z
    3. Herder C
    4. Messias AC
    5. Floegel A
    6. He Y
    7. Heim K
    8. Campillos M
    9. Holzapfel C
    10. Thorand B
    11. Grallert H
    12. Xu T
    13. Bader E
    14. Huth C
    15. Mittelstrass K
    16. Döring A
    17. Meisinger C
    18. Gieger C
    19. Prehn C
    20. Roemisch-Margl W
    21. Carstensen M
    22. Xie L
    23. Yamanaka-Okumura H
    24. Xing G
    25. Ceglarek U
    26. Thiery J
    27. Giani G
    28. Lickert H
    29. Lin X
    30. Li Y
    31. Boeing H
    32. Joost H-G
    33. de Angelis MH
    34. Rathmann W
    35. Suhre K
    36. Prokisch H
    37. Peters A
    38. Meitinger T
    39. Roden M
    40. Wichmann H-E
    41. Pischon T
    42. Adamski J
    43. Illig T

    (2012) Novel biomarkers for pre-diabetes identified by metabolomics

    Molecular Systems Biology 8:615.

    https://doi.org/10.1038/msb.2012.43

    • PubMed
    • Google Scholar
42. 1. Wei Q
    2. Xiao X
    3. Fogle P
    4. Dong Z

    (2014) Changes in metabolic profiles during acute kidney injury and recovery following ischemia/reperfusion

    PLOS ONE 9:e106647.

    https://doi.org/10.1371/journal.pone.0106647

    • PubMed
    • Google Scholar
43. 1. Williams CJ

    (2016) The role of ANKH in pathologic mineralization of cartilage

    Current Opinion in Rheumatology 28:145–151.

    https://doi.org/10.1097/BOR.0000000000000247

    • PubMed
    • Google Scholar
44. 1. Wishart DS
    2. Guo A
    3. Oler E
    4. Wang F
    5. Anjum A
    6. Peters H
    7. Dizon R
    8. Sayeeda Z
    9. Tian S
    10. Lee BL
    11. Berjanskii M
    12. Mah R
    13. Yamamoto M
    14. Jovel J
    15. Torres-Calzada C
    16. Hiebert-Giesbrecht M
    17. Lui VW
    18. Varshavi D
    19. Varshavi D
    20. Allen D
    21. Arndt D
    22. Khetarpal N
    23. Sivakumaran A
    24. Harford K
    25. Sanford S
    26. Yee K
    27. Cao X
    28. Budinski Z
    29. Liigand J
    30. Zhang L
    31. Zheng J
    32. Mandal R
    33. Karu N
    34. Dambrova M
    35. Schiöth HB
    36. Greiner R
    37. Gautam V

    (2022) HMDB 5.0: the human metabolome database for 2022

    Nucleic Acids Research 50:D622–D631.

    https://doi.org/10.1093/nar/gkab1062

    • PubMed
    • Google Scholar
45. 1. Wu Y
    2. Zeng J
    3. Zhang F
    4. Zhu Z
    5. Qi T
    6. Zheng Z
    7. Lloyd-Jones LR
    8. Marioni RE
    9. Martin NG
    10. Montgomery GW
    11. Deary IJ
    12. Wray NR
    13. Visscher PM
    14. McRae AF
    15. Yang J

    (2018) Integrative analysis of omics summary data reveals putative mechanisms underlying complex traits

    Nature Communications 9:918.

    https://doi.org/10.1038/s41467-018-03371-0

    • PubMed
    • Google Scholar
46. 1. Yin X
    2. Bose D
    3. Kwon A
    4. Hanks SC
    5. Jackson AU
    6. Stringham HM
    7. Welch R
    8. Oravilahti A
    9. Fernandes Silva L
    10. FinnGen
    11. Locke AE
    12. Fuchsberger C
    13. Service SK
    14. Erdos MR
    15. Bonnycastle LL
    16. Kuusisto J
    17. Stitziel NO
    18. Hall IM
    19. Morrison J
    20. Ripatti S
    21. Palotie A
    22. Freimer NB
    23. Collins FS
    24. Mohlke KL
    25. Scott LJ
    26. Fauman EB
    27. Burant C
    28. Boehnke M
    29. Laakso M
    30. Wen X

    (2022) Integrating transcriptomics, metabolomics, and GWAS helps reveal molecular mechanisms for metabolite levels and disease risk

    American Journal of Human Genetics 109:1727–1741.

    https://doi.org/10.1016/j.ajhg.2022.08.007

    • PubMed
    • Google Scholar
47. 1. Zhu Z
    2. Zhang F
    3. Hu H
    4. Bakshi A
    5. Robinson MR
    6. Powell JE
    7. Montgomery GW
    8. Goddard ME
    9. Wray NR
    10. Visscher PM
    11. Yang J

    (2016) Integration of summary data from GWAS and eQTL studies predicts complex trait gene targets

    Nature Genetics 48:481–487.

    https://doi.org/10.1038/ng.3538

    • PubMed
    • Google Scholar

Decision letter after peer review:

Thank you for submitting your article "Exploiting the mediating role of the metabolome to unravel transcript-to-phenotype associations" for consideration by eLife. Your article has been reviewed by 2 peer reviewers, including Evan Graehl Williams as Reviewing Editor and Reviewer #1, and the evaluation has been overseen by David James as the Senior Editor.

The reviewers have discussed their reviews with one another, and the Reviewing Editor has drafted this to help you prepare a revised submission.

In brief, besides the reviewer-specific comments, there was a common thread regarding significances of the findings (both in the statistical sense and in the colloquial sense) and description of the method. I also apologize for the delay in the reviews; it was a huge challenge to get reviewers in July and August, and we potentially may consult a third reviewer for the revision.

Reviewer #1 (Recommendations for the authors):

Auwerx et al. have taken a new (or recent?) approach to mine large existing datasets of intermediary molecular data between GWAS and phenotype, with the aim of uncovering novel insight into the molecular mechanisms which lead a GWAS hit to have a phenotypic effect. The authors show that you can get additional insight by integrating multiple omics layers rather than analyzing only a single molecular type, including a handful of specific examples, e.g. that the effect of SNPs in ANKH on calcium are mediated by citrate. Such additional data is necessary because, as the authors' point out, while we have thousands of SNPs with significant impact on phenotypes of interest, we often don't know at all the mechanism, given that the majority of significant SNPs found through GWAS are in non-coding (and often intergenic) regions.

Scientific comments

(1) It's not quite clear to me whether the authors have designed a new method or have applied a recently-developed method. I follow that the general approach uses MR across multiple data types and then merges the results. I did look at the Github repo ("Gene_Metab_Pheno") but I'm not enough of a mathematical bioinformatician to really judge the novelty there, as to whether this is a new method or application of relatively recent other methods.

(2) The sample sizes are huge here since it is a meta-analysis, which means it's not quite clear to me how this method would scale to studies where there are not tens of thousands of individuals measured. It would be good to already have some general idea of what are the requirements in terms of power (sample size/target effect size / etc) for this approach – either added as a figure panel or two in this paper or just as text if the necessary parameters have already been explained elsewhere. I guess this comes back to the p >> n issue. I see for instance that the authors have 21k transcripts but they narrow it down to 7883 with {greater than or equal to} 3 IVs, used against their sample size of n = 7824 patients. Is this just convenient that n ~= p, or were the required parameters tweaked exactly so that you could reduce the number of transcripts down to more or less the sample size? You have a section on "Power Analysis" already which I appreciate, but reading through it again I can't tell what the minimum N would be. Or can I use an N = 500 I just have to use some other way to select only 500 (or whatever) transcripts? The power calculations go up to N = 90'000 which indeed is "state-of-the-art sample sizes" and would be unattainably large for all but the largest of international human cohorts, and perhaps a few large yeast or cell line studies.

(3) Related, I wasn't quite clear if the data from TwinsUK and KORA were analyzed separately, or if they were somehow merged.

(4) The authors mention "overall, 83% of the involved genes were causally influencing the level of a single metabolite, while TMEM258 and FADS2 affected 12 metabolites". Firstly I would mention that they affect 12 metabolites each just to be clear – although it's clearly annotated in Table S1. Second, and more relevantly, it would be nice to know in this table how many of these hits already are known in the literature. For instance I see the most significant hit is MAF1 with 5-oxoproline which doesn't immediately turn up anything (albeit I did not look very long). Hit #7, for SLC22A4 and isovalerylcarnitine does appear to be known in the literature ( http://www.metabolomix.com/rs272889-slc22a4-octn1-with-acylcarnitines/ ). Doing a decent literature analysis on all 257 of these hits in Table S1 might take a while (and/or require someone who knows how to do automated literature mining, which isn't trivial for this type of task especially given metabolite names like "X-12696") but even just some cherry-picking and/or (better) systematic analysis of the top 20-50 would be nice to see. It greatly improved my confidence in the findings when I just looked into two of your hits just now and found that one of them appears to be reasonably well known by literature.

(5) The authors mention that only 33% of the transcript-metabolite-phenotype triplets were observed by TWMR, but that 91% of those observed by all are in concordant directionality. I guess that means for those with concordant directionality, if no significant association was observed in one layer, there was just no cross-layer comparison? i.e. that the 91% "concordant directionality" are only checking among the 33% of hits by TWMR that were also observed both within and across layer?

(6) The authors mention no availability of trans-eQTL data; is that because the sample size is not sufficient to do a proper trans-eQTL analysis, or is it because of data privacy and an inability to get the precomputed trans-eQTLs from the eQTLGen Consortium due to some potential conflict (GDPR?).

(7) Back to Table S1: There are 257 hits here, with 191 genes (out of 7883) and 154 metabolites (out of 453), meaning that 299 metabolites didn't even have one significant association?

Reviewer #2 (Recommendations for the authors):

1. In general, the methodological approach should be explained in more detail. The study presents a novel MR-based approach and thus, the discussion of the method basics should be elaborated. This pertains to the following points:

– In the definition of causal effects from molecular traits on the outcome (eq. (1)), the symbol β' is not explained. The method used for standardization of effect sizes is not elaborated.

– For the multivariable MR analysis (equation not numbered), no causal effect statistic is specified. As discussed in the public review, no significance testing has been performed, but we strongly urge the authors to do so or give reasons for this lack of testing in the case of causal triplets.

– In the description of the simulation analyses, it should be made more clear how the definition of βGWAS as a function of βmQTL is reached.

– The authors should include a thorough discussion and a description of the method employed to evaluate the comparison between total and direct transcript-phenotype effects in Figure 2B.

2. There is some apparent inconsistency in the interpretation of results from the empirical study: The abstract states that 67 transcript-phenotype effects present within the 206 causal triplets were missed in a TWMR-only analyses, while the paragraph "Power Analysis" states that only 67 triplets showed a significant effect estimated by TWMR. The text relating to Figure 2B states that 79 % (~163) of the triplets were estimated to be due to direct effects by the transcript only. This should be presented in a more consistent and clear manner.

3. Plots for the distribution of the number of IVs per transcript/metabolite should be included (as supplemental figures).

4. It would be interesting to discuss potential modes of mediation between transcripts, metabolites and phenotypes. For example, individual metabolites might integrate the effects of several transcripts. This would not result in triplets, but ultimately, networks of mediation as shown for the TMEM258, FADS1, FADS2 case. How would such cases be treated in a statistically sound manner?

Reviewer #1 (Recommendations for the authors):

Auwerx et al. have taken a new (or recent?) approach to mine large existing datasets of intermediary molecular data between GWAS and phenotype, with the aim of uncovering novel insight into the molecular mechanisms which lead a GWAS hit to have a phenotypic effect. The authors show that you can get additional insight by integrating multiple omics layers rather than analyzing only a single molecular type, including a handful of specific examples, e.g. that the effect of SNPs in ANKH on calcium are mediated by citrate. Such additional data is necessary because, as the authors' point out, while we have thousands of SNPs with significant impact on phenotypes of interest, we often don't know at all the mechanism, given that the majority of significant SNPs found through GWAS are in non-coding (and often intergenic) regions.

Scientific comments

(1) It's not quite clear to me whether the authors have designed a new method or have applied a recently-developed method. I follow that the general approach uses MR across multiple data types and then merges the results. I did look at the Github repo ("Gene_Metab_Pheno") but I'm not enough of a mathematical bioinformatician to really judge the novelty there, as to whether this is a new method or application of relatively recent other methods.

While the core Mendelian randomization (MR) and mediation analysis methodology used in this study is well-established, the novelty of the manuscript lies in (i) the strategy that sequential application of univariable and multivariable MR analyses across the different omics layers (i.e., genome, transcriptome, metabolome, phenome) can improve discovery power and (ii) the resulting key message that by incorporating metabolomic data to genetic studies, it becomes possible to gain functional and targetable insights into the pathways linking genetic variation to phenotypic diversity. To clarify this point, we have modified the Introduction of our revised manuscript:

“Here, we develop a framework based on established MR methodology that hypothesizes a mediating role of the metabolome in the transcript-to-phenotype axis. Specifically, our integrative MR analysis combines summary-level multi-omics (i.e., GWAS, eQTL, and mQTL) data to compute the indirect effect of gene expression on complex traits mediated by metabolites in three steps (Figure 1).”

(2) The sample sizes are huge here since it is a meta-analysis, which means it's not quite clear to me how this method would scale to studies where there are not tens of thousands of individuals measured. It would be good to already have some general idea of what are the requirements in terms of power (sample size/target effect size / etc) for this approach – either added as a figure panel or two in this paper or just as text if the necessary parameters have already been explained elsewhere. I guess this comes back to the p >> n issue. I see for instance that the authors have 21k transcripts but they narrow it down to 7883 with {greater than or equal to} 3 IVs, used against their sample size of n = 7824 patients. Is this just convenient that n ~= p, or were the required parameters tweaked exactly so that you could reduce the number of transcripts down to more or less the sample size?

Requiring at least 3 IVs is common practice in the field of MR as it allows to estimate the causal estimate more robustly, as previously described (PMID: 31341166). In general, the more instruments, the less outliers can bias results. If the InSIDE assumption holds the MR estimates are asymptotically unbiased, hence more instruments reduce bias. Finally, more instruments also increase discovery power. Hence, the fact that the number of transcripts is reduced to the sample size of the study is pure coincidence. Of note, while this reduces the multiple testing burden, it also means that we potentially miss out on true causal effects of transcripts that cannot be instrumented properly.

You have a section on "Power Analysis" already which I appreciate, but reading through it again I can't tell what the minimum N would be. Or can I use an N = 500 I just have to use some other way to select only 500 (or whatever) transcripts? The power calculations go up to N = 90'000 which indeed is "state-of-the-art sample sizes" and would be unattainably large for all but the largest of international human cohorts, and perhaps a few large yeast or cell line studies.

While we are hopeful that in the future, ever larger multi-omics datasets will become available to the research community (e.g., initiatives similar to the NMR metabolite measurement in >120,000 UK Biobank participants; https://biobank.ctsu.ox.ac.uk/crystal/label.cgi?id=220), we agree that the majority of currently available molecular datasets have smaller sample sizes than the one used in this study. To address this comment, we ran power analyses simulating mQTLs studies performed on 1000, 2000, and 4000 individuals. Following the previously observed trend, the gain in power of mediation analysis over TWMR decreases with decreasing $N$ and is quasi inexistant with sample sizes smaller than 2,000. Hence, the reviewer’s point is well taken, the mediator sample size is crucial for our main conclusions. These additional results are detailed in the revised manuscript, along with Figure 4 —figure supplement 2:

“Repeating the simulations with a mQTL sample size of 90,000, nearing state-of-the-art sample sizes [7], we observe a strong shift in the above-described trends (Figure 4B). Specifically, when the effect of the transcript on the phenotype is dominated by the effect mediated by the metabolite (${\displaystyle \rho <0.3}$ and ${\displaystyle \rho >1.7}$), mediation analysis has more power than TWMR when ${\displaystyle \sigma >0.2}$. For larger proportions of direct effect, TWMR has increased power the more ${\displaystyle \sigma }$ differs from 1. In line with the increased power of mediation analysis with larger mQTL datasets, the gain in power of mediation analysis over TWMR decreases with decreasing mQTL dataset sample sizes (ranging between ${\displaystyle N=1,000}$ and ${\displaystyle N=4,000}$; Figure 4 —figure supplement 2), indicating that our approach is dependent on large sample sizes to reach its full potential.”

(3) Related, I wasn't quite clear if the data from TwinsUK and KORA were analyzed separately, or if they were somehow merged.

Data from KORA and TwinsUK have previously been analyzed, as described in Shin et al., 2014 (PMID: 24816252). Specifically, genome-wide association analyses were performed separately in KORA and TwinsUK for all 486 metabolite concentrations present in both studies after quality control steps. Association results were then combined using inverse variance meta-analysis based on effect size estimates and standard errors, adjusting for genomic control. Our study solely uses the publicly available meta-analysis summary statistics provided by Shin et al., 2014 for 453 of these metabolites.

We would also like to highlight that based on comments by Reviewer #1 and #2, the revised manuscript now focuses on a restricted set of 242 metabolites for which we were able to identify a Human Metabolome Database (HMDB) identifier through manual annotation. This approach is motivated by only analyzing metabolites whose potential associations can be followed up.

The current version of the manuscript reads as follows:

Results: “Summary statistics for cis-eQTLs stem from the eQTLGen Consortium meta-analysis of 19,942 transcripts in 31,684 individuals [3], while summary statistics for mQTLs originate from a meta-analysis of 453 metabolites in 7,824 individuals from two independent European cohorts: TwinsUK (N = 6,056) and KORA (N = 1,768). After selecting SNPs included in both the eQTL and mQTL studies, our analysis was restricted to 7,884 transcripts with ≥ 3 instrumental variables (IVs) (see Methods, Figure 1 —figure supplement 1A) and 242 metabolites with an identifier in The Human Metabolome Database (HMDB) (see Methods, Supplementary File 1a).”

Methods: “mQTL data originate from Shin et al., 2014 which used ultra-high performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) to measure 486 whole blood metabolites in 7,824 European individuals. Association analyses were carried out on ~2.1 million SNPs and are available for 453 metabolites at the Metabolomics GWAS Server (http://metabolomics.helmholtz-muenchen.de/gwas/). Among these metabolites, 242 were manually annotated with Human Metabolome Database (HMDB) identifiers (Supplementary File 1a) and used in this study.”

(4) The authors mention "overall, 83% of the involved genes were causally influencing the level of a single metabolite, while TMEM258 and FADS2 affected 12 metabolites". Firstly I would mention that they affect 12 metabolites each just to be clear – although it's clearly annotated in Table S1.

This has been clarified in the revised manuscript:

“Most involved genes (86%; 83/96) were causally influencing the level of a single metabolite, with some notable exceptions acting as mQTL hubs, such as TMEM258 and FADS2 both affecting the same 11 metabolites, followed by FADS1 affecting a subset of 6 metabolites.”

Second, and more relevantly, it would be nice to know in this table how many of these hits already are known in the literature. For instance I see the most significant hit is MAF1 with 5-oxoproline which doesn't immediately turn up anything (albeit I did not look very long). Hit #7, for SLC22A4 and isovalerylcarnitine does appear to be known in the literature (http://www.metabolomix.com/rs272889-slc22a4-octn1-with-acylcarnitines/ ). Doing a decent literature analysis on all 257 of these hits in Table S1 might take a while (and/or require someone who knows how to do automated literature mining, which isn't trivial for this type of task especially given metabolite names like "X-12696") but even just some cherry-picking and/or (better) systematic analysis of the top 20-50 would be nice to see. It greatly improved my confidence in the findings when I just looked into two of your hits just now and found that one of them appears to be reasonably well known by literature.

This is a great suggestion which led us to conduct a more detailed follow-up investigation. Having a rough estimation of the number of associations that have previously been described in the literature would indeed be of great interest. Due to the large number of identified associations, performing a thorough, manual investigation of each of the identified associations is hardly manageable. As an alternative, we used an automated literature review strategy that we briefly describe below.

First, we reduced the number of investigated metabolites to 242 interpretable compounds for which we could identify an HMDB identifier through manual annotation (see reply to scientific comment #3). Second, using this annotation, we identified transcript-to-metabolite effects for which an interaction between the encoded protein and metabolite was reported in the HMDB. In parallel, we performed an automated literature review in PubMed for all identified transcript-to-metabolite effects. Next, for the 37 causal transcript-metabolite-phenotype triplets for which either the HMDB or the PubMed search for the transcript-to-metabolite yielded at least one results, we (i) added a third term to the PubMed search relating to the associated phenotype and (ii) replaced the metabolite term by the phenotype term to determine whether the transcript-to-phenotype analysis had previously been reported. Owing to the low specificity of the phenotype search term (e.g., “creatinine” appearing in “X mmol/mol creatinine”), we manually reviewed the abstracts of retained publications, excluding findings in which the search terms were used in an unintended way. We describe the methods and results of these analyses in detail in the revised manuscript – along with a discussion of the limitations of our approach – and provide the PMIDs of our automated search in the Supplementary File 1b:

Methods:

“Omics and traits summary statistics

[…] mQTL data originate from Shin et al. 2014, which used ultra-high performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) to measure 486 whole blood metabolites in 7,824 European individuals. Association analyses were carried out on ~2.1 million SNPs and are available for 453 metabolites at the Metabolomics GWAS Server (http://metabolomics.helmholtz-muenchen.de/gwas/). Among these metabolites, 242 were manually annotated with Human Metabolome Database (HMDB) identifiers (Supplementary File 1a) and used in this study. GWAS summary statistics for 28 outcome traits measured in the UK Biobank (UKB) originate from the Neale Lab (http://www.nealelab.is/uk-biobank/). Protein interactions with metabolites were downloaded from HMDB v5.0 (https://hmdb.ca/downloads/hmdb_proteins.xml) and were used to annotate transcript-metabolites associations.

Automated literature review

An automated literature review of all transcript-metabolites associations (Supplementary File 1b) was conducted in PubMed (20/09/2022) following the scheme:${\displaystyle (Gene[Title/Abstract])AND(({\text{Met}}_{\text{HMDB}}[Title/Abstract])OR({\text{Met}}_{\text{Shin}}[Title/Abstract]))}$ With ${\displaystyle \text{Gene}}$ being the name of the gene whose transcript is involved in the association, ${\displaystyle {\text{Met}}_{\text{HMDB}}}$ the involved metabolite’s common name on HMDB, and ${\displaystyle {\text{Met}}_{Sh\text{in}}}$ the involved metabolite’s name as reported in Shin et al., [6]. Returned PubMed identifiers were retrieved (Supplementary File 1b).

For transcript-metabolite associations involved in a causal triplet and for which the transcript-metabolites returned at least one publication (Supplementary File 1f), the search was extend by (i) adding a further search term for the trait (i.e., ${\displaystyle \text{AND}(\text{trait}[\text{Title}/\text{Abstract}])}$) and (ii) substituting the metabolite term for the trait term. Returned PubMed identifiers were retrieved and corresponding abstracts were manually curated to exclude abstracts in which the search terms were used in a meaning other than the intended one (Supplementary File 1f).”

Results:

“Mapping the transcriptome onto the metabolome

[…] While only 5 (3.8%) of the 133 associations were reported in HMDB, an automated literature review (see Methods) identified a match for 22 (16.5%) of the identified transcript-metabolite pairs (Supplementary File 1b).

Mapping metabolome-mediated effects of the transcriptome onto complex phenotypes

[…] Among the 37 triplets for which the transcript and metabolite had previously been linked through automated literature review, none remained after incorporating a third term for the phenotype in the search and manually removing abstracts for which the search terms were used in an erroneous context. Relaxing the search criteria by omitting the metabolite term, 13/37 (35%) triplets returned at least one match for the gene-trait association.”

Discussion:

“Our automated literature review indicates that while some detected associations were previously reported, a large fraction, especially among the triplets, appears to be novel. It should be noted that the number of previously reported associations is likely underestimated as our approach does not account for all synonyms of a given feature and requires the terms to appear in the title or abstract of the publication. This makes it more likely for hypothesis-driven studies, inherently biased towards well-studied genes and metabolites, to be identified. Conversely, high-throughput, hypothesis-free studies that report the given association in a supplemental table are likely to be missed. Furthermore, due to its automated nature, our search is context-blind, so that some of the identified studies might report negative results, associations only under specific conditions (e.g., different organisms, experimental settings), or usage of the search term with a different meaning. To attenuate the latter, we also performed manual review of the retained abstracts for transcript-to-phenotype searches. While flawed, this rough estimate of the amount of existing evidence supporting our findings can be interpreted in combination with other lines of evidence. For instance, […].”

(5) The authors mention that only 33% of the transcript-metabolite-phenotype triplets were observed by TWMR, but that 91% of those observed by all are in concordant directionality. I guess that means for those with concordant directionality, if no significant association was observed in one layer, there was just no cross-layer comparison? i.e. that the 91% "concordant directionality" are only checking among the 33% of hits by TWMR that were also observed both within and across layer?

In response to Reviewer #2’s comments, the current version of the manuscript uses more stringent significant thresholds. At FDR 5%, we now observe 216 transcript-metabolite-phenotype triplets. Among these 216 triplets, 90 (42%) had a significant transcript-to-phenotype TWMR effect. Among these 90 triplets, 84 (93%) showed a directionally concordant effect between the transcript-to-phenotype, transcript-to-metabolite, and metabolite-to-phenotype estimates, meaning that the sign of the product of the transcript-to-metabolite and metabolite-to-phenotype effects equals the sign of the transcript-to-phenotype effect. This is now clarified in the Discussion:

“[…] among the 90 signals that were also identified through TWMR, 93% showed a directionally concordant effect between the transcript-to-phenotype, transcript-to-metabolite, and metabolite-to-phenotype estimates (i.e., sign of the product of the transcript-to-metabolite and metabolite-to-phenotype effects agrees with the sign of the transcript-to-phenotype effect).”

(6) The authors mention no availability of trans-eQTL data; is that because the sample size is not sufficient to do a proper trans-eQTL analysis, or is it because of data privacy and an inability to get the precomputed trans-eQTLs from the eQTLGen Consortium due to some potential conflict (GDPR?).

The trans-eQTLs data from eQTLGen Consortium are restricted to ~10,000 SNPs reported to be associated with complex traits in GWAS catalog due to the massive size of this type of data. There are current efforts within the eQTLGen Consortium to meta-analyse trans-eQTL summary statistics genome-wide, but these data will not be available before Q3 of 2023. Currently available data is unfortunately not sufficient to investigate reverse causality on metabolites and gene expression, as specified in the manuscript’s Discussion:

“Due to the lack of genome-wide trans-eQTL association summary statistics, our method does not investigate reverse causality on metabolites and gene expression […].”

(7) Back to Table S1: There are 257 hits here, with 191 genes (out of 7883) and 154 metabolites (out of 453), meaning that 299 metabolites didn't even have one significant association?

This is correct. In the current version of the manuscript and as explained in the reply to previous points, we restricted our analysis to 242 metabolites manually annotated with HMDB identifiers. In Supplementary File 1b, we report 133 significant (FDR < 5%) transcript-to-metabolite effects involving 75 unique metabolites, so that in case of 167 metabolites we do not have evidence for them being impacted by any of the 7,824 tested genes. It is interesting to notice that despite changes in our analysis (i.e., restrict to interpretable metabolites, more stringent criterion for significance), the proportion of metabolites that are not affected by any transcripts remains similar (previously: 299/453 = 66%; currently: 167/242 = 69%). One potential explanation is lack of power, given the relatively small sample size of the mQTL dataset (< 8,000), as compared to GWASs for other complex traits (> hundreds of thousands), that prevents detection of associations for a large fraction of metabolites. Alternatively, these metabolites might be influenced by transcripts that cannot be instrumented, making it impossible to detect association.

Reviewer #2 (Recommendations for the authors):

1. In general, the methodological approach should be explained in more detail. The study presents a novel MR-based approach and thus, the discussion of the method basics should be elaborated. This pertains to the following points:

– In the definition of causal effects from molecular traits on the outcome (equation (1)), the symbol β' is not explained. The method used for standardization of effect sizes is not elaborated.

We clarified this in the current version of the manuscript, which now reads as follows:

“Here, ${\displaystyle \beta }$ is a vector of length n containing the standardized effect size of n independent SNPs on the gene/metabolite, derived from eQTL/mQTL studies, with ${\displaystyle {\beta }^{\prime }}$ being the transpose of ${\displaystyle \beta }$. ${\displaystyle \gamma }$ is a vector of length n containing the standardized effect size of each SNP on the outcome. C is the pairwise LD matrix between the n SNPs. The standardized effect sizes for molecular and outcome GWASs were obtained from Z-score of summary statistics standardized by the square root of the sample size to be on the same standard deviation scale.”

– For the multivariable MR analysis (equation not numbered), no causal effect statistic is specified. As discussed in the public review, no significance testing has been performed, but we strongly urge the authors to do so or give reasons for this lack of testing in the case of causal triplets.

The equation has now been labeled (6). As specified in the Methods, we used multivariable MR (MVMR) only to estimate the direct effect of gene expression on the phenotype. We do not test the significance of this estimate as we only compare it with the total effect estimated by TWMR regardless its significance. Regarding significance testing for other sections of the manuscript, please see extensive reply to major concern #1 and #2.

– In the description of the simulation analyses, it should be made more clear how the definition of βGWAS as a function of βmQTL is reached.

As described in the Methods, to estimate the effect of the metabolite on the phenotype, we considered a metabolite with heritability ${\displaystyle h_M^2=0.04}$ (corresponding to the median of ${\displaystyle h^2}$ in the KORA+TwinsUK mQTL data) and ${\displaystyle N_{\text{IVs}}=5}$ (corresponding to the median number of IVs used in MWMR analyses). Effect size ${\displaystyle {\beta }_{\text{mQTL}}}$ are from a normal distribution ${\displaystyle {\beta }_{\text{mQTL}}\sim N(0,\frac{h_M^2}{N_{\text{IVs}}})}$. These quantities allowed to define ${\displaystyle {\beta }_{\text{GWAS}}}$ as ${\displaystyle {\beta }_{\text{GWAS}}=\sqrt{{\alpha }_{\text{TP}}\ast \left(1-\rho \right)/\sigma }\ast {\beta }_{\text{mQTL}}+{\epsilon }_P}$, where ${\displaystyle {\epsilon }_P\sim N(0,\frac{1}{N_{\text{GWAS}}})}$.

– The authors should include a thorough discussion and a description of the method employed to evaluate the comparison between total and direct transcript-phenotype effects in Figure 2B.

We have now added a description of the comparison of the direct and total effects in the Methods:

“The proportion of direct effect ${\displaystyle (\rho )}$ is calculated by regressing direct effects (${\displaystyle {\alpha }_d}$) on total effects (${\displaystyle {\alpha }_{\text{TP}}}$) and then correcting for regression dilution bias:

${\displaystyle {\rho }_{\text{corrected}}=\frac{\rho }{\sqrt{1-\frac{\sum _{}^{}{(SE\left({\alpha }_{\text{TP}}\right))}^2}{\sum _{}^{}{\alpha }_{\text{TP}}^2}}}}$ ”

It is crucial that due to lack of power, we do not compare individual direct-vs-total effects but only assess overall direct-vs-total effects across all transcript-phenotype pairs. The central message of the paper is not the difference between the two, but the difference in detection power between the two strategies, as we do illustrate in Figure 4.

2. There is some apparent inconsistency in the interpretation of results from the empirical study: The abstract states that 67 transcript-phenotype effects present within the 206 causal triplets were missed in a TWMR-only analyses, while the paragraph "Power Analysis" states that only 67 triplets showed a significant effect estimated by TWMR.

We apologize for this inconsistency, which has been clarified. The revised sections now read:

Abstract: “We identified 216 transcript-metabolite-trait causal triplets involving 26 medically relevant phenotypes. Among these associations, 58% were missed by classical transcriptome-wide MR […].”

Results: “Importantly, only 42% (90/216) of the causal triplets showed a significant total transcript-to-phenotype effect (i.e., estimated by TWMR), suggesting that the method lacks power under current settings.”

The text relating to Figure 2B states that 79 % (~163) of the triplets were estimated to be due to direct effects by the transcript only. This should be presented in a more consistent and clear manner.

In Figure 2B, we plot the direct (${\displaystyle {\alpha }_d}$; y-axis) and total (${\displaystyle {\alpha }_{\text{TP}}}$; x-axis) effects for the 216 transcript-metabolite-trait causal triplets. Regressing direct effects on total effects, we observed that 77% of the transcript effect on the phenotype is direct and thus not mediated by the metabolites. This is explained in the Results and Methods of the manuscript:

Results:

“Regressing direct effects (${\displaystyle {\alpha }_d}$) on total effects (${\displaystyle {\alpha }_{\text{TP}}}$) and accounting for regression dilution bias (see Methods; Figure 2A), it was estimated that 77% [95% CI: 70%-85%] of the transcript effect on the phenotype was direct and thus not mediated by the metabolites (Figure 2B).”

Methods:

“The proportion of direct effect ${\displaystyle (\rho )}$ is calculated by regressing direct effects (${\displaystyle {\alpha }_d}$) on total effects (${\displaystyle {\alpha }_{\text{TP}}}$) and then correcting for regression dilution bias:

${\displaystyle {\rho }_{\text{corrected}}=\frac{\rho }{\sqrt{1-\frac{\sum _{}^{}{(SE\left({\alpha }_{\text{TP}}\right))}^2}{\sum _{}^{}{\alpha }_{\text{TP}}^2}}}}$ ”

3. Plots for the distribution of the number of IVs per transcript/metabolite should be included (as supplemental figures).

We now include Figure 1 —figure supplement 1 to illustrate this.

We furthermore refer to these plots in the Results section:

“Mapping the transcriptome onto the metabolome

[…] After selecting SNPs included in both the eQTL and mQTL studies, our analysis was restricted to 7,884 transcripts with ≥ 3 instrumental variables (IVs) (see Methods, Figure 1 —figure supplement 1A) […].

Mapping the metabolome onto complex phenotypes

Univariable metabolome-wide MR (MWMR; Figure 1A) was used to identify causal relationships between 48 metabolites with ≥ 3 IVs and 28 complex phenotypes (Figure 1 —figure supplement 1B). […]

Mapping the transcriptome onto complex phenotypes

We applied univariable transcriptome-wide MR (TWMR [12]; Figure 1B) to identify associations between expression levels of 10,435 transcripts from the eQTLGen Consortium with ≥ 3 IVs (Figure 1 —figure supplement 1C) measured in both exposure and outcome datasets and the same 28 UKB phenotypes described in the previous section (Supplementary File 1c)[…].”

4. It would be interesting to discuss potential modes of mediation between transcripts, metabolites and phenotypes. For example, individual metabolites might integrate the effects of several transcripts. This would not result in triplets, but ultimately, networks of mediation as shown for the TMEM258, FADS1, FADS2 case. How would such cases be treated in a statistically sound manner?

This is a very interesting idea and we believe that the field is increasingly interested in modelling such complex biological network. The situation where a single metabolite mediates the effect of multiple transcripts only needs a systematic search for metabolites repeatedly occurring in triplets. We identified 12 such metabolites (see Table integrated in revised manuscript as Supplementary File 1g). Nine out of 12 metabolites integrate effect of multiple transcripts to in turn affects multiple phenotypes, as illustrated in the last column of the Table. Of note 8/12 cases, the metabolites integrate the effect of transcript whose genes are in close genomic proximity (shown as bold transcripts in the same color). We already made this observation in the original manuscript and discussed this limitation:

“Owing to linkage disequilibrium and regulatory variants affecting multiple genes, transcripts from adjacent genes might appear to be involved in the same signals, as exemplified with the TMEM258/FADS1/FADS2 locus (Figure 3). While literature supports the role of the FADS genes, one cannot exclude a role for TMEM258, nor disentangle the specific function of FADS1 and FADS2.”

Investigating more complex scenarios where multiple correlated metabolites mediate the effect of a single transcript would require more advanced modelling techniques that go beyond the scope of this study. We strengthened the Discussion on the topic in the revised manuscript which now reads as follows:

“Another consideration is that complex phenotypes can have a stronger impact on gene expression than the opposite [15]. Due to the lack of genome-wide trans-eQTL association summary statistics, our method does not investigate reverse causality on metabolites and gene expression, nor the role of metabolites as regulators of gene expression. Metabolites might also integrate the effect of several transcripts (i.e., multiple transcripts causally impact the levels of the same metabolite) before affecting complex phenotypes (Supplementary File 1g) or multiple metabolites may jointly mediate the impact of a single transcript. Modelling the latter phenomenon, which is beyond the scope of our current work, requires the development of Structural Equation Models accounting for such effects and will eventually lead to a more comprehensive modelling of causal relations in complex biological networks, nuancing the interpretation of the molecular mechanisms shaping complex traits.”

Author details

1. Chiara Auwerx

   1. Center for Integrative Genomics, University of Lausanne, Lausanne, Switzerland
   2. Swiss Institute of Bioinformatics, Lausanne, Switzerland
   3. University Center for Primary Care and Public Health, Lausanne, Switzerland
   4. Department of Computational Biology, University of Lausanne, Lausanne, Switzerland

   Contribution

   Writing – original draft, Investigation, Visualization

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0003-3613-8450

2. Marie C Sadler

   1. Swiss Institute of Bioinformatics, Lausanne, Switzerland
   2. University Center for Primary Care and Public Health, Lausanne, Switzerland
   3. Department of Computational Biology, University of Lausanne, Lausanne, Switzerland

   Contribution

   Formal analysis, Writing – review and editing

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0002-2599-9207

3. Tristan Woh

   Department of Computational Biology, University of Lausanne, Lausanne, Switzerland

   Contribution

   Formal analysis, Writing – review and editing

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0001-6916-0174

4. Alexandre Reymond

   Center for Integrative Genomics, University of Lausanne, Lausanne, Switzerland

   Contribution

   Writing – review and editing

   Competing interests

   No competing interests declared

5. Zoltán Kutalik

   1. Swiss Institute of Bioinformatics, Lausanne, Switzerland
   2. University Center for Primary Care and Public Health, Lausanne, Switzerland
   3. Department of Computational Biology, University of Lausanne, Lausanne, Switzerland

   Contribution

   Conceptualization, Writing – original draft

   For correspondence

   zoltan.kutalik@unil.ch

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0001-8285-7523

6. Eleonora Porcu

   1. Center for Integrative Genomics, University of Lausanne, Lausanne, Switzerland
   2. Swiss Institute of Bioinformatics, Lausanne, Switzerland
   3. University Center for Primary Care and Public Health, Lausanne, Switzerland

   Contribution

   Conceptualization, Formal analysis, Writing – original draft

   For correspondence

   eleonora.porcu@unil.ch

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0003-2878-7485

• 2,357

  Page views

• 244

  Downloads

• 4

  Citations

Article citation count generated by polling the highest count across the following sources: PubMed Central, Crossref, Scopus.