rMVP包 (GWAS)

2023-08-03  本文已影响0人  Felix_xpz

1. 处理数据

####初次质控
plink --bfile mvp --maf 0.05 --geno 0.1 --mind 0.1 --recode vcf --out mvp_qc1
####Beagle填充
java -Xmx2g -jar beagle.28Sep18.793.jar nthreads=4 gt=mvp_qc1.vcf out=imput.mvp.geno 
####对填充后文件进行再质控
plink --vcf imput.mvp.geno.vcf --maf 0.05 --geno 0.1 --mind 0.1 --make-bed --out mvp_qc2

2. 运行R包rMVP

#library("devtools") 
#devtools::install_github("xiaolei-lab/rMVP")
install.packages("rMVP")
library(rMVP)

# Step1: 数据格式转换
MVP.Data(filePhe = "mvp.phe", fileBed = "mvp_qc2", out = "demo")

# Step2: 读取转换好的数据
pheno <- read.table("demo.phe", header = TRUE)       # 表型文件中第一列为个体 ID,第二列是要分析的性状
geno  <- attach.big.matrix("demo.geno.desc")         # rMVP 使用的基因型文件是 .geno.desc 和 .geno.bin 两个文件一组,前者储存的是元数据,后者是以二进制格式储存的基因型数据。
map   <- read.table("demo.geno.map", header = TRUE)  # map 文件中存储的是标记的位置信息

# Step3: 进行分析
imMVP <- MVP(
    phe=pheno,    #输入表型文件
    geno=geno,    #输入基因型文件
    map=map,      #输入map文件
    nPC.GLM=5,    #该参数设置的是GLM模型中加入的PC数量,此处设置为前5个PC
    nPC.MLM=3,    #该参数设置的是MLM模型中加入的PC数量,此处设置为前3个PC
    nPC.FarmCPU=3,#该参数设置的是FarmCPU模型中加入的PC数量,此处设置为前3个PC
    threshold=0.05,#该参数设置的是显著水平线的阈值,采用bonfferoni校正,threshold = 0.05/n, n为marker数量
    method=c("GLM", "MLM", "FarmCPU") # 模型可以设置多个同时进行分析。
)

# 多表型情况下,可如下运行
for (i in 2:ncol(pheno)){
    imMVP <- MVP(
        phe=pheno,    #输入表型文件
        geno=geno,    #输入基因型文件
        map=map,      #输入map文件
        nPC.GLM=5,    #该参数设置的是GLM模型中加入的PC数量,此处设置为前5个PC
        nPC.MLM=3,    #该参数设置的是MLM模型中加入的PC数量,此处设置为前3个PC
        nPC.FarmCPU=3,#该参数设置的是FarmCPU模型中加入的PC数量,此处设置为前3个PC
        threshold=0.05,#该参数设置的是显著水平线的阈值,采用bonfferoni校正,threshold = 0.05/n, n为marker数量
        method=c("GLM", "MLM", "FarmCPU") # 模型可以设置多个同时进行分析。
    )
}
上一篇下一篇

猜你喜欢

热点阅读