当我们拿到某个研究的
GWAS summary data时，因为大多是二代测序的缘故，它的SNP ID可能不是正常的rsid，但是，我们下游有些分析却要用到rsid，这时我们会想到去dbSNP库去下载对应版本参考基因组的rsID。
当下载完毕，并且整理成只剩下SNP CHR POS这三列时，又发现文件太大了，如果我们没有基因型数据支持我们利用PLINK软件注释SNP ID的话，就会很难办，因为用R或Python读取大型文件又是个问题（
GRCh37大概15G），这时候，当然是分染色体去进行注释咯。
第一步
使用awk对SNP CHR POS三列文件的GRCh37 参考SNP ID分染色体保存
awk '{print > $1".txt"}' dbSNP_GRCh37
mkdir GRCh37 && mv *.txt GRCh37
第二步
要注释的文件需包括CHR, POS, SNP三列
library(data.table)
setwd("~/new_run")
file = "reformatMETAL.gz"
df = fread(file)[, c(CHR, POS, SNP, A1, A2, freq, beta, SE, p, N)]
out = data.table()
for (chr in c(1:22, "X")){
setwd("~/dbSNP/GRCh37")
db = fread(paste0(chr, ".txt"), header=FALSE, col.names=c("CHR", "POS", "SNP"))
df_merge = merge(df, db, by=c("CHR", "POS"), all.x=TRUE, suffix=c("df", "db"))
df_merge$SNP = ifelse(is.na(df_merge$SNP.db), df_merge$SNP.df, df_merge$SNP.db)
out = rbind(out, df_merge)
}
fwrite(out, file="my_test.gz", row.names=FALSE, quote=FALSE, sep="\t")
就是这么简单，散会！！！
