​​​​R语言学习笔记(五)——曼哈顿图

iJournal

学术期刊信息查询
386篇原创内容
Official Account

 ↑ ↑  关注iJournal,选刊快人一步  ↑ ↑ 

iJournal后台回复“2021学科”即可获得20种

不同学科的最新影响因子文件(Excel版本)

导 语

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

01

GWAS简介

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

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

02

qqman

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

> install.packages('qqman')> library(qqman)           #加载qqman包> library(RColorBrewer)    #用于颜色变化> str(gwasResults)> head(gwasResults)SNP CHR BP         P1 rs1   1  1 0.91480602 rs2   1  2 0.93707543 rs3   1  3 0.28613954 rs4   1  4 0.83044765 rs5   1  5 0.64174556 rs6   1  6 0.5190959

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

> manhattan(gwasResults)

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

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

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

这种翘尾巴的形式也是最理想的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)

03

CMplot

这个包实际上只有一个函数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 trait31 ALGA0000009 1 52297 0.7738187 0.51194318 0.511943182 ALGA0000014 1 79763 0.7738187 0.51194318 0.511943183 ALGA0000021 1 209568 0.7583016 0.98405289 0.984052894 ALGA0000022 1 292758 0.7200305 0.48887140 0.488871405 ALGA0000046 1 747831 0.9736840 0.22096836 0.220968366 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)

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)

上面的图片是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)

可以看到由于不同的信号叠加到一起显得非常拥挤,所以常规的曼哈顿图不适合展示多性状的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结果。

3.4 qq图

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

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

04

小 结

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

(0)

相关推荐

  • 全了!!曼哈顿图样样式、方法大汇总(Python R)~

    最近小编在后台看到有的小伙伴留言咨询曼哈顿图(Manhattan Plot) 的绘制方法,小编一开始也是比较不了解,奈何我又是一个宠读者的小编,这就汇总了曼哈顿图(Manhattan Plot) R和 ...

  • R语言学习笔记

    R语言学习笔记

  • R语言学习笔记(四):apply,sapply,lapply,tapply,vapply以及mapply的用法

    apply() apply(m,dimcode,f,fargs) m 是一个矩阵. dimcode是维度编号,取1则为对行应用函数,取2则为对列运用函数. f是函数 fargs是f的可选参数集 > ...

  • 那些年倒腾的R语言学习笔记,全都在这里了~

    今天这一篇整理一下我以往推送过的所有R语言相关文章,一来是方便大家的检索,二来也是阶段性学习的一次总结. 关于内容分类,我分成了学习心得篇.R语言基础.数据可视化.网络数据爬取,然后各自类别进行详细的 ...

  • R语言学习系列之“多变的热图”

    咱公众号也不能只做一个系列,所以经过深思熟虑,打算将来慢慢增加一些内容,主要有以下几个系列TCGA数据分析系列GEO数据分析系列"老板给一个基因,我该怎么办"系列文献阅读系列R语言 ...

  • 跟着Nature Genetics 学画图:R语言ggplot2画箱线图(boxplot)展示D s...

    简介:R语言统计与绘图公众号目前致力于分享医学统计与R绘图知识,手把手教你使用R语言绘制基线特征表.KM生存曲线.森林图.ROC曲线等.每天一篇精彩R语言推文教程,手把手带你入门R语言绘图. 今天推文 ...

  • 手把手教你使用R语言绘制交互效应的森林图

    交互效应分析在SCI中是一个重要的组成,可以表示特殊的人群(也就是亚组)中自变量和因变量的结果关系,可以说是高分SCI必做的项目,在SCI中,这一部分结果很多都是以森林图来展示,这样结果又直观有方便分 ...

  • C# LINQ学习笔记五:LINQ to XML

    本笔记摘抄自:https://www.cnblogs.com/yaozhenfa/p/CSharp_Linq_For_Xml.html,记录一下学习过程以备后续查用. 一.生成xml 1.1创建简单的 ...

  • 零基础R语言学习路线

    其实相对于常见的编程语言,R语言还是非常容易上手.学习1年多时间,就可以找一份不错的工作了. 前言 我当初学习R的时候在网上搜到一则流传很广的R语言学习路线图(R语言学习由浅入深路线图),我在微信圈, ...

  • 技术贴 | R语言pheatmap聚类分析和热图

    本文由阿童木根据实践经验而整理,希望对大家有帮助. 原创微文,欢迎转发转载. 导读 pheatmap默认会对输入矩阵数据的行和列同时进行聚类,但是也可以通过布尔型参数cluster_rows和clus ...