DAVID_functional_annotation.R

In Emask7/ZIKVtranscriptomics: ZIKV Transcriptomics Analyses

library(limma)
library(openxlsx)
library(org.Hs.eg.db)
library(RDAVIDWebService)
library(stringi)
library(dplyr)
library(GOplot)


multi_DAVID <- function(DE_data, min_lfc, list_name, FA_cats, split_at){
  DE_data <- subset(DE_data, !grepl("ENSPANG", rownames(DE_data)))

  t3v0_res <- run_DAVID(
    david, DE_data[, 1:2], min_lfc,
    stri_join(list_name, "_t3v0"), FA_cats, split_at
  )
  t15v0_res <- run_DAVID(
    david, DE_data[, 3:4], min_lfc,
    stri_join(list_name, "_t15v0"), FA_cats, split_at
  )
  t15v3_res <- run_DAVID(
    david, DE_data[, 5:6], min_lfc,
    stri_join(list_name, "_t15v3"), FA_cats, split_at
  )

  if (is.null(t3v0_res$DAVID_FA)) {
    t3v0_res$DAVID_FA <- data.frame(NA, NA, NA, NA, NA, NA, NA, NA, NA)
  }
  if (is.null(t15v0_res$DAVID_FA)) {
    t15v0_res$DAVID_FA <- data.frame(NA, NA, NA, NA, NA, NA, NA, NA, NA)
  }
  if (is.null(t15v3_res$DAVID_FA)) {
    t15v3_res$DAVID_FA <- data.frame(NA, NA, NA, NA, NA, NA, NA, NA, NA)
  }

  colnames(t3v0_res$DAVID_FA)  <- c(
    "Annotation_Category", "ID", "Term", "Count_t3v0", "%_of_input_t3v0",
    "Genes_t3v0", "Fold.Enrichment_t3v0", "FDR_t3v0", "Z-score_t3v0"
  )
  colnames(t15v0_res$DAVID_FA)  <- c(
    "Annotation_Category", "ID", "Term", "Count_t15v0", "%_of_input_t15v0",
    "Genes_t15v0", "Fold.Enrichment_t15v0", "FDR_t15v0", "Z-score_t15v0"
  )
  colnames(t15v3_res$DAVID_FA)  <- c(
    "Annotation_Category", "ID", "Term", "Count_t15v3", "%_of_input_t15v3",
    "Genes_t15v3", "Fold.Enrichment_t15v3", "FDR_t15v3", "Z-score_t15v3"
  )

  all_res <- full_join(t3v0_res$DAVID_FA, t15v0_res$DAVID_FA)
  all_res <- full_join(all_res, t15v3_res$DAVID_FA)
  all_res <- subset(all_res, !is.na(all_res[, 1]))

  list(t3v0 = t3v0_res, t15v0 = t15v0_res, t15v3 = t15v3_res, all = all_res)
}



david <- DAVIDWebService(
  email = "emask@txbiomed.org",
  url = "https://david.ncifcrf.gov/webservice/services/DAVIDWebService.DAVIDWebServiceHttpSoap12Endpoint/"
)
is.connected(david)
show(david)
getIdTypes(david)

DAVID_results_human <- list(
  GO = multi_DAVID(
    DESeq2_LRT_res_human, 0.8, "DESeq2_LRT_human",
    c("GOTERM_BP_FAT", "GOTERM_MF_FAT"), "~"
  ),
  KEGG = multi_DAVID(
    DESeq2_LRT_res_human, 0.8, "DESeq2_LRT_human", c("KEGG_PATHWAY"), ":"
  ),
  Reactome = multi_DAVID(
    DESeq2_LRT_res_human, 0.8, "DESeq2_LRT_human", c("REACTOME_PATHWAY"), ":"
  )
)

wb <- createWorkbook("ZIKV DAVID Functional Annotation - Human Alignment.xlsx")
addWorksheet(wb, "GO Analysis")
addWorksheet(wb, "KEGG Pathways")
addWorksheet(wb, "Reactome Pathways")
writeData(wb, "GO Analysis", DAVID_results_human$GO$all)
writeData(wb, "KEGG Pathways", DAVID_results_human$KEGG$all)
writeData(wb, "Reactome Pathways", DAVID_results_human$Reactome$all)
saveWorkbook(
  wb, "ZIKV DAVID Functional Annotation - Human Alignment.xlsx", overwrite = TRUE
)



DAVID_results_baboon <- list(
  GO = multi_DAVID(
    DESeq2_LRT_res_baboon, 0.8, "DESeq2_LRT_baboon_08",
    c("GOTERM_BP_FAT", "GOTERM_MF_FAT"), "~"
  # ),
  # KEGG = multi_DAVID(
  #   DESeq2_LRT_res_baboon, 0.8, "DESeq2_LRT_baboon_08", c("KEGG_PATHWAY"), ":"
  ),
  Reactome = multi_DAVID(
    DESeq2_LRT_res_baboon, 0.8, "DESeq2_LRT_baboon_08", c("REACTOME_PATHWAY"), ":"
  )
)

DAVID_results_baboon <- list(
  GO = multi_DAVID(
    DESeq2_LRT_res_baboon, 1, "DESeq2_LRT_baboon",
    c("GOTERM_BP_FAT", "GOTERM_MF_FAT"), "~"
  ),
  KEGG = multi_DAVID(
    DESeq2_LRT_res_baboon, 1, "DESeq2_LRT_baboon", c("KEGG_PATHWAY"), ":"
  ),
  Reactome = multi_DAVID(
    DESeq2_LRT_res_baboon, 1, "DESeq2_LRT_baboon", c("REACTOME_PATHWAY"), ":"
  )
)


# wb <- createWorkbook("ZIKV DAVID Functional Annotation - Baboon Alignment.xlsx")
wb <- createWorkbook("ZIKV DAVID Functional Annotation - Panubis1 Alignment.xlsx")

addWorksheet(wb, "GO Analysis")
addWorksheet(wb, "KEGG Pathways")
addWorksheet(wb, "Reactome Pathways")
writeData(wb, "GO Analysis", DAVID_results_baboon$GO$all)
writeData(wb, "KEGG Pathways", DAVID_results_baboon$KEGG$all)
writeData(wb, "Reactome Pathways", DAVID_results_baboon$Reactome$all)

# saveWorkbook(
#   wb, "ZIKV DAVID Functional Annotation - Baboon Alignment.xlsx", overwrite = TRUE
# )
saveWorkbook(
  wb, "ZIKV DAVID Functional Annotation - Panubis1 Alignment.xlsx", overwrite = TRUE
)