1
$\begingroup$

There are some questions about lifting between reference builds, e.g. this one. But there doesn't appear to be a question about lifting a GWAS results file to a new reference build (except off-site). GWAS results files are not VCFs or plink files, so they aren't covered by the usual lifting tools (such as CrossMap or LiftOver). How do I lift a GWAS results file to hg38?

For instance, the Lee et al 2018 GWAS on education has this file:

user@computer:/home/user/Lee_2018$ head GWAS_EA.to10K.txt
MarkerName  CHR POS A1  A2  EAF Beta    SE  Pval
rs9859556   3   49455986    G   T   0.6905  -0.02901    0.00151 4.61e-82
rs7623659   3   49414791    T   C   0.3095  0.02899 0.00151 6.05e-82
rs11917431  3   49644012    C   T   0.6973  -0.0292 0.00152 7.40e-82
rs1873625   3   49666964    C   A   0.6973  -0.02916    0.00152 8.10e-82
rs11921590  3   49644193    T   C   0.6973  -0.02916    0.00152 1.17e-81
rs2352974   3   49890613    T   C   0.5 -0.02569    0.00141 6.35e-74
rs7029718   9   23358495    G   A   0.5646  -0.0244 0.00143 7.82e-65
rs7868984   9   23357826    C   T   0.4388  0.02437 0.00143 8.21e-65
rs11793831  9   23362311    T   G   0.4354  0.02418 0.00143 9.34e-64
$\endgroup$

3 Answers 3

2
$\begingroup$

Yes. LiftOver tools, if they work on VCF data (some only take BED as input) work on GWAS summary statistics. You could use https://pypi.org/project/liftover/, for example.

$\endgroup$
1
  • $\begingroup$ Yeah, I figure out the simpler solution is to make a small VCF with an empty case, and then lift that to hg38 using liftover, then merge back the column of new positions. $\endgroup$ Jun 15, 2023 at 23:27
1
$\begingroup$

It turned out to be easier than expected. You figure out which dbSNP version they used. dbSNP is the classification system for snps based on the rsid format. Pretty much every GWAS publishes their results with the rsids for the variants, as done above as well (MarkerName). So you download one of the dbSNP files from the official site. In this case, I am using version 142 since this was used by the 1000 genomes dataset.

After this, it is a matter of joining the tables. This can be done easily in R (you need about 15GB free memory), we load the dbSNP database into memory. Since we don't need all the columns, we skip some of them to save memory:

#library
library(tidyverse)

#cols
vcf_cols = cols_only(
  `#CHROM` = col_character(),
  `POS` = col_number(),
  `ID` = col_character(),
  `REF` = col_character(),
  `ALT` = col_character()
  )

#this is from the 1kg
#ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/technical/reference/GRCh38_reference_genome/other_mapping_resources/
dbsnp142 = readr::read_tsv(
  "ALL_20141222.dbSNP142_human_GRCh38.snps.vcf"),
  comment = "##",
  col_types = vcf_cols
  )

GWAS_EA = read_tsv(GWAS_EA.to10K.txt,
                   col_types = cols(chrpos37 = col_character(),
                                    chrpos38 = col_character()))

#make chrpos
dbsnp142 %<>% mutate(
  chrpos38 = str_c((CHROM %>% str_replace("chr", "")), ":", POS)
)

#count overlap
(GWAS_EA$MarkerName %in% dbsnp142$ID) %>% table2()

#mutate and merge
GWAS_EA = GWAS_EA %>%
  mutate(
    chrpos37 = str_c(CHR, ":", POS)
  ) %>%
  left_join(dbsnp142 %>% select(ID, chrpos38), by = c("MarkerName" = "ID"))

#update file on disk
GWAS_EA %>% write_tsv("GWAS_EA.to10K.txt")

The two new columns, chrpos37 and chrpos38, have the chr:pos variant names in both versions.

$\endgroup$
1
$\begingroup$

I would advise against using dbSNP rsIDs to perform a liftover as they get updated over time. I recently wrote a software to liftover summary statistics and the documentation includes instructions that show you how to easily generate educational attainment summary statistics for GRCh38

$\endgroup$
2
  • $\begingroup$ Your answer could be improved with additional supporting information. Please edit to add further details, such as citations or documentation, so that others can confirm that your answer is correct. You can find more information on how to write good answers in the help center. $\endgroup$
    – Community Bot
    Jan 11, 2023 at 20:27
  • $\begingroup$ Once the plugins have been correctly installed, you could lift your GWAS results with a command like this: bcftools +munge -f human_g1k_v37.fasta -C colheaders.tsv GWAS_EA.to10K.txt | bcftools +liftover -- -s human_g1k_v37.fasta -f Homo_sapiens_assembly38.fasta -c hg19ToHg38.over.chain.gz $\endgroup$ Jan 12, 2023 at 21:59

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.