-
2018-12-07 21:36:24
转录组GO富集与微生物相关性分析
原始数据格式
我相信做生物分析的都知道很多时候难点往往不在于代码怎么写,而在于数据格式很多时候与各种包的要求不相符合,因此我先列一下本次分析中用到的数据格式
表1:转录组数据使用原始矩阵
A_1 A_2 A_3 B_1 B_2 B_3 C_1 C_2 C_3 D_1 …… Gene0000001 …… …… …… …… …… …… …… …… …… …… …… Gene0000002 …… …… …… …… …… …… …… …… …… …… …… Gene0000003 …… …… …… …… …… …… …… …… …… …… …… Gene0000004 …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… …… 表2:分组数据
sample_ID group A_1 A A_2 A A_3 A B_1 B B_2 B B_2 B C_1 C …… …… 表3:细菌丰度表
genus1 genus2 genus3 genus4 …… A_1 …… …… …… …… …… A_2 …… …… …… …… …… A_3 …… …… …… …… …… B_1 …… …… …… …… …… B_2 …… …… …… …… …… …… …… …… …… …… …… 表4:基因注释文件
由于我做的是非模式生物,所以进行GO富集时还需要一个基因注释文件,这个文件的后缀是.map(用文本文档按下列格式弄好改后缀就行),gene_id列与GO_number列之间用\t制表符分隔gene_id GO_number Gene0000001 GO:0005515 Gene0000001 GO:0051259, GO:0046872, GO:0030176, GO:0005783, GO:0030968, GO:0016874, GO:0016021, GO:0016020, GO:0008270, GO:0000299, GO:0006928, GO:0006512, GO:0006511, GO:0007165, GO:0005515, GO:0030433, GO:0005789, GO:0004842, GO:0000209, GO:0004872 Gene0000003 GO:0008076, GO:0005242, GO:0030955, GO:0005244 …… …… 使用TCC包进行差异基因分析
TCC是什么包?为什么要用TCC包?目前常用的差异基因分析包为DEseq2或者edgeR,但是这两个包在寻找差异基因时只能够找出两组之间的差异基因,即使你的数据里有A、B、C、D四个组,DESeq2分析后也只会返回A vs B,A vs C,A vs D,B vs C……的结果,而不会找出在四个组中有显著性差异的基因。这对于做临床研究、药理研究的同志们来说可能够了,但是对于植物学方面的就不够用了(不同器官、不同生长阶段等)。那么有没有可以多组同时比较的R包呢?经过检索我发现了这个TCC包可以满足要求1。其实这个包是在DESeq,DESeq2和edgeR三个包基础上构建了一个适用于多组分析的pipeline。这里我使用的是基于DESeq2的管道,更多内容可以参考帮助文档。分析代码如下:
#gene_count_matrix,为基因表达矩阵,读入过程略 #group为分组,读入过程略 library(TCC) #构建TCC类 tcc <- new("TCC", gene_count_matrix, group) #下面使用TCC包基于DESeq2的管道方案 ssss <- calcNormFactors(tcc, norm.method = "deseq2", test.method = "deseq2", iteration = 3) DE <- estimateDE(ssss, test.method = "deseq2", full = ~ group, reduced = ~ 1) res <- getResult(DE)
res就是分析后的结果了,根据names可以很清楚这个结果包含哪些信息
使用topGO包进行GO富集分析
这里没什么说的,就是一个典型的topGO分析,有不明白的可以参考其他博主的内容
#准备基因列表,注意,一定要因子化!!! genelist<-res$estimatedDEG names(genelist)<-res$gene_id genelist<-as.factor(genelist) #读入基因注释文件 gene2GO_map<-readMappings("D:/gene2GO.map") #构建topGO类,注意MF。CC和BP是分开构建三个。 Godata_BP <- new( "topGOdata", ontology = "BP", allGenes = genelist, annot = annFUN.gene2GO, gene2GO = gene2GO_map) #寻找差异GO resultFis_BP <- runTest(Godata_BP, algorithm = "classic", statistic = "fisher") #找出排名前100的GO sig.tab <- GenTable( Godata_BP, Fis = resultFis_BP, topNodes = 100)
计算转录组与微生物组相关性
这边由于微生物和转录组在两个不同的矩阵中,所以自己写了个函数进行相关性计算,比较懒直接用了循环嵌套,速度有些慢,用apply族写的话速度应该会快很多。需要注意的是,循环数小的应该放在外层,循环数大的应该放在内层
get_matrix_cor<-function(x,y) { cormatrix<-as.data.frame(matrix(numeric(0),ncol=length(x),nrow = length(y))) for (i in (1:length(x))){ for (j in (1:length(y))){ cormatrix[j,i]<-cor(x[,i],y[,j]) } } }
计算相关性矩阵
trans_bac_cor_matrix<-get_matrix_cor(bac_genus_matrix,t(gene_count_matrix)) %>% as.data.frame()
根据GO号提取相关性矩阵并绘图
一个转录组数据里面的gene数以万计,画在一张图上就什么都看不见了,因此根据我感兴趣的GO号提取出来画图会好的多。GO号的选取可以根据富集结果的sig.tab选取,这里重点讲述画图和输出结果的过程
提取相关性矩阵
#输入一个GO号 GO="0005985" #获取感兴趣的GO号相对应的gene名字 a<-genesInTerm(Godata_BP,paste0("GO:",GO))[[1]] %>% as.vector() #从相关系数矩阵中提取对应基因的行,建议用apply族函数改写,简化代码 b<-c() for (i in 1:length(a) ){ b[i]<-which(row.names(trans_bac_cor_matrix)==a[i]) } matrix_BF<-trans_bac_cor_matrix[b,]
这里的matrix_BF就是转录组与微生物组的相关性矩阵了,后期用来进行热图绘制。但是用pheatmap包进行热图绘制前需要先去除NA
#删除NA matrix_BF<-matrix_BF[-which(is.na(matrix_BF)),] na_flag <- apply(is.na(matrix_BF), 2, sum) matrix_BF<-matrix_BF[,na_flag==0] #注意,我这里删除NA的策略是将某种存在NA的菌直接删去,所以如果有某个基因全是 #NA的话,那matrix_BF就会成为一个空数据框,这个时候需要将这行删去。
绘制热图
#绘制热图 annotation_col = data.frame( ClassBac = factor(paste0('Cluster',cutree(hclust_BF,10))) ) rownames(annotation_col) = colnames(matrix_BF) pdf(paste0("D:/GO_",GO,".pdf")) heatmap_BF<-pheatmap::pheatmap(matrix_BF,show_colnames = F,annotation_col = annotation_col) dev.off()
这里展示一张(涉及数据因此基因名盖住了)
输出与热图对应的相关性矩阵
pheatmap会自动聚类,这样绘制完成的热图就与之前提取的相关性矩阵不对应了,因此在这里进行一个处理,使得我们输出的矩阵与热图顺序一致。该部分内容参考了一位博主的方法,但是后来找不到页面了,如有侵权,与我联系
hclust_BF=hclust(dist(t(matrix_BF))) hclust_BF_row=hclust(dist(matrix_BF)) col_cluster=cutree(hclust_BF,k=10) row_cluster=cutree(hclust_BF_row,k=2) newOrder=matrix_BF[hclust_BF_row$order,] newOrder[,ncol(newOrder)+1]=row_cluster[match(row.names(newOrder),names(row_cluster))] names(newOrder)[ncol(newOrder)]="Cluster" names2=names(matrix_BF[,hclust_BF$order]) newOrder[nrow(newOrder)+1,]=col_cluster[match(c(names2,"cluster"),names(col_cluster))] row.names(newOrder)[nrow(newOrder)]="Cluster" #输出对应矩阵 write.csv(newOrder,paste0("D:/GO_",GO,".csv"))
这样输出的矩阵就与热图上的位置一致了。
更多相关内容 -
Go语言使用组合的方式实现多继承的方法
2020-09-22 10:00:58主要介绍了Go语言使用组合的方式实现多继承的方法,实例分析了多继承的原理与使用组合方式来实现多继承的技巧,需要的朋友可以参考下 -
Go语言实现的排列组合问题实例(n个数中取m个)
2020-09-21 10:22:44主要介绍了Go语言实现的排列组合问题,结合实例形式分析了Go语言实现排列组合数学运算的原理与具体操作技巧,需要的朋友可以参考下 -
Go语言中接口组合的实现方法
2020-09-22 10:01:29主要介绍了Go语言中接口组合的实现方法,实例分析了接口中包含接口的实现技巧,需要的朋友可以参考下 -
一张图展示多组富集分析结果
2022-03-30 08:39:15单细胞之富集分析-3:GO和KEGG富集分析及绘图 单细胞之富集分析-4:分组水平GSVA 文献里经常看到这样不同组的通路选择性的一起展示的点图,只要有不同分组的通路富集结果其实就可以画了。简单做一下复现吧~ 还是...单细胞富集分析系列:
单细胞之富集分析-1:单细胞GSEA分析流程
单细胞之富集分析-2:批量GSEA和GSVA分析
单细胞之富集分析-3:GO和KEGG富集分析及绘图
单细胞之富集分析-4:分组水平GSVA
文献里经常看到这样不同组的通路选择性的一起展示的点图,只要有不同分组的通路富集结果其实就可以画了。简单做一下复现吧~
还是使用我们熟悉的pbmc3k的数据做演示。因为这个数据集只有一个样品,所以展示的还是不同细胞群间的通路。多样本的数据集比较不同样本/分组间的差异只需要改一下Idents就可以了。
library(Seurat)
library(patchwork)
library(clusterProfiler)
library(org.Mm.eg.db) ##加载小鼠
library(org.Hs.eg.db) ##加载人类
library(tidyverse)导入注释好的pbmc数据集
pbmc <- readRDS(“pbmc.rds”)
寻找每个群的marker基因
Idents(pbmc) <- ‘cell_type’#因为这个是单样品数据集,在做多样品的时候,把Idents设置为不同样品即可计算不同样品的差异基因
table(pbmc$cell_type)Naive CD4 T Memory CD4 T CD14+ Mono B CD8 T
711 480 472 344 279
FCGR3A+ Mono NK DC Platelet
162 144 32 14
markers <- FindAllMarkers(pbmc)
sig_dge.all <- subset(markers, p_val_adj<0.05&abs(avg_log2FC)>0.15) #所有差异基因
#saveRDS(sig_dge.all, file = “deg.rds”)
View(sig_dge.all)
##CD4T的上调基因
sig_dge.CD4T_up <- subset(sig_dge.all, avg_log2FC>0.15&cluster==‘Naive CD4 T’) #CD4T细胞的上调差异基因
#富集分析
ego_CD4 <- enrichGO(gene = row.names(sig_dge.CD4T_up),
#universe = row.names(dge.celltype),
OrgDb = ‘org.Hs.eg.db’,
keyType = ‘SYMBOL’,
ont = “ALL”, #设置为ALL时BP, CC, MF都计算
pAdjustMethod = “BH”,
pvalueCutoff = 0.01,
qvalueCutoff = 0.05)
ego_cd4 <- data.frame(ego_CD4)
#write.csv(ego_CD4,‘enrichGO_CD4.csv’)
View(ego_cd4)
ego_cd4$Group=‘Naive CD4 T’一样的方法计算其他细胞群的富集结果,可以写成循环
sig_dge.B_up <- subset(sig_dge.all, avg_log2FC>0.15&cluster==‘B’) #CD4T细胞的上调差异基因
#富集分析
ego_B <- enrichGO(gene = row.names(sig_dge.B_up),
#universe = row.names(dge.celltype),
OrgDb = ‘org.Hs.eg.db’,
keyType = ‘SYMBOL’,
ont = “ALL”, #设置为ALL时BP, CC, MF都计算
pAdjustMethod = “BH”,
pvalueCutoff = 0.01,
qvalueCutoff = 0.05)
ego_b <- data.frame(ego_B)
#write.csv(ego_B,‘enrichGO_B.csv’)
View(ego_b)
ego_b$Group=‘B’sig_dge.MCD4T_up <- subset(sig_dge.all, avg_log2FC>0.15&cluster==‘Memory CD4 T’) #CD4T细胞的上调差异基因
#富集分析
ego_MCD4T <- enrichGO(gene = row.names(sig_dge.MCD4T_up),
#universe = row.names(dge.celltype),
OrgDb = ‘org.Hs.eg.db’,
keyType = ‘SYMBOL’,
ont = “ALL”, #设置为ALL时BP, CC, MF都计算
pAdjustMethod = “BH”,
pvalueCutoff = 0.01,
qvalueCutoff = 0.05)
ego_mcd4t <- data.frame(ego_MCD4T)
#write.csv(ego_MCD4T,‘enrichGO_MCD4T.csv’)
View(ego_mcd4t)
ego_ego_mcd4t$Group=‘Memory CD4 T’
每个群选了前3个通路,这里的通路可以随便选CD4T=ego_cd4[1:3,]
B=ego_b[1:3,]
MCD4T=ego_mcd4t[1:3,]
all <- rbind(CD4T,B,MCD4T)
all G e n e R a t i o < − c ( 0.3798883 , 0.3798883 , 0.3854749 , 0.2105263 , 0.1184211 , 0.09210526 , 0.1502591 , 0.1450777 , 0.0880829 ) l i b r a r y ( f o r c a t s ) a l l GeneRatio <- c(0.3798883,0.3798883,0.3854749,0.2105263,0.1184211,0.09210526,0.1502591,0.1450777,0.0880829) library(forcats) all GeneRatio<−c(0.3798883,0.3798883,0.3854749,0.2105263,0.1184211,0.09210526,0.1502591,0.1450777,0.0880829)library(forcats)allDescription <- as.factor(all D e s c r i p t i o n ) a l l Description) all Description)allDescription <- fct_inorder(all$Description)
最终的目的是得到这样的表格,挑选想展示的通路,标记清楚分组即可ggplot(all, aes(Group, Description)) +
geom_point(aes(color=p.adjust, size=GeneRatio))+theme_bw()+
theme(panel.grid = element_blank(),
axis.text.x=element_text(angle=90,hjust = 1,vjust=0.5))+
scale_color_gradient(low=’#6699CC’,high=’#CC3333’)+
labs(x=NULL,y=NULL)+guides(size=guide_legend(order=1))作者:Hayley笔记
链接:https://www.jianshu.com/p/9acaa566f04b
来源:简书
著作权归作者所有。商业转载请联系作者获得授权,非商业转载请注明出处。 -
R 实战 | 使用clusterProfiler进行多组基因富集分析
2021-11-09 01:17:28R 实战 | 使用clusterProfiler进行多组基因富集分析clusterProfiler这个包我就不再介绍了,网上关于用这个包做的基础的富集分析的教程已经非常多了,今天主要介绍一...R 实战 | 使用clusterProfiler进行多组基因富集分析
clusterProfiler
这个包我就不再介绍了,网上关于用这个包做的基础的富集分析的教程已经非常多了,今天主要介绍一下使用compareCluster
功能进行多组基因的富集分析。示例数据
富集分析
多个基因集的富集分析
多个分组的富集分析
基本用法
参数
References
往期
示例数据
library(clusterProfiler) data(gcSample) #载入 str(gcSample) #数据集格式 lapply(gcSample, head) #查看数据集
> str(gcSample) List of 8 $ X1: chr [1:216] "4597" "7111" "5266" "2175" ... $ X2: chr [1:805] "23450" "5160" "7126" "26118" ... $ X3: chr [1:392] "894" "7057" "22906" "3339" ... $ X4: chr [1:838] "5573" "7453" "5245" "23450" ... $ X5: chr [1:929] "5982" "7318" "6352" "2101" ... $ X6: chr [1:585] "5337" "9295" "4035" "811" ... $ X7: chr [1:582] "2621" "2665" "5690" "3608" ... $ X8: chr [1:237] "2665" "4735" "1327" "3192" ... > lapply(gcSample, head) $X1 [1] "4597" "7111" "5266" "2175" "755" "23046" $X2 [1] "23450" "5160" "7126" "26118" "8452" "3675" $X3 [1] "894" "7057" "22906" "3339" "10449" "6566" $X4 [1] "5573" "7453" "5245" "23450" "6500" "4926" $X5 [1] "5982" "7318" "6352" "2101" "8882" "7803" $X6 [1] "5337" "9295" "4035" "811" "23365" "4629" $X7 [1] "2621" "2665" "5690" "3608" "3550" "533" $X8 [1] "2665" "4735" "1327" "3192" "5573" "9528"
富集分析
多个基因集的富集分析
xx <- compareCluster(gcSample, fun="enrichKEGG", organism="hsa", pvalueCutoff=0.05) table(xx@compareClusterResult$Cluster) #每个基因集富集个数 head(as.data.frame(xx)) #查看完整结果
> table(xx@compareClusterResult$Cluster) X1 X2 X3 X4 X5 X6 X7 X8 0 3 3 18 10 1 15 19 > head(as.data.frame(xx)) Cluster ID Description GeneRatio 1 X2 hsa04110 Cell cycle 18/384 2 X2 hsa05169 Epstein-Barr virus infection 23/384 3 X2 hsa05340 Primary immunodeficiency 8/384 4 X3 hsa04512 ECM-receptor interaction 9/187 5 X3 hsa04060 Cytokine-cytokine receptor interaction 17/187 6 X3 hsa04151 PI3K-Akt signaling pathway 19/187 BgRatio pvalue p.adjust qvalue 1 126/8105 2.441211e-05 0.007470105 0.006989572 2 202/8105 7.911793e-05 0.012105043 0.011326356 3 38/8105 3.297441e-04 0.033633898 0.031470314 4 88/8105 1.815667e-04 0.045098637 0.042158192 5 295/8105 4.490651e-04 0.045098637 0.042158192 6 354/8105 5.048355e-04 0.045098637 0.042158192 geneID 1 991/1869/890/1871/701/990/10926/9088/8317/9700/9134/1029/2810/699/11200/23594/8555/4173 2 4067/3383/7128/1869/890/1871/578/864/637/9641/6891/355/9134/5971/916/956/6850/7187/3551/919/4734/958/6772 3 100/6891/3932/973/916/925/958/64421 4 7057/3339/1299/3695/1101/3679/3910/3696/3693 5 2919/4982/3977/6375/8200/608/8792/3568/2057/1438/8718/655/652/10220/50615/51561/7042 6 894/7057/6794/2247/1299/3695/2252/2066/1101/8817/1021/5105/3679/3082/2057/3910/3551/3696/3693 Count 1 18 2 23 3 8 4 9 5 17 6 19
dotplot(xx) #气泡图
多个分组的富集分析
示例数据
data(geneList, package="DOSE") #载入DOSE包中的数据 head(geneList) #查看数据 包含基因名及其foldchange mydf <- data.frame(Entrez=names(geneList), FC=geneList) # 以下内容目的是构建分组 也可以用别的分组 # 将FC大于1的标注为上调 反之为下调 mydf <- mydf[abs(mydf$FC) > 1,] mydf$group <- "upregulated" mydf$group[mydf$FC < 0] <- "downregulated" # 将FC绝对值大于2 的标注为B 反之为A mydf$othergroup <- "A" mydf$othergroup[abs(mydf$FC) > 2] <- "B" head(mydf) # 查看示例数据(两个分组)
> head(mydf) Entrez FC group othergroup 4312 4312 4.572613 upregulated B 8318 8318 4.514594 upregulated B 10874 10874 4.418218 upregulated B 55143 55143 4.144075 upregulated B 55388 55388 3.876258 upregulated B 991 991 3.677857 upregulated B
分析及其可视化
# 可以进行单组或多组分析,在 + 号后添加即可 formula_res <- compareCluster(Entrez~group+othergroup, data=mydf, fun='enrichKEGG') dotplot(formula_res)
# 同样可以进行分面操作 dotplot(formula_res, x=~group) + ggplot2::facet_grid(~othergroup) #对于只用一次该包的功能可以这么写
基本用法
最后附上用法参数。
compareCluster(geneClusters, fun = "enrichGO", data = "", ...)
参数
geneClusters
entrez基因list. 或者, Entrez~分组 形式的列表 fun
"groupGO", "enrichGO", "enrichKEGG", "enrichDO" or "enrichPathway" . data
数据集 References
Chapter 11 Biological theme comparison | clusterProfiler: universal enrichment tool for functional and comparative study (hiplot.com.cn)
往期
-
Protocols-4pub:Lyu缔结的多组学分析协议
2021-04-27 15:22:14(北京大学成力实验室, )提出的多组学分析方案。 声明:该项目例外学习用途,请共同,禁止进行商业盈利。 上存在集成的分步渐晕。 这些文件包括三个主要部分: 预处理管道(等bulkRNA,bulkATAC,ChIP,WGBS等... -
宏基因组分析-基于组装
2022-05-26 16:59:03宏基因组的分析目前主要包括三种方法:基于组装分析、基于reads分析、基于bin分析。 下面我们介绍基于组装的分析方法。 二、分析流程介绍 数据分析从下机原始序列开始,首先对原始序列进行去接头、 质量剪切以及...一、介绍
宏基因组 ( Metagenome) 指特定环境下所有生物遗传物质的总和。它包含了可培养的和未可培养的微生物的基因。一般从环境样品中提取基因组DNA, 进行高通量测序,从而分析微生物多样性、种群结构、功能信息、与环境之间的关系等。
宏基因组的分析目前主要包括三种方法:基于组装分析、基于reads分析、基于bin分析。
下面我们介绍基于组装的分析方法。
二、分析流程介绍
数据分析从下机原始序列开始,首先对原始序列进行去接头、 质量剪切以及去除污染等优化处理。然后使用优质序列进行拼接组装和基因预测,将各样本预测得到的基因集合并在一起去冗余,得到非冗余基因集;对得到的非冗余基因集进行物种和功能上的注释,并使用BWA软件将优化序列比对到非冗余基因集,计算得到各基因在各样品中的丰度信息(RPKM); 对物种和功能注释结果进行统计分析。
三、详细流程
- 使用fastp软件使用划框方法去除低质量碱基,同时去除接头序列;如果样品来源于宿主(比如人或动物的粪便),而且该宿主本身的基因组已被发表, 则通过软件Bowtie2将reads比对宿主DNA序列,并去除比对相似性高的污染reads;
- MEGAHIT是一个二代测序从头组装工具,尤其在土壤等复杂环境样本组装、大量样本混合组装方面优势明显,同时提供更好的完整性和连续性,为行业的主流组装软件。使用Megahit软件通过设置不同kmer参数,对优化序列进行组装得到Contigs;可以通过N50判断组装结果的质量。
样品名称
序列数(条)
碱基数(bp)
Mean(bp)
Max(bp)
N50(bp)
N70(bp)
N90(bp)
Sample1
118571
2.47E+08
2079
736220
13365
2240
598
Sample2
160636
2.74E+08
1703
764995
7111
1569
513
Sample3
120704
2.23E+08
1845
627035
10082
1827
536
- Prodigal是一种采用无监督的机器学习算法,用于原核基因组的蛋白质编码基因预测软件工具。我们使用Prodigal软件对Contigs进行基因预测(ORF);
- 使用CD-Hit软件对基因集去冗余,得到非冗余基因集 (non-redundant gene catalog);
- 基因丰度统计:使⽤BWA软件,将质控优化序列⽐对到非冗余基因集上,从比对上的 序列数目及基因长度出发,计算得到各基因在各样品中的丰度信息(RPKM);
- 使用DIAMOND软件将非冗余基因集与 NCBI 的 NR序列进行blastp比对,得到基因序列对应的物种注释信息;
- 使用DIAMOND软件将非冗余基因集与ARDB、CAZy、KEGG、eggNOG、VFDB等数据库进行blastp比对,得到基因序列对应的功能信息。
四、主要结果介绍
1、基于基因的分析
1.1 样品间相关性分析
生物学重复是任何生物学实验所必须的,高通量测序技术也不例外。 样品间基因丰度相关性是检验实验可靠性和样本选择是否合理性的重要指标。 相关系数越接近1,表明样品之间基因丰度模式的相似度越高。
注:图中不同颜色代表 spearman 相关系数的高低; 相关系数与颜色间的关系见右侧图例说明; 颜色越深代表样品间相关系数的绝对值越大; 椭圆向左偏表明相关系数为正,右偏为负; 椭圆越扁说明相关系数的绝对值越大。
2、基于物种的分析
从不同分类层级(门纲目科属种)的相对丰度表出发,展示相对丰度前25的物种,其余的合并为Others,绘制出各样品对应的物种注释结果在不同分类层级上的相对丰度柱形图。
注:横坐标为样本,纵坐标为各级别的相对丰度;其它级别的图请在下面目录查看。
3、功能分析
抗生素抗性基因数据库(ARDB,ARDB-Antibiotic Resistance Genes Database) 收集了不同环境来源的(如肠道、生活废水、河流等)细菌抗药性 基因及其抗性谱、作用机制、本体论、COG和CDD等注释信息,为研 究药物作用、环境治理提供研究依据。
抗性基因类型统计(Type 左) (Class 右)
The Comprehensive Antibiotic Resistance Database (http://arpcard.mcmaster.ca)CARD 数据库核心是 ARO(Antibiotic Resistance Ontology),ARO 包含了与抗生素抗性基因,抗性机制,抗生素和靶相关的term。CARD 数据库已成为目前最受欢迎的耐药基因研究工具之一。
抗性基因统计(Drug_Class左) (Resistance_Mechanism右)
抗性基因统计(Best_Hit_ARO左) (AMR_Gene_Family 右)
3.2 碳水化合物酶功能注释(CAZY)
碳水化合物活性酶(Carbohydrate-active enzymes,CAZyme) 对地球上所有碳水化合物的合成、降解与修饰起重要作用, 因此深入研究CAZyme,对于了解微生物碳水化合物的代谢机制非常重要。
注: 横坐标6大功能类,纵坐标基因丰度,颜色表示功能。
KEGG数据库(Kyoto Encyclopedia of Genes and Genomes, 京都基因和基因组百科全书,KEGG: Kyoto Encyclopedia of Genes and Genomes) 是系统分析基因功能,联系基因组信息和功能信息的大型知识库。 KEGG GENES数据库提供关于在基因组计划中发现的基因和蛋白质的序列信息; KEGG PATHWAY数据库包括各种代谢通路、合成通路、膜转运、信号传递、 细胞周期以及疾病相关通路等。
注: 横坐标每个样本,纵坐标基因丰度,颜色表示名称。
3.4 eggNOG功能注释
eggNOG 数据库:是国际上普遍认可的同源聚类基因群的专业注释数据库,包括来自原始COG/KOG的功能分类,以及基于分类学的功能注释。
注: 横坐标每个功能,纵坐标基因丰度,颜色表示每个功能。
GO数据库分别从功能、参与的生物途径及细胞中的定位对基因产物进行了标准化描述,所谓的 GO,是生物学功能注释的一个标准词汇表术语(GO term),将基因的功能分为三部分: 基因执行的分子功能( Molecular Function), 基因参与的生物学过程(Biological Process), 基因所处的细胞组分( Cellular Component)。
注: 横坐标表示 Level2 水平上 GO term,不同颜色代表不同 GO 类,分为三大类,纵坐标表示 注释到该 term 上面的 Unigene 数目
4 Alpha多样性:
Alpha 多样性( Alphadiversity)包括 :
chao 指数和 ACE 指数反映样品中群落的丰富度( speciesrichness),即简单指群落中物种的数量,而不考虑群落中每个物种的丰度情况。
shannon 指数以及 simpson 指数反映群落的多样性( speciesdiversity),受样品群落中物种丰富度( speciesrichness)和物种均匀度( speciesevenness)的影响。
相同物种丰富度的情况下,群落中各物种具有越大的均匀度,则认为群落具有越大的多样性。
5 Beta多样性分析
5.1 基于基因分析
Beta多样性(Beta diversity)分析是用来比较一对样品在物种多样性方面存在的差异大小。分析各类群在样品中的含量,进而计算出不同样品间的Beta多样性值。多种指数可以衡量Beta多样性,常用的为Bray-Curtis,Jaccard。我们将以这两种方法得到的矩阵,做NMDS,PCoA及聚类树分析。
5.2 基于物种分析
左边是样品间基于群落组成的层次聚类分析(Bray-curtis 算法),右边是样品的群落结构柱状图。
注:左边是样品间基于群落组成的层次聚类分析(Bray-curtis 算法),右边是样品的群落结构柱状图。
6差异分析
如果老师提供了分组信息,我们可以从物种、功能、alpha、beta多样性不同的角度进行差异分析。包括参数检验、非参数检验、lefse分析等。
相关性热图
7 相关性分析
功能和物种之间的相关性、物种和临床指标之间的相关性,只要是同样的样本检测了两种不同的指标,都可以进行两个指标之间的关联分析。展示方式主要是spearman热图和网络图的形式。
图中
(1)颜色代表相关性系数:蓝色为负相关,且颜色越深,相关性系数越大,粉红色为正相关,且颜色越深,相关性系数越大,具体的颜色与相关性系数的对应关系见图中右上角的图注。
(2)横轴代表的是临床因子,纵轴代表的是代谢物。
(3)图中的*代表P值,*为0.05>P>0.01,**为0.01>P>0.001…..,只要是图中标注了*的,都是有显著相关的。
(4)左侧和上面的树都是根据相关性系数的相似性情况进行聚类的。
图中:
(1)图中形状代表不同的检测类型:本图是细菌和真菌,也可以是物种和不同功能之间的相关性。
(2)spearman的结果中p小于0.05,r绝对值大于0.6.
(3)节点的大小代表相对丰度的均值或是相关性程度。
(4)图中节点的颜色代表不同的门。注意:网络图中展示的信息需要根据老师的要求来设计。
-
Majorbio Cloud:一站式多组学数据分析平台
2022-04-07 12:32:44点击蓝字 关注我们Majorbio Cloud:一站式多组学数据分析平台https://doi.org/10.1002/imt2.12全文解读前 言随着高通量检测技术的发展,高通量组学数据生信平台应运而生,MG-RAST,Qiita,BIGDdb,TRAPR,imageGP,... -
利用R语言绘制热图、韦恩图、GO富集分析图(有了转录组数据不知道该怎么写文章,看我就对了!)
2020-12-30 06:18:52GO富集分析的R包和网页工具很多,我习惯用clusterProfiler包,因为它提供了很多画图函数,比较方便。genes intersect(RIP$`Gene name`,RNA_seq$`gene name`) #1756个基因library(clusterProfiler)library(org... -
【转载】无参转录组GO、KEGG富集分析——diamond+idmapping+GOstats
2021-08-30 15:08:02(1) 用无参转录组分析软件得到unigene fasta文件,命名为my_unigenes.fa,格式如下表所示: >MSTRG.5.1 gene=MSTRG.5 TGATGTCATCGATCCGTGACGTTTAGTATTCAACCAATAGGAATCAACCACGTAGGATTGGCGATCCTCG ... -
go语言介绍及应用场景分析
2022-07-15 16:26:13golang特定、应用场景及简单的使用介绍 -
理解差异表达与GO分析-Go语言中文社区
2021-06-04 09:51:08基因表达差异的显著性分析简称表达差异分析,其目的是比较两个条件(包括种属、表型等)下的基因表达差异,通过一定的统计学方法,从中识别出与条件相关的特异性基因,然后进一步分析这些特异性基因的生物学意义。... -
R语言GO富集分析
2022-06-12 10:26:17GO富集分析 —— 差异圈图、聚类圈图一、数据链接 链接: https://pan.baidu.com/s/1PPUW5YyJHjwxvkjmMi0Vtwpwd=c9s8 提取码: c9s8二、安装包介绍 clusterProfiler包:功能也比较强大,主要是做GO和KEGG的功能富集... -
无参转录组GO、KEGG富集分析——diamond+idmapping+GOstats
2018-03-29 18:54:19(1)用无参转录组分析软件得到unigene fasta文件,命名为my_unigenes.fa,格式如下表所示: >MSTRG.5.1 gene=MSTRG.5 TGATGTCATCGATCCGTGACGTTTAGTATTCAACCAATAGGAATCAACCACGTAGGATTGGCGATCCTCG ... -
GO富集分析(转载)-Go语言中文社区
2021-05-17 18:32:48GO富集介绍每个基因都会对应有一个或多个GO term(也就是GO功能)。富集涉及到两个概念:前景基因和背景基因。前景基因就是你关注的要重点研究的基因集,背景基因就是所有的基因集。比如做两个样本对照组和处理组的... -
Crosshub:癌症基因组图谱(TCGA)项目数据集的多向分析-开源
2021-04-29 14:25:55Crosshub可以对癌症基因组图谱(TCGA)项目的RNA-Seq,miRNA-Seq和甲基化组数据进行多路分析:1.差异表达分析(基因,替代转录本和miRNA)2.调控性miRNA预测(TargetScan,DIANA microT ,mirSVR,PicTar,... -
Go pprof内存指标含义备忘录及案例分析
2020-12-20 18:05:14最近组内一些Go服务碰到内存相关的问题,所以今天抽时间看了下Go pprof内存指标的含义,为后续查问题做准备。 内容主要来自于Go代码中对这些字段的注释,加自己的理解。理解不对的地方欢迎指正。 // ... -
转录组分析技术
2020-03-19 11:54:30@TOC1.转录组 @1.转录组 -
an吟:时程转录组学分析
2021-02-08 04:10:11an吟:时程转录组学分析 -
论文研究 - 泌乳期和干燥期奶牛瘤胃上皮细胞转录组变化的组装与分析
2020-05-27 23:42:40本报告介绍了在干燥和哺乳期使用新一代测序技术组装和分析牛牛瘤胃上皮组织转录组的过程。 转录组学分析和比较显示,由于对泌乳的适应性,与瘤胃上皮组织代谢相关的基因表达发生了广泛变化。 泌乳期间瘤胃上皮适应... -
一些GO及KEGG分析的知识
2021-03-07 22:01:01什么是GO分析?Gene Ontology(简称GO)是一个国际标准化的基因功能分类体系,提供了一套动态更新的标准词汇表(controlled vocabulary)来全面描述生物体中基因和基因产物的属性。GO总共有三个ontology(本体),分别描述... -
Go-Gonum是Go编程语言的一组数字库。它包含用于矩阵统计优化等的库
2019-08-13 11:20:16Gonum是Go编程语言的一组数字库。 它包含用于矩阵,统计,优化等的库 -
对转录组测序数据进行分析以及注释
2021-05-15 03:25:54随着二代测序技术的高速发展,人们获得了大量的转录组数据序列如何从数据中挖掘具有生物意义的信息已经成为很多研究的关键,对未知基因的功能进行预测和注释就是一个重要问题这篇文章主要是跟着 刘粉香,杨文国,孙... -
转录组-蛋白组-代谢组关联分析
2021-01-25 16:57:10对生物体内生命过程中产生的一系列代谢产物做全面的分析有助于揭示基因型和表型之间的联系,整合多组学分析是目前综合分析代谢产物的最有效的方法。转录组与蛋白质组数据依据 mRNA 与蛋白之间的翻译关系彼此关联,... -
基于GO分析的全基因组芯片筛选小鼠牙囊发育早期差异表达基因的研究.pdf
2021-07-26 16:17:41基于GO分析的全基因组芯片筛选小鼠牙囊发育早期差异表达基因的研究.pdf -
Golang算法之田忌赛马问题实现方法分析
2020-12-26 08:32:12输入有多组测试数据。 每组测试数据包括3行: 第一行输入N(1≤N≤1000),表示马的数量。 第二行有N个整型数字,即渊子的N匹马的速度(数字大表示速度快)。 第三行有N个整型数字,即对手的N匹马的速度。 当N为0时... -
GO和KEGG富集分析详细步骤
2022-05-20 12:38:42GO和KEGG富集分析 文章目录GO和KEGG富集分析@[toc]1. 将差异表达结果的基因名称转化为id2. GO富集分析3. GO圈图绘制4. KEGG富集分析5. KEGG圈图绘制 1. 将差异表达结果的基因名称转化为id 因为GO和KEGG分析需要用到... -
IPA-蛋白质组、代谢组、转录组分析利器
2020-05-10 13:05:05做为一个超过5年蛋白质组学、代谢组学行业从业以及实验研究者,虽然在几年前早已经听过Ingenuity Pathway Analysis (IPA)的大名,但是对其了解可是真的不多,直到去年接触了IPA分析,才明白原来自己以前做的常规go...