咨询热线:400-065-6886   天昊基因

中文 / English

主页 > 技术支持 > 科研进展 >

【绘图进阶】之lefse定制化绘图(五)


生信团 上海天昊生物

   lefse分析在微生物多样性中经常使用,但是由于数据库分类问题,我们会对其中一些没有意义的菌进行过滤,或挑选LDA值较高且P值显著的菌,重新绘制条形图。这里,我们展示了几种lefse的条形图。


 

0一

 

数据预处理


In [1]:
df = read.table('taxon.lefse.result.xls',sep = '\t',header = T)
head(df,4)

Out[1]:
Taxonomy Group LDA Pvalue
k__Archaea A 6.879678 0.03899022
k__Archaea.p__Euryarchaeota B 5.548066 0.03899022
k__Archaea.p__Euryarchaeota.c__Thermoplasmata B 5.497864 0.03899022
k__Archaea.p__Euryarchaeota.c__Thermoplasmata.o__Methanomassiliicoccales B 5.514124 0.02732372


   如上图所示,我们拿到的lefse的文件格式一般是这个样子的,为了减少文字长度,我们通常只显示最后一个分类层级的名称,需要用“.” 作为分隔符,取最后一列。linux上使用awk -F '.' '{print $NF}' taxon.lefse.result.xls 获取。结果如下所示:


数据读取

In [2]:
df = read.table('lefse.result.xls',sep = '\t',header = T)
head(df)

Out[2]:
Taxonomy Group LDA Pvalue
k__Archaea A 6.879678 0.0389902
2
p__Euryarchaeota B 5.548066 0.03899022
c__Thermoplasmata B 5.497864 0.03899022
o__Methanomassiliicoccales B 5.514124 0.02732372
f__Methanomassiliicoccaceae B 5.501945 0.02732372
g__Methanomassiliicoccus  B  5.520169  0.02732372

 

 

Lefse绘图


In [3]:
library(ggplot2)


1.按照LDA值进行排序

In [4]:
ggplot(df,aes(x = reorder(Taxonomy,LDA),y=LDA,fill = Group)) +
geom_bar(stat = 'identity',colour = 'black',width = 0.9,position = position_dodge(0.7))+
xlab('') + ylab('LDA') + coord_flip() +
theme_bw() +
theme(panel.grid.major.y = element_blank()) + theme(axis.text.y = element_blank(),
            axis.ticks = element_blank(),
axis.text.x = element_text(face = 'bold'),axis.title.x = element_text(face='bold')) +
theme(panel.border = element_blank())+
geom_text(aes(y = ifelse(df$LDA >0,-0.1,0.1),label=Taxonomy),fontface=4,size=1,hjust = ifelse(df$LDA>0,1,0))

Out[4]:



2. 按照分组和LDA值进行排序


In [5]:
df=df[order(df[,2],df[,3]),]
head(df)

Out[5]:
Taxonomy Group LDA Pvalue
161 g__Oceanibaculum A 3.942752 0.02413390
147 s__Methylobacterium_oxalidis A 4.063120 0.03794154
71 s__Actinocorallia_longicatena A 4.134189 0.04593813
116 s__Paenibacillus_edaphicus A 4.138317 0.03899022
26 s__Aquihabitans_daechungensis A 4.149702 0.04964780
88  g__Segetibacter  A  4.174600  0.03473526


In [6]:
df$Taxonomy = factor(df$Taxonomy,levels = as.character(df$Taxonomy))

In [7]:
df$Taxonomy = factor(df$Taxonomy,levels = as.character(df$Taxonomy))
ggplot(df,aes(x = Taxonomy,y=LDA,fill = Group)) +
  geom_bar(stat = 'identity',colour = 'black',width = 0.9,position = position_dodge(0.7))+
  xlab('') + ylab('LDA') + coord_flip() +
  theme_bw() +
  theme(panel.grid.major.y = element_blank()) +
  theme(axis.text.y = element_blank(),axis.ticks = element_blank(),
      axis.text.x = element_text(face = 'bold'),axis.title.x = element_text(face='bold')) +
  theme(panel.border = element_blank())+
  geom_text(aes(y = ifelse(df$LDA >0,-0.1,0.1),label=Taxonomy),fontface=4,size=1,hjust = ifelse(df$LDA>0,1,0))

Out[7]:


 


3. 带正负值的绘图

将B组的LDA值改成负数


In [8]:
df[df$Group=='A',3] = 0 - df[df$Group=='A',3]
head(df)

Out[8]:
Taxonomy Group LDA Pvalue
161 g__Oceanibaculum A -3.942752 0.02413390
147 s__Methylobacterium_oxalidis A -4.063120 0.03794154
71 s__Actinocorallia_longicatena A -4.134189 0.04593813
116 s__Paenibacillus_edaphicus A -4.138317 0.03899022
26 s__Aquihabitans_daechungensis A -4.149702 0.04964780
88  g__Segetibacter  A  -4.174600  0.03473526

In [9]:
df$Taxonomy = factor(df$Taxonomy,levels = as.character(df$Taxonomy))

In [10]:
ggplot(df,aes(x = Taxonomy,y=LDA,fill = Group)) +
  geom_bar(stat = 'identity',colour = 'black',width = 0.9,position = position_dodge(0.7))+
  xlab('') + ylab('LDA') + coord_flip() +
  theme_bw() +
  theme(panel.grid.major.y = element_blank()) +
  theme(axis.text.y = element_blank(),axis.ticks = element_blank(),
      axis.text.x = element_text(face = 'bold'),axis.title.x = element_text(face='bold')) +
  theme(panel.border = element_blank()) +
  geom_text(aes(y = ifelse(df$LDA >0,-0.1,0.1),label=Taxonomy),fontface=4,size=1,hjust = ifelse(df$LDA>0,1,0))

Out[10]:



4.绘制top10的正负区分条形图


In [11]:
df = rbind(head(df[df$Group=='A',],10),head(df[df$Group=='B',],10),head(df[df$Group=='C',],10))head(df)

Out[11]:
Taxonomy Group LDA Pvalue
161 g__Oceanibaculum A -3.942752 0.02413390
147 s__Methylobacterium_oxalidis A -4.063120 0.03794154
71 s__Actinocorallia_longicatena A -4.134189 0.04593813
116 s__Paenibacillus_edaphicus A -4.138317 0.03899022
26 s__Aquihabitans_daechungensis A -4.149702 0.04964780
88  g__Segetibacter  A  -4.174600  0.03473526


通过geom_text(size=4) 调节字体大小

In [12]:
df$Taxonomy = factor(df$Taxonomy,levels = as.character(df$Taxonomy))
ggplot(df,aes(x = Taxonomy,y=LDA,fill = Group)) +
  geom_bar(stat = 'identity',colour = 'black',width = 0.9,position = position_dodge(0.7))+
  xlab('') + ylab('LDA') + coord_flip() +
  theme_bw() +
  theme(axis.text.y = element_blank(),axis.ticks = element_blank(),
      axis.text.x = element_text(face = 'bold'),axis.title.x = element_text(face='bold')) +
  theme(panel.border = element_blank()) +
  geom_text(aes(y = ifelse(df$LDA >0,-0.1,0.1),label=Taxonomy),fontface=4,size=4,hjust = ifelse(df$LDA>0,1,0))

Out[12]:


 

关注“上海天昊生物”公众号,

回复关键字“lefse”,

获取文章链接,数据和代码的百度云盘链接。

 

往期相关链接:

1、R基础篇

excel不熟练怎么办,R来帮您(一)数据分类汇总
如何使用Rstudio练习R基础教程
R相关软件及R包安装
【零基础学绘图】之气泡图绘制(六)
【零基础学绘图】之绘制venn图(五)
【零基础学绘图】之绘制barplot柱状图图(四)
【零基础学绘图】之绘制heatmap图(三)
【零基础学绘图】之绘制PCA图(二)
【零基础学绘图】之alpha指数箱体图绘制(一)

2、R进阶

【绘图进阶】之六种带中心点的PCA 图和三维PCA图绘制(四)

【绘图进阶】之交互式可删减分组和显示样品名的PCA 图(三)

【绘图进阶】之绘制PCA biplot图(二)

【进阶篇绘图】之带P值的箱体图、小提琴图绘制(一)

3、数据提交

3分钟学会微生物多样性云平台数据分析;

3分钟学会CHIP-seq类实验测序数据可视化 —IGV的使用手册

10分钟搞定多样性数据提交,最快半天内获取登录号,史上最全的多样性原始数据提交教程

20分钟搞定GEO上传,史上最简单、最详细的GEO数据上传攻略

4、表达谱分析

表达谱分析(二)通路富集分析和基因互作网络图绘制
如何对GEO数据进行差异分析
miRNA靶基因预测软件__miRWalk 3.0

5、医学数据分析

KING: 样本亲缘关系鉴定工具
【WGS服务升级】人工智能软件SpliceAI助力解读罕见和未确诊疾病中的非编码突变
隐性疾病trio家系别忽视单亲二倍体现象——天昊数据分析助力临床疾病诊断新添UPD(单亲二倍体)可视化分析工具
【昊工具】Oh My God! 太好用了吧!疾病或表型的关键基因查询数据库,我不允许你不知道Phenolyzer
 
天昊客户服务中心
手机/微信号:18964693703
 
 

【本群将为大家提供】

分享生信分析方案

提供数据素材及分析软件支持

定期开展生信分析线上讲座

QQ号:1040471849

 
 

作者:大熊

审核:有才

来源:天昊生信团

 

微信扫一扫
关注该公众号

 

 




上海天昊生物科技有限公司 版权所有 沪ICP备17008908号
地址:上海市浦东新区康桥路787号9号楼 邮箱:techsupport@geneskies.com 电话:400-065-6886