Skip to contents

Given harmonised data, this function conducts a two-sample MR analysis.

Usage

run_TwoSampleMR(TwoSampleMRinput, mr_plot = "None", prefix = "")

Arguments

TwoSampleMRinput

Harmonised data.

mr_plot

one of "None", "TwoSampleMR", "pQTLtools" for no, the original and the revised plots, respectively.

prefix

a prefix for output files.

Value

No value is returned but several files.

Details

As TwoSampleMR faces seemingly perplexing options, this function intends to simplify various steps in a two-sample MR as in Dimou and Tsilidis (2018) . It is particularly useful when a large numbher of MRs are necessary, e.g., multiple proteins and their cis/trans regions need to be examined, in which case prefix could direct the output to various directories.

References

Dimou NL, Tsilidis KK (2018). “A Primer in Mendelian Randomization Methodology with a Focus on Utilizing Published Summary Association Data.” In Evangelou E (ed.), Genetic Epidemiology: Methods and Protocols, chapter 13, 211-230. Springer New York, New York, NY. ISBN 978-1-4939-7868-7, doi:10.1007/978-1-4939-7868-7_13 .

Examples

library(TwoSampleMR)
library(pQTLtools)
outcomes <- "ebi-a-GCST007432"
prot <- "MMP.10"
type <- "cis"
f <- paste0(prot,"-",type,".mrx")
d <- read.table(file.path(find.package("pQTLtools",lib.loc=.libPaths()),"tests",f),
                header=TRUE)
exposure <- TwoSampleMR::format_data(within(d,{P=10^logP}), phenotype_col="prot", snp_col="rsid",
                                     chr_col="Chromosome", pos_col="Posistion",
                                     effect_allele_col="Allele1", other_allele_col="Allele2",
                                     eaf_col="Freq1", beta_col="Effect", se_col="StdErr",
                                     pval_col="P", log_pval=FALSE,
                                     samplesize_col="N")
clump <- TwoSampleMR::clump_data(exposure)
#> Please look at vignettes for options on running this locally if you need to run many instances of this command.
#> Clumping vmQoON, 1106 variants, using EUR population reference
#> Server code: 502; Server is possibly experiencing traffic, trying again...
#> Retry succeeded!
#> Removing 1102 of 1106 variants due to LD with other variants or absence from LD reference panel
outcome <- TwoSampleMR::extract_outcome_data(snps=exposure$SNP,outcomes=outcomes)
#> Extracting data for 1106 SNP(s) from 1 GWAS(s)
#> Finding proxies for 155 SNPs in outcome ebi-a-GCST007432
#> Extracting data for 155 SNP(s) from 1 GWAS(s)
#> Server code: 502; Server is possibly experiencing traffic, trying again...
#> Retry succeeded!
harmonise <- TwoSampleMR::harmonise_data(clump,outcome)
#> Harmonising MMP.10 (vmQoON) and FEV1 || id:ebi-a-GCST007432 (ebi-a-GCST007432)
prefix <- paste(outcomes,prot,type,sep="-")
run_TwoSampleMR(harmonise, mr_plot="pQTLtools", prefix=prefix)
#> Analysing 'vmQoON' on 'ebi-a-GCST007432'

#> Warning: Removed 1 rows containing missing values (`geom_errorbarh()`).
#> Warning: Removed 1 rows containing missing values (`geom_point()`).


#> Warning: Removed 1 rows containing missing values (`geom_errorbarh()`).
#> Warning: Removed 1 rows containing missing values (`geom_point()`).

caption <- "Table. MMP.10 variants and FEV1"
knitr::kable(read.delim(paste0(prefix,"-result.txt"),header=TRUE),
             caption=paste(caption, "(result)"))
#> 
#> 
#> Table: Table. MMP.10 variants and FEV1 (result)
#> 
#> |id.exposure |id.outcome       |outcome                               |exposure |method                    | nsnp|         b|      se|  pval|
#> |:-----------|:----------------|:-------------------------------------|:--------|:-------------------------|----:|---------:|-------:|-----:|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |MR Egger                  |    4| -0.000577| 0.01001| 0.959|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |Weighted median           |    4|  0.002002| 0.00625| 0.749|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |Inverse variance weighted |    4|  0.002067| 0.00594| 0.728|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |Simple mode               |    4|  0.004152| 0.00945| 0.690|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |Weighted mode             |    4|  0.001134| 0.00666| 0.876|
knitr::kable(read.delim(paste0(prefix,"-heterogeneity.txt"),header=TRUE),
             caption=paste(caption,"(heterogeneity)"))
#> 
#> 
#> Table: Table. MMP.10 variants and FEV1 (heterogeneity)
#> 
#> |id.exposure |id.outcome       |outcome                               |exposure |method                    |     Q| Q_df| Q_pval|
#> |:-----------|:----------------|:-------------------------------------|:--------|:-------------------------|-----:|----:|------:|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |MR Egger                  | 0.114|    2|  0.944|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |Inverse variance weighted | 0.222|    3|  0.974|
knitr::kable(read.delim(paste0(prefix,"-pleiotropy.txt"),header=TRUE),
             caption=paste(caption,"(pleiotropy)"))
#> 
#> 
#> Table: Table. MMP.10 variants and FEV1 (pleiotropy)
#> 
#> |id.exposure |id.outcome       |outcome                               |exposure | egger_intercept|      se|  pval|
#> |:-----------|:----------------|:-------------------------------------|:--------|---------------:|-------:|-----:|
#> |vmQoON      |ebi-a-GCST007432 |FEV1 &#124;&#124; id:ebi-a-GCST007432 |MMP.10   |         0.00143| 0.00434| 0.774|
knitr::kable(read.delim(paste0(prefix,"-single.txt"),header=TRUE),
             caption=paste(caption,"(single)"))
#> 
#> 
#> Table: Table. MMP.10 variants and FEV1 (single)
#> 
#> |exposure |outcome                               |id.exposure |id.outcome       | samplesize|SNP                             |         b|      se|     p|
#> |:--------|:-------------------------------------|:-----------|:----------------|----------:|:-------------------------------|---------:|-------:|-----:|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs11225354                      |  0.006425| 0.01217| 0.598|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs142915220                     | -0.003993| 0.02856| 0.889|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs17099622                      |  0.004360| 0.02229| 0.845|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs17860955                      |  0.000616| 0.00739| 0.934|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|All - Inverse variance weighted |  0.002067| 0.00594| 0.728|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|All - MR Egger                  | -0.000577| 0.01001| 0.959|
knitr::kable(read.delim(paste0(prefix,"-loo.txt"),header=TRUE),
             caption=paste(caption,"(loo)"))
#> 
#> 
#> Table: Table. MMP.10 variants and FEV1 (loo)
#> 
#> |exposure |outcome                               |id.exposure |id.outcome       | samplesize|SNP         |        b|      se|     p|
#> |:--------|:-------------------------------------|:-----------|:----------------|----------:|:-----------|--------:|-------:|-----:|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs11225354  | 0.000703| 0.00681| 0.918|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs142915220 | 0.002341| 0.00608| 0.700|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs17099622  | 0.001891| 0.00617| 0.759|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|rs17860955  | 0.004730| 0.01001| 0.636|
#> |MMP.10   |FEV1 &#124;&#124; id:ebi-a-GCST007432 |vmQoON      |ebi-a-GCST007432 |     321047|All         | 0.002067| 0.00594| 0.728|
for (x in c("result","heterogeneity","pleiotropy","single","loo"))
    unlink(paste0(prefix,"-",x,".txt"))