######################The ORA method that assessing candidate gene set relevance to RA.######################

#Draw the intersection graph of the gold standard gene and calculate the corresponding hypergeometric test p value
enrichment_test <- function(pathway_genes, test_genes, background_genes,savepath) {
  # Calculation parameters
  m <- length(background_genes)           # Total number of background genes
  n <- length(pathway_genes)              # Number of pathway genes
  k <- length(test_genes)                  # Number of genes to be tested
  x <- length(intersect(pathway_genes, test_genes)) # Number of intersection genes
  
  # Hypergeometric test
  p_value <- phyper(x - 1, n, m - n, k, lower.tail = FALSE)
  adjusted_p_values <- p.adjust(p_value, method = "BH") 
  ##Draw a Venn diagram
  library(ggvenn)
  
  venndata<- list( `Gold Standard Gene Set`=  pathway_genes, `Your candidate gene` = test_genes)
  set_sizes <- sapply(venndata, length)
  venn.plot <- ggvenn(
    data = venndata,
    fill_color = c("lightblue", "lightpink"), 
    stroke_size = 0.5,
    set_name_size = 5,
    text_size = 4)+ 
    annotate("text", 
              x = 0, y =-1.5, 
              label = paste("P-value(hypergeometric test):",sprintf("%.3e", p_value)), 
              size = 5,fontface = "bold", color = "black")
  ggsave(str_c(savepath,"venneplot.png"),venn.plot,width = 5, height = 5,bg="white")

}

  library(data.table)
  library(tidyverse)
  library(org.Hs.eg.db)
  #Gold standard gene
  all_result <- fread(str_c(path,"\\","RA_gene_gold_matrix.txt")) 
  all_result$entrez_id <- as.character(all_result$entrez_id )
  
  savepath <- path
  #candidate gene sets
  genes <- data.frame(Symbol=sample(all_result$Symbol,50),entrez_id=sample(all_result$entrez_id,50))
 
  all_sig<- all_result %>% dplyr::filter(P.value<=0.05)
  pathway_genes <- all_sig$entrez_id
  gene_sets <-genes$entrez_id
  background_genes <- all_result$entrez_id
  enrichment_test(pathway_genes, gene_sets, background_genes,savepath)


