本篇推文的示例数据来源于参考书 《Genome-Wide Association Studies》的第十章 A Practical Guide to Using Structural Variants for Genome-Wide Association Studies。
植物里做GWAS分析通常是选择某个群体做二代基因组测序（有的已经研究比较多的物种比如 水稻、玉米可以能已经发表过很多数据，），测序数据与参考基因组进行比对鉴定变异位点，然后用变异位点和表型特征去做关联分析。变异位点包括snp indel SV SRR STR TIP 等等。表型数据通常常规的植物表型，比如株高、产量、花瓣数量等。也可以是分子层面的代谢物，基因表达量等。也可以是环境数据，鉴定哪些变异位点和环境有关。
变异位点的数据通常是用vcf文件存储（当然也有其他格式）。我这篇推文介绍用vcf文件去做GWAS。
首先介绍一下vcf文件的格式
vcf文件是文本文件，我们在自己电脑上直接用记事本打开就可以查看文件里的内容。通常样本很多的话，对应的文件也会非常大。自己电脑打开还挺费劲的。
vcf 文本里的内容按照特定的模式排列
vcf简单的可以划分为三个部分
1、两个#号开头的行
2、一个#号开头的行
3、零个#号开头的行
表型数据
两列，第一列是样本名字，第二列是表型的值
如果是用rMVP这个R包来做GWAS的话表型数据的样本顺序和vcf文件的样本顺序不一致也可以，但是其他软件有的会要求样本顺序一致
rMVP 这个R包的github主页
https://github.com/xiaolei-lab/rMVP
只使用混合线性模型的代码
代码语言：javascript复制setwd("D:/Jupyter/practice/rMVP_GWAS/")
library(rMVP)
MVP.Data(fileVCF="smoove_filtered.vcf",
filePhe="phenotype.txt",
sep.phe = "t",
fileKin=TRUE,
filePC=TRUE,
out="smoove"
)
genotype <- attach.big.matrix("smoove.geno.desc")
phenotype <- read.table("smoove.phe",head=TRUE)
map <- read.table("smoove.geno.map" , head = TRUE)
Kinship <- attach.big.matrix("smoove.kin.desc")
Covariates_PC <- bigmemory::as.matrix(attach.big.matrix("smoove.pc.desc"))
imMVP <- MVP(
phe=phenotype,
geno=genotype,
map=map,
K=Kinship,
CV.MLM=Covariates_PC,
nPC.MLM=10,
priority="speed", ##for Kinship construction
ncpus=2,
vc.method="BRENT", ##only works for MLM
maxLoop=10,
threshold=0.05,
method=c("MLM"),
file.output=c("pmap", "pmap.signal", "plot", "log")
)
