If you use the MungeSumstats package, please cite
The MungeSumstats package is designed to facilitate the standardisation of GWAS summary statistics as utilised in our Nature Genetics paper1.
The package is designed to handle the lack of standardisation of output files by the GWAS community. There is a group who have now manually standardised many GWAS: R interface to the IEU GWAS database API • ieugwasr and gwasvcf but because a lot of GWAS remain closed access, these repositories are not all encompassing.
The GWAS-Download project has collated summary statistics from 200+ GWAS. This repository has been utilsed to identify the most common formats, all of which can be standardised with MungeSumstats.
Moreover, there is an emerging standard of VCF format for summary statistics files with multiple, useful, associated R packages such as vcfR. However, there is currently no method to convert VCF formats to a standardised format that matches older approaches.
The MungeSumstats package standardises both VCF and the most common summary statistic file formats to enable downstream integration and analysis.
MungeSumstats also offers comprehensive Quality Control (QC) steps which are important prerequisites for downstream analysis like Linkage disequilibrium score regression (LDSC) and MAGMA.
Moreover, MungeSumstats is efficiently written resulting in all reformatting
and quality control checks completing in minutes for GWAS summary statistics
with 500k SNPs on a standard desktop machine. This speed can be increased
further by increasing the number of threads (nThread) for data.table
to use.
Currently MungeSumstats only works on data from humans, as it uses human-based genome references.
MungeSumstats will ensure that the all essential columns for analysis are present and syntactically correct. Generally, summary statistic files include (but are not limited to) the columns:
MungeSumstats uses a mapping file to infer the inputted column names (run
data("sumstatsColHeaders")
to view these). This mapping file is far more
comprehensive than any other publicly available munging tool containing more
than 200 unique mappings at the time of writing this vignette. However, if your
column headers are missing or if you want to change the mapping, you can do so
by passing your own mapping file (see format_sumstats(mapping_file)
).
MungeSumstats offers unmatched levels of quality control to ensure, for example, consistency of allele assignment and direction of effects. Tests run by MungeSumstats include:
Users can specify which checks to run on their data. A note on the allele flipping check: MungeSumstats infers the effect allele will always be the A2 allele, this is the approach done for IEU GWAS VCF and has such also been adopted here. This inference is first from the inputted file’s column headers however, the allele flipping check ensures this by comparing A1, what should be the reference allele, to the reference genome. If a SNP’s A1 DNA base doesn’t match the reference genome but it’s A2 (what should be the alternative allele) does, the alleles will be flipped along with the effect information (e.g. Beta, Odds Ratio, signed summary statistics, FRQ, Z-score*).
*-by default the Z-score is assumed to be calculated off the effect size not the P-value and so will be flipped. This can be changed by a user.
If a test is failed, the user will be notified and if possible, the
input will be corrected. The QC steps from the checks above can also be
adjusted to suit the user’s analysis, see
MungeSumstats::format_sumstats
.
MungeSumstats can handle VCF, txt, tsv, csv file types or .gz/.bgz versions of these file types. The package also gives the user the flexibility to export the reformatted file as tab-delimited, VCF or R native objects such as data.table, GRanges or VRanges objects. The output can also be outputted in an LDSC ready format which means the file can be fed directly into LDSC without the need for additional munging.
The MungeSumstats package contains small subsets of GWAS summary statistics files. Firstly, on Educational Attainment by Okbay et al 2016: PMID: 27898078 PMCID: PMC5509058 DOI: 10.1038/ng1216-1587b.
Secondly, a VCF file (VCFv4.2) relating to the GWAS Amyotrophic lateral sclerosis from ieu open GWAS project. Dataset: ebi-a-GCST005647: https://gwas.mrcieu.ac.uk/datasets/ebi-a-GCST005647/
These datasets will be used to showcase MungeSumstats functionality.
MungeSumstats is available on Bioconductor. To install the package on Bioconductor run the following lines of code:
if (!require("BiocManager"))
install.packages("BiocManager")
BiocManager::install("MungeSumstats")
Once installed, load the package:
library(MungeSumstats)
To standardise the summary statistics’ file format, simply call
format_sumstats()
passing in the path to your summary statistics file or
directly pass the summary statistics as a dataframe or datatable.
You can specify which genome build was used in the GWAS(GRCh37 or
GRCh38) or, as default, infer the genome build from the data.The
reference genome is used for multiple checks like deriving missing data
such SNP/BP/CHR/A1/A2 and for QC steps like removing non-biallelic SNPs,
strand-ambiguous SNPs or ensuring correct allele and direction of SNP effects.
The path to the reformatted summary statistics file can be returned by the
function call, the user can specify a location to save the file or the user can
return an R native object for the data: data.table, VRanges or GRanges
object.
Note that for a number of the checks implored by MungeSumstats a
reference genome is used. If your GWAS summary statistics file of
interest relates to GRCh38, you will need to install
SNPlocs.Hsapiens.dbSNP144.GRCh38
and BSgenome.Hsapiens.NCBI.GRCh38
from Bioconductor as follows:
BiocManager::install("SNPlocs.Hsapiens.dbSNP144.GRCh38")
BiocManager::install("BSgenome.Hsapiens.NCBI.GRCh38")
If your GWAS summary statistics file of interest relates to GRCh37,
you will need to install SNPlocs.Hsapiens.dbSNP144.GRCh37
and
BSgenome.Hsapiens.1000genomes.hs37d5
from Bioconductor as follows:
BiocManager::install("SNPlocs.Hsapiens.dbSNP144.GRCh37")
BiocManager::install("BSgenome.Hsapiens.1000genomes.hs37d5")
These may take some time to install and are not included in the package as some users may only need one of GRCh37/GRCh38.
The Educational Attainment by Okbay GWAS summary statistics file is saved as a text document in the package’s external data folder so we can just pass the file path to it straight to MungeSumstats.
NOTE - By default, Formatted results will be saved to tempdir()
. This
means all formatted summary stats will be deleted upon ending the R session if
not copied to a local file path. Otherwise, to keep formatted summary stats,
change save_path
( e.g.file.path('./formatted',basename(path))
), or make
sure to copy files elsewhere after processing
( e.g.file.copy(save_path, './formatted/' )
.
eduAttainOkbayPth <- system.file("extdata","eduAttainOkbay.txt",
package="MungeSumstats")
reformatted <-
MungeSumstats::format_sumstats(path=eduAttainOkbayPth,
ref_genome="GRCh37")
The arguments format_sumstats
in that control the level of QC
conducted by MungeSumstats are:
p-values < 5e-324
be converted
to 0? Small p-values pass the R limit and can cause errors with
LDSC/MAGMA and should be converted. Default is TRUE.data.table
, GRanges
or VRanges
directly to user.
Otherwise, return the path to the save data. Default is FALSE.data(sumstatsColHeaders)
for default mapping and necessary format.See ?MungeSumstats::format_sumstats()
for the full list of parameters to
control MungeSumstats QC and standardisation steps.
VCF files can also be standardised to the same format as other summary statistic files. A subset of the Amyotrophic lateral sclerosis GWAS from the ieu open GWAS project (a .vcf file) has been added to MungeSumstats to demonstrate this functionality.Simply pass the path to the file in the same manner you would for other summary statistic files:
#save ALS GWAS from the ieu open GWAS project to a temp directory
ALSvcfPth <- system.file("extdata","ALSvcf.vcf", package="MungeSumstats")
#set a low INFO filter so we get a return
reformatted_vcf <-
MungeSumstats::format_sumstats(path=ALSvcfPth,
ref_genome="GRCh37",
INFO_filter=0.1)
You can also get more information on the SNPs which have had data imputed or
have been filtered out by MungeSumstats by using the imputation_ind
and
log_folder_ind
parameters respectively. For example:
#set
reformatted_vcf_2 <-
MungeSumstats::format_sumstats(path=ALSvcfPth,
ref_genome="GRCh37",
INFO_filter = 0.1,
log_folder_ind=TRUE,
imputation_ind=TRUE,
log_mungesumstats_msgs=TRUE)
Note that using log_folder_ind
returns a list from format_sumstats
which
includes the file locations of the differing classes of removed SNPs. Using
log_mungesumstats_msgs
saves the messages to the console to a file which is
returned in the same list. Note that not all the messages will also print to
screen anymore when you set log_mungesumstats_msgs
:
names(reformatted_vcf_2)
## [1] "sumstats" "log_files"
A user can load a file to view the excluded SNPs:
data.table::fread(reformatted_vcf_2$log_files$info_filter)[1:10]
And we can see all these SNPs had lower INFO values that our set threshold.
The different types of exclusion which lead to the names are explained below:
Note to export to another type such as R native objects including data.table,
GRanges, VRanges or save as a VCF file, set return_format
:
#set
reformatted_vcf_2 <-
MungeSumstats::format_sumstats(path=ALSvcfPth,
ref_genome="GRCh37",
INFO_filter = 0.1,
log_folder_ind=TRUE,
imputation_ind=TRUE,
log_mungesumstats_msgs=TRUE,
return_format="GRanges")
See our publication for further discussion of these checks and options:
MungeSumstats also contains a function to quickly infer the genome build of
multiple summary statistic files. This can be called separately to
format_sumstats()
which is useful if you want to just quickly check the genome
build:
# Pass path to Educational Attainment Okbay sumstat file to a temp directory
eduAttainOkbayPth <- system.file("extdata", "eduAttainOkbay.txt",
package = "MungeSumstats")
sumstats_list <- list(ss1 = eduAttainOkbayPth, ss2 = eduAttainOkbayPth)
ref_genomes <- MungeSumstats::get_genome_builds(sumstats_list = sumstats_list)
The MungeSumstats package aims to be able to handle the most common summary statistic file formats including VCF. If your file can not be formatted by MungeSumstats feel free to report the bug on github: https://github.com/neurogenomics/MungeSumstats along with your summary statistic file header.
We also encourage people to edit the code to resolve their particular issues too and are happy to incorporate these through pull requests on github. If your summary statistic file headers are not recognised by MungeSumstats but correspond to one of:
SNP, BP, CHR, A1, A2, P, Z, OR, BETA, LOG_ODDS,
SIGNED_SUMSTAT, N, N_CAS, N_CON, NSTUDY, INFO or FRQ
feel free to update the MungeSumstats::sumstatsColHeaders
following
the approach in the data.R file and add your mapping. Then use a pull
request on github and we will incorporate this change into the package.
A note on MungeSumstats::sumstatsColHeaders
for summary statistic files with
A0/A1. The mapping in MungeSumstats::sumstatsColHeaders
converts A0 to A*,
this is a special case so that the code knows to map A0/A1 to A1/A2 (ref/alt).
The special case is needed since ordinarily A1 refers to the reference not the
alternative allele.
See the Open GWAS vignette for how MungeSumstats can be used along with data from the MRC IEU Open GWAS Project and also Mungesumstats’ functionality to handle lists of summary statistics files.
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.4 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.14-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.14-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] MungeSumstats_1.2.4 BiocStyle_2.22.0
##
## loaded via a namespace (and not attached):
## [1] fs_1.5.2
## [2] bitops_1.0-7
## [3] matrixStats_0.61.0
## [4] bit64_4.0.5
## [5] filelock_1.0.2
## [6] progress_1.2.2
## [7] httr_1.4.2
## [8] GenomeInfoDb_1.30.1
## [9] googleAuthR_2.0.0
## [10] tools_4.1.3
## [11] bslib_0.3.1
## [12] utf8_1.2.2
## [13] R6_2.5.1
## [14] DBI_1.1.2
## [15] BiocGenerics_0.40.0
## [16] tidyselect_1.1.2
## [17] prettyunits_1.1.1
## [18] bit_4.0.4
## [19] curl_4.3.2
## [20] compiler_4.1.3
## [21] cli_3.2.0
## [22] Biobase_2.54.0
## [23] xml2_1.3.3
## [24] DelayedArray_0.20.0
## [25] rtracklayer_1.54.0
## [26] bookdown_0.25
## [27] sass_0.4.1
## [28] rappdirs_0.3.3
## [29] stringr_1.4.0
## [30] digest_0.6.29
## [31] Rsamtools_2.10.0
## [32] rmarkdown_2.13
## [33] R.utils_2.11.0
## [34] XVector_0.34.0
## [35] BSgenome.Hsapiens.1000genomes.hs37d5_0.99.1
## [36] pkgconfig_2.0.3
## [37] htmltools_0.5.2
## [38] MatrixGenerics_1.6.0
## [39] dbplyr_2.1.1
## [40] fastmap_1.1.0
## [41] BSgenome_1.62.0
## [42] rlang_1.0.2
## [43] RSQLite_2.2.11
## [44] jquerylib_0.1.4
## [45] BiocIO_1.4.0
## [46] generics_0.1.2
## [47] jsonlite_1.8.0
## [48] BiocParallel_1.28.3
## [49] dplyr_1.0.8
## [50] R.oo_1.24.0
## [51] VariantAnnotation_1.40.0
## [52] RCurl_1.98-1.6
## [53] magrittr_2.0.2
## [54] GenomeInfoDbData_1.2.7
## [55] Matrix_1.4-1
## [56] Rcpp_1.0.8.3
## [57] S4Vectors_0.32.4
## [58] fansi_1.0.2
## [59] lifecycle_1.0.1
## [60] R.methodsS3_1.8.1
## [61] stringi_1.7.6
## [62] yaml_2.3.5
## [63] SummarizedExperiment_1.24.0
## [64] zlibbioc_1.40.0
## [65] BiocFileCache_2.2.1
## [66] grid_4.1.3
## [67] blob_1.2.2
## [68] parallel_4.1.3
## [69] crayon_1.5.0
## [70] lattice_0.20-45
## [71] Biostrings_2.62.0
## [72] GenomicFeatures_1.46.5
## [73] hms_1.1.1
## [74] KEGGREST_1.34.0
## [75] knitr_1.37
## [76] pillar_1.7.0
## [77] GenomicRanges_1.46.1
## [78] rjson_0.2.21
## [79] biomaRt_2.50.3
## [80] stats4_4.1.3
## [81] XML_3.99-0.9
## [82] glue_1.6.2
## [83] evaluate_0.15
## [84] SNPlocs.Hsapiens.dbSNP144.GRCh37_0.99.20
## [85] data.table_1.14.2
## [86] BiocManager_1.30.16
## [87] png_0.1-7
## [88] vctrs_0.3.8
## [89] purrr_0.3.4
## [90] assertthat_0.2.1
## [91] cachem_1.0.6
## [92] xfun_0.30
## [93] restfulr_0.0.13
## [94] gargle_1.2.0
## [95] tibble_3.1.6
## [96] GenomicAlignments_1.30.0
## [97] AnnotationDbi_1.56.2
## [98] memoise_2.0.1
## [99] IRanges_2.28.0
## [100] ellipsis_0.3.2
1. Nathan G. Skene, T. E. B., Julien Bryois. Genetic identification of brain cell types underlying schizophrenia. Nature Genetics (2018). doi:10.1038/s41588-018-0129-5