library(clusterProfiler)
library(org.Dm.eg.db)
library(dplyr)
library(ggplot2)
library(VennDiagram)

plotTheme = theme(
    axis.title = element_text(family = "Arial", face = "bold", size = 16),
    legend.title = element_text(family = "Arial", face = "bold", size = 16), 
    text = element_text(family = "Arial",size = 14))

载入需要的程辑包：grid

载入需要的程辑包：futile.logger

anaGO = function(geneids,prefix){
    res = enrichGO(geneids,OrgDb = 'org.Dm.eg.db',keyType = 'ENTREZID',ont = 'ALL',readable = TRUE,maxGSSize = 300)
    if(!is.null(res)){
        write.csv(res@result,file = paste0(prefix,'_GO.csv'))
        resTable = res@result[order(res@result$p.adjust),]
        resTable = resTable %>% group_by(ONTOLOGY) %>% top_n(-p.adjust,n = 10)
        p = ggplot(data = resTable,
                   aes(
                       x = reorder(Description,-p.adjust),
                       y = -log10(p.adjust),
                       fill = ONTOLOGY,
                       label = reorder(Count,-p.adjust)
                    )
            )+
            geom_bar(stat = 'identity')+
            geom_text(stat = 'identity',hjust = 1.2)+
            coord_flip()+
            facet_grid(ONTOLOGY~.,scales = 'free',space = 'free')+
            labs(x = '',fill = '')+
            theme_classic()+
            scale_fill_manual(
                breaks = c('BP','CC','MF'),
                labels = c('BP','CC','MF'),
                values = c('#ffc64b','#66CC33','#33ccff')
            )+
            theme(strip.text = element_blank(),
                  axis.line.y = element_blank(),
                  axis.ticks.y = element_blank(),
                  axis.title.x = element_text(size = 14),
                  text = element_text(size = 12),
                  axis.text = element_text(size = 12))+plotTheme
        ggsave(filename = paste0(prefix,'_GO.png'),plot = p,dpi = 300,width = 10,height = 8,bg = 'white')
        p
    }else{
        print(paste0('No GO terms enriched in ',prefix,'!!!!'))
    }
}

pa14DEG = read.csv('pa14_DEG.csv',header = TRUE,row.names = 1)
pa14DEG

down = pa14DEG$names[pa14DEG$status=='down']
geneTmp = bitr(down,fromType = 'SYMBOL',toType = 'ENTREZID',OrgDb = 'org.Dm.eg.db',drop = FALSE)
numTmp = sum(is.na(geneTmp[,2]))
print(paste0(numTmp,' genes failed to map.'))
geneTmp[is.na(geneTmp[,2]),]

'select()' returned 1:1 mapping between keys and columns

Warning message in bitr(down, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = "org.Dm.eg.db", :
“5.61% of input gene IDs are fail to map...”


[1] "72 genes failed to map."

anaGO(na.omit(geneTmp[,2]),prefix = 'pa14_down')

agingDEG = read.csv('aging_DEG.csv',header = TRUE,row.names = 1)
agingDown = agingDEG$names[agingDEG$status=='down']

venn.diagram(x=list(agingDown,down),category.names = c("aging Down","PA14 Down"),scaled = FALSE,alpha = 0.3,filename = "aging_PA14_Down_venn.png",imagetype = "png",resolution = "720")

1

tmp = intersect(agingDown,down)
geneTmp = bitr(tmp,fromType = 'SYMBOL',toType = 'ENTREZID',OrgDb = 'org.Dm.eg.db',drop = FALSE)
numTmp = sum(is.na(geneTmp[,2]))
print(paste0(numTmp,' genes failed to map.'))
geneTmp[is.na(geneTmp[,2]),]

'select()' returned 1:1 mapping between keys and columns

Warning message in bitr(tmp, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = "org.Dm.eg.db", :
“7.76% of input gene IDs are fail to map...”


[1] "45 genes failed to map."

anaGO(na.omit(geneTmp[,2]),prefix = 'aging_pa14_intersect_down')

tmp = intersect(agingDEG$names[agingDEG$status=='up'],pa14DEG$names[pa14DEG$status=='up'])
tmp[grep('^CG',tmp,invert = TRUE)]

1. 'Nazo'
2. 'Tsp42Ee'
3. 'Jupiter'
4. 'Nrt'
5. 'Trx-2'
6. 'hdly'
7. 'pnt'
8. 'Lmpt'
9. 'GstD2'
10. 'PGRP-LB'
11. 'mnb'
12. 'wun2'
13. 'JhI-26'
14. 'krimp'
15. 'rst'
16. 'msi'
17. 'h'
18. 'Lpin'
19. 'lncRNA:CR46064'
20. 'sna'
21. 'Mec2'
22. 'Ldh'
23. 'B4'
24. 'Fkbp14'
25. 'cnn'
26. 'GstS1'
27. 'Pde6'
28. 'wb'
29. 'inaE'
30. 'C1GalTA'
31. 'Myo61F'
32. 'Cpr51A'
33. 'AcCoAS'
34. 'Swip-1'
35. 'l(3)02640'
36. 'Shrm'
37. 'Gk2'
38. 'Pen'
39. 'Dl'
40. 'del'
41. 'polo'
42. 'Rgl'
43. 'kuk'
44. 'Rcd2'
45. 'chn'
46. 'dve'
47. 'dpr21'
48. 'Ant2'
49. 'arr'
50. 'fry'
51. 'dpn'
52. 'Tsp2A'
53. 'Cks30A'
54. 'Btk29A'
55. 'GstE8'
56. 'Marf1'
57. 'cu'
58. 'pdgy'
59. 'lncRNA:CR43144'
60. 'mbc'
61. 'sog'
62. 'bves'
63. 'pch2'
64. 'Nnf1a'
65. 'CdGAPr'
66. 'Kr-h1'
67. 'msd1'
68. 'MFS3'
69. 'GLS'
70. 'smp-30'
71. 'Wnt5'
72. 'kirre'
73. 'E(spl)m7-HLH'
74. 'Ntan1'
75. 'ImpL2'
76. 'nahoda'
77. 'GstD5'
78. 'lncRNA:roX2'
79. 'Sema2b'

venn.diagram(x=list(agingDEG$names[agingDEG$status=='up'],pa14DEG$names[pa14DEG$status=='up']),category.names = c("aging Up","PA14 Up"),scaled = FALSE,alpha = 0.3,filename = "aging_PA14_Up_venn.png",imagetype = "png",resolution = "720")

1

geneTmp = bitr(tmp,fromType = 'SYMBOL',toType = 'ENTREZID',OrgDb = 'org.Dm.eg.db',drop = FALSE)
numTmp = sum(is.na(geneTmp[,2]))
print(paste0(numTmp,' genes failed to map.'))
geneTmp[is.na(geneTmp[,2]),]
anaGO(na.omit(geneTmp[,2]),prefix = 'aging_pa14_intersect_up')

'select()' returned 1:1 mapping between keys and columns

Warning message in bitr(tmp, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = "org.Dm.eg.db", :
“4.79% of input gene IDs are fail to map...”


[1] "7 genes failed to map."