GWAS分析 一般线性模型GLM中meta、se、T、p值的计算
作者：互联网
001、plink
root@PC1:/home/test# ls gwas_test.map gwas_test.ped root@PC1:/home/test# plink --file gwas_test --assoc 1> /dev/null root@PC1:/home/test# ls gwas_test.map gwas_test.ped plink.log plink.qassoc root@PC1:/home/test# head plink.qassoc CHR SNP BP NMISS BETA SE R2 T P 1 snp1 2802 541 -8.911 8.344 0.002111 -1.068 0.286 1 snp2 2823 541 7.754 10.04 0.001104 0.772 0.4405 1 snp3 4512 541 9.264 9.814 0.00165 0.9439 0.3456 1 snp4 16529 541 -18.49 10.81 0.005401 -1.711 0.08769 1 snp5 16578 541 5.661 9.93 0.0006026 0.5701 0.5689 1 snp6 16579 541 -5.577 5.657 0.0018 -0.9858 0.3247 1 snp7 16635 541 2.985 8.717 0.0002176 0.3425 0.7321 1 snp8 20879 541 9.053 7.29 0.002853 1.242 0.2148 1 snp9 20908 541 9.278 6.695 0.00355 1.386 0.1664
002、R计算
root@PC1:/home/test# ls gwas_test.map gwas_test.ped root@PC1:/home/test# plink --file gwas_test --recode A 1> /dev/null root@PC1:/home/test# ls gwas_test.map gwas_test.ped plink.log plink.raw root@PC1:/home/test# cut -d " " -f 1,3-5 --complement plink.raw > num.raw root@PC1:/home/test# ls gwas_test.map gwas_test.ped num.raw plink.log plink.raw
library(data.table)
dat <- fread("num.raw", data.table = F)
dat[1:3,1:8]
result <- data.frame()
for (i in 3:12) {
reg = lm(PHENOTYPE ~ 1 + dat[,i], data=dat)
result <- rbind(result, summary(reg)$coefficients[2,])
}
names(result) <- c("beta","se","t","p")
result
标签：541,GWAS,GLM,PC1,gwas,meta,test,home,plink 来源： https://www.cnblogs.com/liujiaxin2018/p/16523231.html