R/tcgaEnrichment.R

In dming1024/TCGApackages0226: Analysis of TCGA data

#' Carry out the Enrichment analysis on the KEGG and GO.
#'
#' use the tcgaEnrichment() fuction to get the enrichment results of DEGs.
#' @rawNamespace import(AnnotationDbi,except=select)
#' @import clusterProfiler org.Hs.eg.db patchwork Cairo stringr enrichplot
#' @param dataPath the direction of download Data
#' @param pvalue P Vlaue threshold for the significant results
#' @param fw the width of PDF files
#' @param fh the height of PDF files
#' @export
tcgaEnrichment<-function(dataPath,pvalue,fw=15,fh=10){
  setwd(dataPath)
  message("对差异表达基因进行在线富集分析...",appendLF = F)
  pvalue=pvalue
  load(file=".//output//step07.RData")
  DEGs_entrez_id <- AnnotationDbi::mapIds(
    x = org.Hs.eg.db::org.Hs.eg.db,
    column = 'ENTREZID',
    keys = as.character(TCGA_result_genesymbol$hgnc_symbol),
    keytype = 'SYMBOL'
  )

  enrich_go_BP <- clusterProfiler::enrichGO(
    gene = DEGs_entrez_id,
    OrgDb = 'org.Hs.eg.db',
    ont = 'BP',
    pvalueCutoff = pvalue
    )

  enrich_go_CC <- clusterProfiler::enrichGO(
    gene = DEGs_entrez_id,
    OrgDb = 'org.Hs.eg.db',
    ont = 'CC',
    pvalueCutoff = pvalue
  )
  enrich_go_MF <- clusterProfiler::enrichGO(
    gene = DEGs_entrez_id,
    OrgDb = 'org.Hs.eg.db',
    ont = 'MF',
    pvalueCutoff = pvalue
  )
  enrich_go_KEGG <- clusterProfiler::enrichKEGG(
    gene = DEGs_entrez_id,
    organism = 'hsa',
    keyType = 'kegg',
    pvalueCutoff = pvalue
  )

  Cairo::CairoPDF(".//output//step08_goEnrichment02.pdf",family="Arial",width = fw,height = fh)
  #pdf(file=".//output//step08_goEnrichment02.pdf",width = 12,height = 10)
  p1=enrichplot::dotplot(enrich_go_BP,title='Biological Process of DEGs', showCategory = 8,orderBy="GeneRatio")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))+
    ggplot2::xlim(NA,max(unlist(lapply(enrich_go_BP@result$GeneRatio, function(x){round(eval(parse(text=x)),digits = 4)})))*1.05)
  #print(p1)
  p2=clusterProfiler::dotplot(enrich_go_CC,title='Cell component of DEGs', showCategory = 8,orderBy="GeneRatio")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))+
    ggplot2::xlim(NA,max(unlist(lapply(enrich_go_CC@result$GeneRatio, function(x){round(eval(parse(text=x)),digits = 4)})))*1.05)
  #print(p1+p2)
  p3=enrichplot::dotplot(enrich_go_MF,title='Molecular Function of DEGs', showCategory = 8,orderBy="GeneRatio")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))+
    ggplot2::xlim(NA,max(unlist(lapply(enrich_go_MF@result$GeneRatio, function(x){round(eval(parse(text=x)),digits = 4)})))*1.05)
  #print(p3)
  p4=enrichplot::dotplot(enrich_go_KEGG,title='KEGG pathway of DEGs', showCategory = 8,orderBy="GeneRatio")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))+
    ggplot2::xlim(NA,max(unlist(lapply(enrich_go_KEGG@result$GeneRatio, function(x){round(eval(parse(text=x)),digits = 4)})))*1.05)
  p=patchwork::wrap_plots(list(p1,p2,p3,p4),ncol=2)
  print(p)
  dev.off()



  Cairo::CairoPDF(file = ".//output//step08_goEnrichment01.pdf",family="Arial",width = fw,height = fh)
  p1=barplot(enrich_go_BP,title='Biological Process of DEGs',orderBy="Count")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))
  p2=barplot(enrich_go_CC,title='Cell component of DEGs',orderBy="Count")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))
  p3=barplot(enrich_go_MF,title='Molecular Function of DEGs',orderBy="Count")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))
  p4=barplot(enrich_go_KEGG,title='KEGG pathway of DEGs',orderBy="Count")+
    ggplot2::scale_y_discrete(labels=function(y)stringr::str_wrap(y,width = 30))
  p=patchwork::wrap_plots(list(p1,p2,p3,p4),ncol=2)
  print(p)
  #print((p1+p2)/(p3+p4))
  dev.off()


  #输出富集分析结果
  write.csv(enrich_go_BP@result,".//output//BiologicProgress.csv")
  write.csv(enrich_go_CC@result,".//output//CellComponent.csv")
  write.csv(enrich_go_MF@result,".//output//MolecularFunction.csv")
  write.csv(enrich_go_KEGG@result,".//output//KEGGResults.csv")
  message("Completed!")
}