快捷搜索:  汽车  科技

r语言怎么做图像分析()

r语言怎么做图像分析()最后看一下曼哈顿图的相关参数> qq(gwasResults$P) #qq图只需要一列p值的数据这种翘尾巴的形式也是最理想的qq图结果,可以看到从横坐标大于2开始,GWAS结果的p值与均匀分布的p值就有了明显的差距,说明表型和基因型之间确实存在显著的相关关系。> install.packages("qqman") > library(qqman) #加载qqman包 > library(RColorBrewer) #用于颜色变化 > str(gwasResults) > head(gwasResults) SNP CHR BP P 1 rs1 1 1 0.9148060 2 rs2 1 2 0.9370754 3 rs3 1 3 0.2861395 4 rs4 1 4

导 语

全基因组关联分析(Genome-wide association study)是生物领域挖掘功能基因的常用方法。GWAS的最重要的结果展示就是曼哈顿图,本期给大家介绍两个画曼哈顿图的R包qqman和CMplot,并对二者做一个简单小结。

01GWAS简介

全基因组关联分析最早是在人类疾病研究中被应用,随后在动物和植物研究中也大放异彩。从2005年第一篇GWAS研究开始到现在已经16年了,很多重要的基因都已经被挖掘出来了,但值得注意的是每年仍然有不少高水平的GWAS相关文章。尤其是近年来转录组、蛋白组、代谢组等多组学的兴起,可能会大大扩展GWAS的研究边界,让GWAS焕发出第二春。

关于GWAS具体的分析流程网上的资料很多,在这里不做更多的介绍。这里只介绍GWAS的结果展示方式,即曼哈顿图和qq图。其中曼哈顿图显示所有SNP位点的p-value,可以理解为每个SNP与表型的关联程度;qq图的纵轴是SNP位点的p-value值,横轴是则是均匀分布的p-value值,显示了每个SNP位点p-value实际值与理论值(假设SNP与表型不相关)的差异。

02qqman

这个包的用法比较简单,这里用包中自带的示例文件gwasResults展示

> install.packages("qqman") > library(qqman) #加载qqman包 > library(RColorBrewer) #用于颜色变化 > str(gwasResults) > head(gwasResults) SNP CHR BP P 1 rs1 1 1 0.9148060 2 rs2 1 2 0.9370754 3 rs3 1 3 0.2861395 4 rs4 1 4 0.8304476 5 rs5 1 5 0.6417455 6 rs6 1 6 0.5190959

可以看到数据一共4列,分别是SNP标记的名称、染色体、染色体位置、P值。我们先不调参数画个最简单的曼哈顿图

> manhattan(gwasResults)

r语言怎么做图像分析()(1)

从曼哈顿图可以看到只有3号染色体上有一个很高的峰,这几乎是最理想的GWAS结果。

同样,我们画出相应的qq图:

> qq(gwasResults$P) #qq图只需要一列p值的数据

r语言怎么做图像分析()(2)

这种翘尾巴的形式也是最理想的qq图结果,可以看到从横坐标大于2开始,GWAS结果的p值与均匀分布的p值就有了明显的差距,说明表型和基因型之间确实存在显著的相关关系。

最后看一下曼哈顿图的相关参数

manhattan(x chr = "CHR" bp = "BP" p = "P" snp = "SNP" col = c("gray10" "gray60") chrlabs = NULL suggestiveline = -log10(1e-05) genomewideline = -log10(5e-08) highlight = NULL logp = TRUE annotatePval = NULL annotateTop = TRUE ...)

这里只介绍几个常用的参数,col是调整颜色,suggestiveline是p的阈值,genomewideline是第二条阈值线,highlight是标记某个或某些SNP。

> manhattan(gwasResults col = c("red" "blue") annotatePval = 0.0001)

r语言怎么做图像分析()(3)

03CMplot

这个包实际上只有一个函数CMplot,但集成了SNP密度图、曼哈顿图和qq图多种图形的画法。这个函数的参数有很多,这里列举一些常用的参数。

CMplot(Pmap col=c("#377EB8" "#4DAF4A" "#984EA3" "#FF7F00") bin.size=1e6 bin.max=NULL pch=19 band=1 cir.band=0.5 H=1.5 ylim=NULL cex.axis=1 plot.type="b" multracks=FALSE cex=c(0.5 1 1) r=0.3 xlab="Chromosome" ylab=expression(-log[10](italic(p))) xaxs="i" yaxs="r" outward=FALSE threshold = NULL threshold.col="red" threshold.lwd=1 threshold.lty=2 amplify= TRUE chr.labels=NULL signal.cex = 1.5 signal.pch = 19 signal.col="red" signal.line=1 cir.chr=TRUE cir.chr.h=1.5 chr.den.col=c("darkgreen" "yellow" "red") cir.legend=TRUE cir.legend.cex=0.6 cir.legend.col="black" LOG10=TRUE box=FALSE conf.int.col="grey" file.output=TRUE file="jpg" dpi=300 memo="")

col 设置颜色cex/pch 设置点的大小/形状

bin.size 设置SNP密度图中的窗口大小

cex.axis 设置坐标轴字体和标签字体的大小

plot.type 设置不同的绘图类型,可以设定为 "d" "c" "m" "q" or "b",其中d是SNP密度图,c是环形曼哈顿图,m是曼哈顿图,q是qq图,b是同时画环形曼哈顿图、曼哈顿图和qq图。

threshold/ threshold.col/ threshold.lwd/ threshold.lty 设置阈值并添加阈值线/阈值线的颜色/宽度/类型

signal.cex/signal.pch/signal.col 设置显著点的大小/性状/颜色

cir.legend/cir.legend.cex/cir.legend.col 设置是否显示图例/图例字体大小/图例颜色

3.1 SNP密度图

> #install.packages(CMplot) > library(CMplot) > data = pig60K > head(data) SNP Chromosome Position trait1 trait2 trait3 1 ALGA0000009 1 52297 0.7738187 0.51194318 0.51194318 2 ALGA0000014 1 79763 0.7738187 0.51194318 0.51194318 3 ALGA0000021 1 209568 0.7583016 0.98405289 0.98405289 4 ALGA0000022 1 292758 0.7200305 0.48887140 0.48887140 5 ALGA0000046 1 747831 0.9736840 0.22096836 0.22096836 6 ALGA0000047 1 761957 0.9174565 0.05753712 0.05753712

可以看到示例数据pig60K有6列,前3列是SNP信息,后三列是表型数据。我们自己用CMplot包作图时,可直接保留列名,将数据替换成自己的数据即可。

> CMplot(pig60K plot.type = "d" bin.size = 1e5 col = c("blue" "red" "yellow") file.output = F)

r语言怎么做图像分析()(4)

3.2 曼哈顿图

单性状曼哈顿图

> CMplot(pig60K plot.type = "m" threshold = c(0.01 0.05)/nrow(pig60K) amplify = T signal.cex = c(1 1) signal.pch = c(20 20) signal.col = c("red" "blue") multracks = F file.output = F)

r语言怎么做图像分析()(5)

上面的图片是trait1的曼哈顿图,若绘制多性状曼哈顿图则是这样的:

> CMplot(pig60K plot.type = "m" threshold = c(0.01 0.05)/nrow(pig60K) amplify = T signal.cex = c(1 1) signal.pch = c(20 20) signal.col = c("red" "blue") multracks = F file.output = F)

r语言怎么做图像分析()(6)

可以看到由于不同的信号叠加到一起显得非常拥挤,所以常规的曼哈顿图不适合展示多性状的GWAS结果。那用什么展示呢,环形曼哈顿图!

环形曼哈顿图

> CMplot(pig60K plot.type="c" r=0.5 threshold=c(1e6 1e6) cex = 1 threshold.col = c("red" "blue") cir.chr.h = 2 signal.cex = c(2 2) signal.col=c("red" "green") file.output = F)

可以看到环形的曼哈顿图能够同时显示3个性状的GWAS结果。

r语言怎么做图像分析()(7)

3.4 qq图

qq图的做法比较简单,一般只需要把作图类型改成q,设置阈值,调整一下字体大小即可。

> CMplot(pig60K plot.type = "q" threshold = 0.05)

r语言怎么做图像分析()(8)

04小 结

本期给大家介绍了用于展示GWAS结果曼哈顿图和qq图的R包qqman和CMplot,其中qqman的用法比较简单,功能也相对单一。而CMplot除了曼哈顿图和qq图,还能够画SNP密度图,并且能够用环形曼哈顿图同时展示多个GWAS结果,推荐大家使用!

猜您喜欢: