# Function that import all data (take a lot of time)import_raw_data_microarray <-function(path_i, path_info_range_i, sheet_i, begin_end_line){# Import info #### column_global_info <-c("A", "I") range_global_info =paste0(column_global_info[1],begin_end_line[1],":", column_global_info[2],begin_end_line[2])cat("Read global info excel sheet (", range_global_info, ")\n") df_info =read_excel(path_i, sheet = sheet_i, range = range_global_info, col_names = T) %>%as.data.frame() %>% dplyr::rename(id_probe="N°",other_description1 =...7,other_description2 =...9) %>% dplyr::select(-...6)# import the column where the data is located #### compile_raw_data <-as.data.frame(matrix(data =NA,nrow =as.numeric(begin_end_line[2]) -as.numeric(begin_end_line[1]) ,ncol =0)) df_info_range =read_excel(path_info_range_i, col_names = T) n <-nrow(df_info_range) pb <- progress_bar$new(format ="[:bar] :percent | Elapsed time :elapsed | Estimated time :eta",total = n,clear =FALSE,width =80 ) start_time <-Sys.time() df_info_sample_clean_compile <-as.data.frame(matrix(data =NA, nrow =0, ncol =13)) colnames(df_info_sample_clean_compile) <-c("date", "rep", "project_nb", "project", "sample_type", "num_combination", "sample_name", "color", "sample_num", "simplify_condition", "genotype", "sulfur_condition", "sample_id")for (i in1:n) { column_red_green <-c(df_info_range$red_column[i], df_info_range$green_column[i]) num_combination <- df_info_range$num_combination[i] begin_end_line_info <-c("3", "12") range_info <-paste0(column_red_green[1], begin_end_line_info[1], ":", column_red_green[2], begin_end_line_info[2])cat("Read info of the sample", range_info, "\n") df_info_sample <-read_excel(path_i, sheet = sheet_i, range = range_info, col_names =FALSE) %>%as.data.frame() %>% dplyr::rename(col1 = ...1,col2 = ...2 ) df_info_sample_clean <-as.data.frame(matrix(data =NA, nrow =2, ncol =13)) colnames(df_info_sample_clean) <-c("date", "rep", "project_nb", "project", "sample_type", "num_combination", "sample_name", "color", "sample_num", "simplify_condition", "genotype", "sulfur_condition", "sample_id") df_info_sample_clean$date <-c(df_info_sample$col1[1], df_info_sample$col2[1]) df_info_sample_clean$rep <-c(df_info_sample$col1[2], df_info_sample$col2[2]) df_info_sample_clean$project_nb <-c(df_info_sample$col1[5], df_info_sample$col2[5]) df_info_sample_clean$project <-c(df_info_sample$col1[6], df_info_sample$col2[6]) df_info_sample_clean$num_combination <- num_combination df_info_sample_clean$sample_type <-c(df_info_sample$col1[7], df_info_sample$col2[7]) df_info_sample_clean$sample_name <-c(df_info_sample$col1[8], df_info_sample$col2[9]) df_info_sample_clean$color <-c(df_info_sample$col1[10], df_info_sample$col2[10]) df_info_sample_clean <- df_info_sample_clean %>% dplyr::mutate(sample_num =str_extract(sample_name, "^[0-9_.]+"),sample_num =str_replace_all(sample_num, c("\\_"="\\.")), simplify_condition =str_extract(sample_name, "WT[+-]|Mut[+-]"),simplify_condition =str_replace_all(simplify_condition, c("\\+"="_SS", "\\-"="_SD")),sulfur_condition =str_sub(simplify_condition, -2, -1),genotype =ifelse(sample_num %in%c(1.1, 2.2, 3.3, 4.4, 17.1, 18.2, 19.3, 20.4), "WT2", sample_num),# Add genotypegenotype =ifelse(sample_num %in%c(5.1, 6.2, 7.3, 8.4, 21.1, 22.2, 23.3, 24.4), "E568K", genotype),# Add genotypegenotype =ifelse(sample_num %in%c(9.5, 10.6, 11.7, 12.8, 25.5, 26.6, 27.7, 28.8), "WT1", genotype),# Add genotypegenotype =ifelse(sample_num %in%c(13.5, 14.6, 15.7, 16.8, 29.5, 30.6, 31.7, 32.8), "W78*", genotype),# Add genotypesample_id =paste(sep ="_", simplify_condition, sample_num, rep, color, num_combination) ) range_i <-paste0(column_red_green[1], begin_end_line[1], ":", column_red_green[2], begin_end_line[2])# Importer les données rouges/vertes pour les deux échantillonscat("Read raw data range ",range_i, " \n") df_x <-read_excel(path_i, sheet = sheet_i, range = range_i, col_names =TRUE) %>%as.data.frame()# Vérification des noms de colonnesif (colnames(df_x)[1] !="Red"&colnames(df_x)[2] !="Green"&length(colnames(df_x))!=2) {stop("Column name needs to be Red or Green") }colnames(df_x)[1] = df_info_sample_clean$sample_id[1]colnames(df_x)[2] = df_info_sample_clean$sample_id[2] compile_raw_data =cbind(compile_raw_data, df_x) df_info_sample_clean_compile =rbind(df_info_sample_clean_compile, df_info_sample_clean) pb$tick() elapsed_time <-Sys.time() - start_time } global_df <-cbind(df_info , compile_raw_data)return(list(global_df = global_df,df_info_sample_clean_compile = df_info_sample_clean_compile ))}
19.3 GO term microarray
Code
# Corentin Maslard# Date: 20250220# The aim of this scipt is to automatically perform GOs on different groups, and then to compare the different groups with each other. # pkg library(org.Psativum1c.eg.db) # if not working install itlibrary(clusterProfiler)library(GO.db)library(AnnotationDbi)library(dplyr)library(progressr)# functioncat_col <-function(text, color) { colors <-c("red", "green", "yellow", "blue", "magenta", "cyan", "white") color_code <-switch(tolower(color),"red"="31","green"="32","yellow"="33","blue"="34","magenta"="35","cyan"="36","white"="37")# Vérifier si la couleur est valideif (color_code =="") {cat("Couleur non valide.")return(invisible()) }# Afficher le texte avec la couleur spécifiéecat(paste0("\033[", color_code, "m", text, "\033[0m"))}GO_on_different_group <-function(functional_roles ="BP", #MF or CC group_info, # here is a dataframe with in ID (the name of the genes in database) and group (the group (cluster or category)). group, # column of the different group ID, # column of the DF with the name of the different genesminGSSize_i =10, maxGSSize_i =500,pvalueCutoff_i =0.05,pAdjustMethod_i ="BH",top =FALSE){ #keytypes(org.Psativum1c.eg.db) # verification of all inputsif( !functional_roles %in%c("BP", "CC", "MF")){stop("You need to add BP, CC or MF") }# add column in the data frame corresponding to ID and group group_info = group_info %>%mutate(group =!!sym(group), ID =!!sym(ID)) %>% dplyr::select(group, ID)cat_col("- Initialisation of the database \n", "green") all_genes <-keys(x=org.Psativum1c.eg.db, keytype ="GID") all_GO <- AnnotationDbi::select(x = org.Psativum1c.eg.db, keys = all_genes, columns =c("GO","ONTOLOGY") , keytype ="GID")# Select Biological Processes all_GO_type <- all_GO[all_GO$ONTOLOGY == functional_roles,] %>%na.omit() %>% dplyr::select("GO", "GID") %>% dplyr::rename(go_id = GO, gene_id = GID)cat_col("- Performing Gene Ontology enrichment analysis for each group \n", "green") num_groups <-length(unique(group_info$group)) # Total number of iterations pb <-txtProgressBar(min =0, max = num_groups, style =3) # Create progress bar all_results <-list() group_list <-unique(group_info$group) # Store unique groupsfor (i inseq_along(group_list)) { group <- group_list[i] # Get the group namesetTxtProgressBar(pb, i) # Update progress bar# Extract the gene list for the current group gene_list <- group_info %>%filter(group ==!!group) %>%pull(ID)# Perform enrichment analysis enrichment_results <-enricher(gene = gene_list, TERM2GENE = all_GO_type, pvalueCutoff = pvalueCutoff_i,pAdjustMethod = pAdjustMethod_i,minGSSize = minGSSize_i,maxGSSize = maxGSSize_i)# Store results if enrichment is not emptyif (!is.null(enrichment_results) &&nrow(as.data.frame(enrichment_results)) >0) { enrichment_df <-as.data.frame(enrichment_results) enrichment_df$group <- group # Add group ID all_results[[as.character(group)]] <- enrichment_df } }close(pb) # Close progress bar# Combine all cluster results into a single data frame combined_results <-bind_rows(all_results)# Get the GO description go_terms <-keys(GO.db) go_descriptions <-Term(go_terms) go_table <-data.frame(GO_Term = go_terms, Description_GO = go_descriptions)names(go_table)[1] <-"ID" combined_results <-left_join(combined_results, go_table, by ="ID")if(top ==FALSE){cat_col("- Give all go terms \n", "green")return(combined_results) }elseif(top >0){ combined_results.top <- combined_results %>% dplyr::group_by(group) %>%top_n(-top, wt = p.adjust) %>%arrange(group, p.adjust)# Define the list of selected significant GO terms Signif_GO <-unique(combined_results.top$Description)# Filter results based on this list group_selected_GO <-filter(combined_results, ID %in% Signif_GO)# Ensure GeneRatio and BgRatio are character columns group_selected_GO_2 <- group_selected_GO#str(group_selected_GO_2) group_selected_GO_2$GeneRatio <-as.character(group_selected_GO_2$GeneRatio) group_selected_GO_2$BgRatio <-as.character(group_selected_GO_2$BgRatio)# Calculate the fold enrichment group_selected_GO_2$FoldEnrichment <-with(group_selected_GO_2, (sapply(strsplit(GeneRatio, "/"), function(x) as.numeric(x[1]) /as.numeric(x[2]))) / (sapply(strsplit(BgRatio, "/"), function(x) as.numeric(x[1]) /as.numeric(x[2]))))# Calculate the -log10(pvalue) group_selected_GO_2$'|-log10(Pval)|'<--(log10(group_selected_GO_2$p.adjust))# Ensure that `Term` column has factor levels including all possible terms group_selected_GO_2$Description_GO <-factor(group_selected_GO_2$Description_GO, levels =unique(group_selected_GO_2$Description_GO))# Now, create a data frame with all combinations of `Term` and `filename` all_combinations <-expand.grid(Description_GO =levels(group_selected_GO_2$Description_GO), group =unique(group_info$group))# Merge this with your original data to fill missing combinations group_selected_GO_filled <-merge(all_combinations, group_selected_GO_2, by =c("Description_GO", "group"), all.x =TRUE) return(group_selected_GO_filled) }}# test = GO_on_different_group(functional_roles = "CC", # group_info = cluster_info,# group = "cluster", # ID = "ID", # top = 10# )