Skip to content

WuYang-Lab/GWAStoolkit

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

48 Commits
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

GWAStoolkit

A unified, high-performance C++ toolkit for processing GWAS summary statistics.

GWAStoolkit

GitHub tag GitHub release License Monthly commits Stars

GWAStoolkit integrates common GWAS summary-statistic processing steps into one flag-driven command-line tool. There are no subcommands: the program decides what to run from the flags you provide.


Highlights

  • Fast rsID annotation from dbSNP text, gzipped text, or PLINK .bim.
  • CHR/POS lookup with allele-aware matching: GWAS A1/A2 can match dbSNP REF/ALT in either order.
  • OR to beta/SE conversion.
  • Sample-size filling from fixed case/control, per-row case/control columns, or total-N.
  • Optional beta/SE standardization through --std-effect.
  • Final format conversion to gwas, cojo, popcorn, or mrmega.
  • Built-in QC for available beta, SE, freq, P, and N columns.
  • Optional duplicate SNP removal by smallest P-value.
  • Quiet default logging, with --verbose for column mappings, scan progress, and debug examples.

Execution Model

Operations are selected only by flags:

  • --or COL converts OR to beta/SE.
  • --case --control, --case-col --control-col, or numeric --N fills sample size.
  • --std-effect standardizes beta/SE using binary effective N or continuous total N.
  • --dbsnp FILE annotates/fills rsID.
  • --format writes the final output in a supported downstream format.

When several flags are used together, GWAStoolkit runs them in a fixed order:

OR -> sample size / standardized effect -> rsID -> final format conversion

Input columns are required only when the requested operation needs them. For example, CHR and POS are required only when --dbsnp is used. Beta/SE may be absent if --or is used to generate them before final format conversion.

Installation

Dependencies:

  • g++ with C++17 support
  • zlib

Build:

make clean
make

This creates:

./GWAStoolkit

Show help:

./GWAStoolkit --help

Quick Start

OR to beta/SE, fixed total sample size, rsID annotation, and COJO output in one call:

./GWAStoolkit \
  --gwas-summary gwas_or.txt.gz \
  --out final.cojo.rsid.txt.gz \
  --or OR \
  --case 20000 \
  --control 30000 \
  --dbsnp GRCH37.dbSNP157.txt \
  --dbchr CHROM --dbpos POS --dbA1 REF --dbA2 ALT --dbrsid ID \
  --chr CHR --pos POS \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --beta beta --se se --pval P \
  --format cojo \
  --threads 4

The same style works for every workflow: add the flags for the operations you want, and GWAStoolkit will execute them in the correct order.

Workflow Steps

This section follows the real execution order used internally by GWAStoolkit.

1. OR to Beta/SE

./GWAStoolkit \
  --gwas-summary gwas_or.txt.gz \
  --out gwas.beta.txt.gz \
  --or OR \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --pval P --N N \
  --format gwas

If SE is unavailable, GWAStoolkit can estimate SE from P-value. If both SE and P are absent, SE is set to 999.

2. Add Sample Size Without Standardization

Fixed case/control counts:

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out gwas.N.txt.gz \
  --case 20000 \
  --control 30000 \
  --format gwas

Per-row case/control columns:

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out gwas.N.txt.gz \
  --case-col n_case \
  --control-col n_control \
  --format gwas

3. Standardized Effect

Binary trait with fixed case/control counts:

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out final.std.cojo.txt.gz \
  --case 20000 \
  --control 30000 \
  --std-effect \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --beta beta --se se --pval P \
  --format cojo

Continuous trait with fixed total N:

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out final.std.cojo.txt.gz \
  --N 50000 \
  --std-effect \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --beta beta --se se --pval P \
  --format cojo

Continuous trait with an existing total-N column:

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out final.std.cojo.txt.gz \
  --N N \
  --std-effect \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --beta beta --se se --pval P \
  --format cojo

4. rsID Annotation

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out gwas.rsid.txt.gz \
  --dbsnp GRCH37.dbSNP157.txt \
  --dbchr CHROM --dbpos POS --dbA1 REF --dbA2 ALT --dbrsid ID \
  --chr CHR --pos POS \
  --A1 A1 --A2 A2 \
  --format gwas

rsID matching uses:

GWAS CHR/POS == dbSNP CHR/POS
and
{GWAS A1, GWAS A2} == {dbSNP REF, dbSNP ALT}

Allele matching is case-insensitive and order-insensitive. No strand-complement matching is applied.

5. Final Format Conversion

./GWAStoolkit \
  --gwas-summary gwas.txt.gz \
  --out gwas.cojo.txt.gz \
  --SNP SNP --A1 A1 --A2 A2 \
  --freq freq --beta beta --se se --pval P --N N \
  --format cojo

Formulas

OR to Beta

$$ \beta = \ln(\mathrm{OR}) $$

If P-value is available and SE is absent, GWAStoolkit estimates SE from the two-sided P-value:

$$ z = \left|\Phi^{-1}\left(1 - \frac{P}{2}\right)\right| $$

$$ \mathrm{SE} = \frac{|\beta|}{z} $$

If both SE and P are absent, SE is set to 999 as a missing-SE sentinel.

Total Sample Size

Without --std-effect, case/control flags write total sample size:

$$ N_{\mathrm{total}} = N_{\mathrm{case}} + N_{\mathrm{control}} $$

Effective Sample Size for Binary Traits

With --std-effect, binary traits use effective sample size:

$$ N_{\mathrm{eff}} = \frac{4N_{\mathrm{case}}N_{\mathrm{control}}} {N_{\mathrm{case}} + N_{\mathrm{control}}} $$

For per-row case/control columns:

$$ N_{\mathrm{eff}, i} = \frac{4N_{\mathrm{case}, i}N_{\mathrm{control}, i}} {N_{\mathrm{case}, i} + N_{\mathrm{control}, i}} $$

Standardized Effect

GWAStoolkit standardizes beta/SE using:

$$ z = \frac{\beta}{\mathrm{SE}} $$

$$ \mathrm{SE}_{\mathrm{std}} = \frac{1} {\sqrt{2f(1-f)(N_{\mathrm{scale}} + z^2)}} $$

$$ \beta_{\mathrm{std}} = z \times \mathrm{SE}_{\mathrm{std}} $$

For binary traits:

$$ N_{\mathrm{scale}} = N_{\mathrm{eff}} $$

For continuous traits:

$$ N_{\mathrm{scale}} = N_{\mathrm{total}} $$

dbSNP Preparation

GWAStoolkit can read dbSNP as a tab-delimited text file with columns such as:

CHROM POS ID REF ALT

It also supports PLINK .bim / .bim.gz, interpreted as:

CHR RSID CM POS A1 A2

VCF to Text Example

Example for GRCh37 dbSNP 157:

  1. Split multi-allelic sites into biallelic
wget https://ftp.ncbi.nih.gov/snp/latest_release/VCF/GCF_000001405.25.gz

gunzip -c GCF_000001405.25.gz | bgzip -c > GCF_000001405.25.bgz
tabix -p vcf GCF_000001405.25.bgz

bcftools norm -m -any --thread 10 -Oz -o GRCH37.dbsnp157.vcf.gz GCF_000001405.25.bgz
  1. Rename chromosome names (optional, recommended)
cat > chrmap.txt <<EOF
NC_000001.10    1
NC_000002.11    2
NC_000003.11    3
NC_000004.11    4
NC_000005.9     5
NC_000006.11    6
NC_000007.13    7
NC_000008.10    8
NC_000009.11    9
NC_000010.10    10
NC_000011.9     11
NC_000012.11    12
NC_000013.10    13
NC_000014.8     14
NC_000015.9     15
NC_000016.9     16
NC_000017.10    17
NC_000018.9     18
NC_000019.9     19
NC_000020.10    20
NC_000021.8     21
NC_000022.10    22
NC_000023.10    X
NC_000024.9     Y
EOF

bcftools annotate --rename-chrs chrmap.txt GRCH37.dbsnp157.vcf.gz --thread 10 -Oz -o GRCH37.dbsnp157.chr.vcf.gz
bcftools index GRCH37.dbsnp157.chr.vcf.gz
  1. Extract chromosomes 1 to X, then export the needed columns:
bcftools view -r 1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,X \
  GRCH37.dbsnp157.vcf.gz \
  --thread 10 -Oz \
  -o GRCH37.dbsnp157.1toXchr.vcf.gz

(
  echo -e "CHROM\tPOS\tID\tREF\tALT"
  bcftools query -f '%CHROM\t%POS\t%ID\t%REF\t%ALT\n' GRCH37.dbsnp157.1toXchr.vcf.gz
) > GRCH37.dbSNP157.txt

Then use:

--dbsnp GRCH37.dbSNP157.txt \
--dbchr CHROM --dbpos POS --dbA1 REF --dbA2 ALT --dbrsid ID

Chromosome Name Handling

GWAStoolkit canonicalizes common chromosome names internally. Typical accepted forms include:

Input form Internal meaning
1, 01, chr1, CHR1, NC_000001.11 chromosome 1
X, chrX, 23, NC_000023.11 chromosome X
Y, chrY, 24, NC_000024.11 chromosome Y
M, MT, chrM, NC_012920.1 mitochondrial chromosome

Build consistency is still required: GRCh37 GWAS should be matched against GRCh37 dbSNP, and GRCh38 GWAS against GRCh38 dbSNP.

Arguments

Required

Flag Meaning
--gwas-summary FILE Input GWAS summary statistics. Text and gzipped input are supported.
--out FILE Output file. Text and gzipped output are supported.

Common Column Flags

Flag Meaning Default
--SNP COL SNP identifier column SNP
--chr COL --pos COL GWAS chromosome/position columns. Required only with --dbsnp. CHR / POS
--A1 COL --A2 COL GWAS allele columns A1 / A2
--freq COL Allele frequency column freq
--beta COL Beta column to read/write b
--se COL Standard error column to read/write se
--pval COL P-value column p
--N COL_OR_VALUE, --n COL_OR_VALUE Sample-size column name or fixed numeric total N N

Operation Flags

Flag Meaning
--or COL Convert odds ratio in COL to beta/SE.
--case INT --control INT Use fixed case/control counts. Without --std-effect, writes total N. With --std-effect, uses binary effective N.
--case-col COL --control-col COL Use per-row case/control columns.
--std-effect Standardize beta/SE using sample size.
--dbsnp FILE Annotate/fill rsID from dbSNP text/gz or PLINK .bim.
--dbchr COL --dbpos COL dbSNP chromosome/position columns.
--dbA1 COL --dbA2 COL dbSNP allele columns, usually REF and ALT.
--dbrsid COL dbSNP rsID column, usually ID.
--format gwas|cojo|popcorn|mrmega Final output format.

QC and Runtime Flags

Flag Meaning Default
--maf VALUE MAF threshold. Must be between 0 and 0.5. 0.01
--remove-dup-snp Keep one row per SNP, using the smallest P-value when available. off
--threads N, --thread-num N Number of OpenMP threads. 1
--log FILE Write log output to file. none
--verbose Print column mappings, scan progress, temporary files, and rsID mismatch examples. off
--keep-temp Keep intermediate files generated by multi-step calls. off

Format Requirements

For final non-gwas output, the processed GWAS table must contain the columns required by that format. These columns may already exist in the input or be generated by earlier flags in the same call.

Format Required columns
cojo SNP, A1, A2, freq, beta, se, p, N
popcorn SNP, A1, A2, freq, beta, se, N
mrmega SNP, A1, A2, freq, beta, se, p, N

If required fields are missing, add the relevant operation flags first. For example, use --or to create beta/SE from OR, and use --case --control, --case-col --control-col, or numeric --N to create N.

Example COJO output:

SNP       A1  A2  freq   b       se      p       N
rs1000    A   G   0.37   0.145   0.035   1e-5    50000
rs2000    T   C   0.42  -0.080   0.025   2e-3    50000

Logging

Use --verbose when debugging column names, dbSNP scans, or low rsID match rates:

./GWAStoolkit ... --verbose --log run.log

Verbose mode prints column mapping, temporary file paths, dbSNP progress, allele direction counts, and the first rsID allele mismatch examples.

Troubleshooting

Very Low rsID Match Rate

  • Check genome build consistency: GRCh37 vs GRCh38.
  • Confirm --chr, --pos, --A1, and --A2 point to the same columns used in your reference matching script.
  • Confirm dbSNP columns with --dbchr, --dbpos, --dbA1, --dbA2, and --dbrsid.
  • Use --verbose and inspect rsID allele mismatch example lines.
  • Make sure multi-allelic sites were split with bcftools norm -m -any.

COJO Output Missing N

COJO requires an N column. Use one of:

--N N
--N 50000
--case 20000 --control 30000
--case-col n_case --control-col n_control

Continuous Trait Standardization

For continuous traits, use total N through --N. Do not use case/control unless the trait is binary.

Binary Trait Standardization

For binary traits, use case/control counts and --std-effect so GWAStoolkit can compute effective N.

About

A unified, high-performance C++ toolkit for processing GWAS summary statistics.

Resources

Stars

Watchers

Forks

Releases

Packages

Contributors

Languages