# Enrichment analysis using ReactomePA
Based on the DE genes found, try to group genes into pathways to be more interpretable. 


```{r}
library(clusterProfiler)
library(org.Hs.eg.db)
library(gtools)
library(viridis)

# Set working directory
setwd("/home/roger/data_immunotherapy/CD4/DEgenes/")

# List all files matching the pattern and sort them naturally
file_list <- mixedsort(list.files(pattern = "^branch[0-5]+\\.csv$"))

# Initialize the overall lists to store Entrez IDs and original gene markers
enrich_gene_list <- list()
gene_marker_list <- list()

# Loop through each file
for (file in file_list) {
  # Read the CSV file
  markers <- read.csv(file, stringsAsFactors = FALSE)
  
  # Convert gene symbols to Entrez IDs
  sig.gene2pre <- bitr(markers$gene, 
                       fromType = "ALIAS", 
                       toType = "ENTREZID", 
                       OrgDb = org.Hs.eg.db)

  # Extract Entrez IDs
  entrezz2pre <- sig.gene2pre$ENTREZID

  # Extract branch name
  branch_name <- gsub("\\.csv", "", file) # Extract branch name from file

  # Add Entrez IDs to the list
  enrich_gene_list[[branch_name]] <- entrezz2pre

  # Add original gene markers to a separate list
  gene_marker_list[[branch_name]] <- markers$gene
}


```





```{r}


reactome <- compareCluster(geneCluster=enrich_gene_list, fun="enrichPathway", pvalueCutoff=0.05)

svg("/home/roger/data_immunotherapy/CD4/reactomeCD4.svg")
# Create the dot plot
p1 <- dotplot(reactome, label_format = 75, size = "count", showCategory = 5, font.size = 8) #, color="p.adjust")

# Define significance threshold
significance_threshold <- 0.05  # Adjust if needed

# Calculate logarithmic breaks
my_log_breaks <- round(log10(significance_threshold)):round(log10(min(reactome@compareClusterResult$p.adjust)))
my_log_breaks <- my_log_breaks[my_log_breaks %% 3 == 0]
my_breaks <- 10^my_log_breaks

# Apply the 'cividis' color scale with fill instead of color
p5 <- p1 + scale_fill_viridis(
  option = 'cividis',
  direction = -1,
  trans = 'log',
  breaks = my_breaks
)
print(p5)
dev.off()


# Convert the compareClusterResult object to a data frame
reactome_df <- as.data.frame(reactome@compareClusterResult)



# Specify the output file path
output_file <- "/home/roger/data_immunotherapy/CD4/reactome_compareClusterResult.csv"

# Save the data frame as a CSV file
write.csv(reactome_df, file = output_file, row.names = FALSE)



```