16s分析之进化树+差异分析(二)

上回我们说到,R语言16s进化树的构建,单单一个优秀的进化树对于我们不是做进化的人来讲,还是不够,有些华而不实,但是Y叔这个ggtree功能当然不止如此,Y叔也经常宣称其为真正的’支持图形语法’;相比功能和可塑性都会提升,就下来我们就进一步深入了解ggtree

这里我要做的就是将进化树连同OTU丰度信息展示出来,ggtree可以做,我们就来一一展示,套用上次推送的代码做前面处理:(如果没有看,请点击:16s分析之进化树+差异分析(一)


p<-ggtree(tree3,layout="fan", branch.length = "none",ladderize = FALSE,aes(color=group)) %<+% tax1+

geom_tiplab2(size=1.5)+ theme(legend.position= "right")

#查看节点;这里简单说几句,geom_label2是Y树参照ggplot的geom_label在ggtree中重新写的,主要特点就是可以使用subset,直接进行数据取舍,这里的意思就是不显示进化树叶端点

p+geom_label2(aes(subset=!isTip,label=node),size=2)

##加上阴影,这里我们可以根据自己的需求选择合适的端点来添加高亮阴影,就16s数据,我挑选了这么几个大类

p1<-p+

geom_hilight(node=323, fill="gray",alpha=0.5) +

geom_hilight(node=320, fill="pink",alpha=0.5) +

geom_hilight(node=313,fill="beige", alpha=0.5) +

geom_hilight(node=298,fill="yellow", alpha=0.5)+

geom_hilight(node=255, fill="blue",alpha=0.5) +

geom_hilight(node=237, fill="red",alpha=0.5) +

geom_hilight(node=194,fill="green", alpha=0.5)+

geom_hilight(node=269,fill="lightblue", alpha=0.5)

p1#这就是我们上次得到的图形,这里我们再次展示一下:

#这里我仅仅添加了OTU名称信息


添加数据,今天的重头戏在于添加数据,表示差异,这里表示:我接下来的代码参考了宏基因组的推送,进化树,上期我已经转发过来了,大家可以自行查看(抄作业从不脸红)

###下面我们对进化树添加数据

## 读取OTU表

otu_table =read.delim("otu_table.txt", row.names= 1,  header=T, sep="\t")

## 分组文件我这里只选择我这次做的三个组

design =read.table("map_lxdjhg_CK.txt", header=T, row.names= 1,sep="\t")

## 这里我进行筛选

idx=intersect(rownames(design),colnames(otu_table))

sub_design=design[idx,]

## 计的样品顺序重排列

otu_table=otu_table[,idx]

## 将OTU表count转换为百分比

norm =t(t(otu_table)/colSums(otu_table,na=T)) * 100

## 筛选树中OTU对应的数据

tax_tomato =norm[rownames(tax),]

head(tax_tomato)

# 树+ 组均值热图

## 提取实验设计中的分组信息

sampFile =as.data.frame(sub_design$SampleType,row.names = row.names(sub_design))

colnames(sampFile)[1]= "group"

## OTU表转置,让样品名为行

mat_t = t(tax_tomato)

## 合并分组信息至丰度矩阵,并去除样品名列

mat_t2 =merge(sampFile, mat_t, by="row.names")[,-1]

## 按组求均值

mat_mean =aggregate(mat_t2[,-1], by=mat_t2[1], FUN=mean)

## 去除非数据列并转置

mean_tomato =do.call(rbind, mat_mean)[-1,]

## 重命名列名为组名

colnames(mean_tomato)= mat_mean$group

wt2<-sqrt(mean_tomato)#

p4 <-  gheatmap(p3, wt2, low="blue",high="red",offset = 15, width=0.2, font.size=1, hjust=-.1)

plot(p4)

#这个颜色似乎很重口味,换一个呗

p( p1, wt2,offset = 15,width=0.2, low="white", high="black",colnames_position = "top", font.size=2)

# plot

plot(p5)#这个颜色似乎就清爽许多

当然这还不够,达不到我们的要求

###这里调整热图和标签位置

p1<-ggtree(tree3,layout="fan", branch.length = "none",ladderize = FALSE,aes(color=group)) %<+% tax1+

geom_tiplab2(size=1.5, offset = 5)+theme(legend.position = "right")+

geom_hilight(node=323, fill="gray",alpha=0.5) +

geom_hilight(node=320, fill="pink",alpha=0.5) +

geom_hilight(node=313,fill="beige", alpha=0.5) +

geom_hilight(node=298,fill="yellow", alpha=0.5)+

geom_hilight(node=255, fill="blue",alpha=0.5) +

geom_hilight(node=237, fill="red",alpha=0.5) +

geom_hilight(node=194,fill="green", alpha=0.5)+

geom_hilight(node=269,fill="lightblue", alpha=0.5)

p1

p5 <-  gheatmap( p1, wt2, offset = 0,width=0.2,low="white", high="black", colnames_position ="top", font.size=1)

# plot

plot(p5)

##下面我们来对热图添加分类信息

p5+geom_cladelabel(node=323,"Galagoidea", offset=10,angle = 350, fontsize=3,  barsize=1,offset.text=1.5,color="gray") +

geom_cladelabel(320, "Lemuroidea",offset=9,angle = 10, fontsize=3, barsize=1, offset.text=1.5,color="pink") +

geom_cladelabel(313, "Tarsioidea",offset=11,angle = 50, fontsize=3, barsize=1, offset.text=1.5,color="beige") +

geom_cladelabel(298,"Ceboidea",offset=9.5,angle = 80, fontsize=3,  barsize=1,offset.text=1.5,color="yellow") +

geom_cladelabel(255, "Hominoidea",offset=9.5,angle = 305, fontsize=3, barsize=1, offset.text=1.5,color="blue") +

geom_cladelabel(237,"Cercopithecoidea",offset=10.5,angle = 360, fontsize=3,  barsize=1,offset.text=1.5,color="red") +

geom_cladelabel(194,"Cercopithecoidea",offset=11.5,angle = 70, fontsize=3,  barsize=1,offset.text=1.5,color="green")+

geom_cladelabel(269,"Cercopithecoidea",offset=9,angle = 320, fontsize=3,  barsize=1,offset.text=1.5,color="lightblue")

这里角度我没有调整好,所以自己调整,参数如下:

angle:就是设置分类label的角度

offset:距离中心的距离

fontsize:字体大小

barsize:标签柱子宽度

到这里,我们的圈图就到这里了,似乎在圈图在ggtree中只可以添加热图,柱状图,饼图目前也没有看到添加方法;

这里提到Y树ggtree的facet_plot,可以添加各种图形同树图一块展示,但是在layout为fan情况下,要想展示,就必须将其坐标转化为极坐标,而且需要相应调整,因为这种树图和数据结合采用圈图方式的确是逼格满满:

这里展示iTOL的主页:的确是漂亮

于是我问Y树,可不可以做成这种形式,答案在此:

这张代表我的心情:

我的耳旁响起:一千个伤心的理由,一千个伤心的理由,

怎么办呢,等Y叔的扩展,还是自己写一个,我写的了吗?不试试谁知道呢?


学习永无止境,分享永不停歇!

(0)

相关推荐

  • R语言GEO数据处理(七)

    # 6. 可视化展示 ---------------------------------------------------------------- ##6.1 火山图 library(ggplot ...

  • 16s分析之进化树+差异分析(一)

    之前我便说过,16s出图主要是用R,自这两年来Y叔的ggtree出来,好评不断,引用不断攀升 这里当然借用Y树的包来做我们的进化树,当然Y叔公众号讲的更加深入和详细,有需要请关注:biobabble ...

  • 16s分析之不同分类水平差异分析及气泡图绘制

    对otu的差异分析并不是我们唯一的选择,差异往大的做,可以往往门,纲,目,科做. 今天要做一张长的图,我们可以和别的图一起配合使用会好? 比如这篇文章,还是挺好看的: 下面是一份完整的代码,我仅仅只做 ...

  • 16s分析之差异分析(DESeq2)

    今天我们来学习R语言DESeq2包,原理什么的后不说,在操作过程中点缀一下,等四个差异包推送完成后,咱们在对这四个包做差异分析的原理做一个比较分析: 这个包适用于: 高通量数据分析过程中,基于coun ...

  • 16s分析之LEfSe分析学习笔记(数据整理的二三事)

    之前做 LEfSe 分析,使用全部的六个等级分类数据,总观测达到 1500,做出来树图很不好看,因此这次做一下调整: 1.  选择相对丰度在 0.001 以上的观测: 2.  将每个等级没有注释上名称 ...

  • 语文阅读分析训练(十二)

    分析人物形象及作用: (2019三州联考)阅读下面文字,回答1-3题.(12分) 火车侠 火车侠是我见过的最特别的一个人. 我记得是在高中的某节课上,我正被讲台上的老师念叨得昏昏欲睡,旁边一个女生悄悄 ...

  • 行为分析技术(十二)不同调整级别的确定

    2008/10/08 趋势在什么时间会出现多大级别的调整,相信这是人人关心的问题.有时候我们会发现,有些调整行情的级别非常小,可能几天就调整完毕.而有些时候却需要数周甚至数月的调整. 那么,趋势开始出 ...

  • 微生物组-扩增子16S分析和可视化(线上/线下同时开课,2021.7)

    福利公告:为了响应学员的学习需求,经过易生信培训团队的讨论筹备,现决定安排扩增子16S分析.宏基因组.转录组的线上/线下同时开课.报名参加线上直播课的老师可在1年内选择参加同课程的一次线下课 .期待和 ...

  • 八字分析 命中注定要娶二婚男的八字女

    几乎每个女人都对未来另一半充满幻想,希望是个称心如意的白马王子,然而命运是很残酷的,缘分到来的时候,很多事情都会身不由己,由不得自己选择,所以我们可以来看看传说中注定要嫁个二婚男的八字是什么样的,看看 ...

  • Microbiome:16s分析的价格,宏基因分析的结果,这个你可以有!

    俄亥俄州立大学Ann L. Griffen等人于2018年9月5日在<Microbiome>发表了题目为<High-resolution ISR amplicon sequencin ...