# Defining the experiment variable (Comment out one experiment at a time)
#experiment <- "240322_Gly"
#experiment <- "240325_Gly"
experiment <- "240331_Gly"

# Read data from an Excel file into two separate data frames: data_lum from the "csvLum" sheet and data_OD from the "csvOD" sheet
data_lum <- read_excel(file.path(paste0("data/glycogen/", experiment, ".xlsx")), sheet = "csvLum")
data_OD <- read_excel(file.path(paste0("data/glycogen/", experiment, ".xlsx")), sheet = "csvOD")
data_norm <- read_excel(file.path(paste0("data/glycogen/", experiment, ".xlsx")), sheet = "csvNorm")

# Filter out rows where "Reject" column is "F" and where "Type" column is not "Empty" in the data_lum data frame, and remove the "Reject" column for both dataframes. Remove species "Cna", as this strain was not included on the manuscript.
data_lum <- data_lum %>% 
  filter(Reject == "F") %>%
  filter(Type != "Empty") %>%
  filter(Species != "Cna") %>%
  select(-Reject)

data_OD <- data_OD %>% 
  filter(Reject == "F") %>%
  filter(Species != "Cna") %>%
  select(-Reject)

data_norm <- data_norm %>% 
  filter(Species != "Cna")

# Split the data_lum data frame into two separate data frames: data_lum_samples containing rows where "Type" is "Sample" and data_lum_standards containing rows where "Type" is "Standard", and remove unnecessary columns
data_lum_samples <- data_lum %>% 
  filter(Type == "Sample")

data_lum_standards <- data_lum %>% 
  filter(Type == "Standard") %>%
  select(-Species, -Strain)

# Convert the "Amy" column in data_lum_standards to a factor and "Glycogen" column to numeric
data_lum_standards$Amy <- as.factor(data_lum_standards$Amy)
data_lum_standards$Glycogen <- as.numeric(data_lum_standards$Glycogen)

# Function to analyze a single timepoint
analyze_glycogen <- function(time, data_lum_standards, data_lum_samples) {
  # Filter the data_lum_standards data frame for rows where "TimePoint" is equal to the specified time
  time_standard <- data_lum_standards %>% 
    filter(TimePoint == time) %>%
    group_by(Glycogen) %>%
    mutate(Difference = ifelse(Amy == "T", Lum - lag(Lum), NA)) %>%
    filter(!is.na(Difference))

  # Fit a linear model (lm) for the standards with "Difference" as the dependent variable and "Glycogen" as the independent variable
  lmStandard <- lm(Difference ~ Glycogen, time_standard) 

  # Extract coefficients from the linear model
  coefficients <- lmStandard$coefficients

  # Filter the data_lum_samples data frame for rows where "TimePoint" is equal to the specified time
  time_samples <- data_lum_samples %>%
    filter(TimePoint == time) %>%
    mutate(Glycogen = (Lum - coefficients[1]) / coefficients[2])

  # Group the data by Strain and Amy, calculate the mean Lum value within each group
  time_samples_diff <- time_samples %>% 
    group_by(Species, Amy) %>%
    summarise(LumAvg = mean(Lum)) %>%
    mutate(Difference = ifelse(Amy == "T", LumAvg - lag(LumAvg), NA)) %>%
    filter(!is.na(Difference)) %>%
    mutate(Glycogen = (Difference - coefficients[1]) / coefficients[2]) %>%
    select(-Amy, -LumAvg) %>%
    mutate(TimePoint = time)

  # Adds individual IDs to the replicas, with and without Amy.
  time_samples <- time_samples %>%
    arrange(Species, Amy)
  time_samples$Replica <- as.factor(rep(1:6, length.out = nrow(time_samples)))

  return(list(time_samples_diff = time_samples_diff))
}

# Loop through different values of time
all_results <- list()

for (time_point in unique(data_lum_standards$TimePoint)) {
  result <- analyze_glycogen(time_point, data_lum_standards, data_lum_samples)
  all_results[[paste0("Time_", time_point)]] <- result

# Combine all time_samples_diff dataframes into a single dataframe
full_Lum_df <- do.call(rbind, lapply(all_results, function(result) result$time_samples_diff))

full_Lum_df <- left_join(full_Lum_df, data_norm) %>%
  mutate(GlycogenNorm = Glycogen/ODNorm)
}

# Convert the 'Species' column to a factor with specified levels
full_Lum_df$Species <- factor(full_Lum_df$Species, levels = c("Cpr", "Cps", "Cac", "Ctu", "Cgl"))

# Convert the 'TimePoint' column to a factor with default level ordering
full_Lum_df$TimePoint <- factor(full_Lum_df$TimePoint)

# Assign the first 'Date' value from 'data_OD' to all rows in 'full_Lum_df'
full_Lum_df$Date <- data_OD$Date[1]

# Write the 'full_Lum_df' dataframe to a CSV file, excluding row names
#write.csv(full_Lum_df, paste0("data/glycogen/Lum/", experiment,"_Lum.csv"), row.names = F)

# Write the 'data_OD' dataframe to a CSV file, excluding row names
#write.csv(data_OD, paste0("data/glycogen/OD/", experiment,"_OD.csv"), row.names = F)

# Define the directory path for OD data
folder_OD <- "data/glycogen/OD"

# Create a list of full file names for all '_OD.csv' files in the directory
file_list_OD <- list.files(folder_OD, pattern = "_OD.csv", full.names = TRUE)

# Read each file into a list of dataframes
df_list_OD <- lapply(file_list_OD, read.csv)

# Combine all dataframes in the list into a single dataframe
df_OD <- bind_rows(df_list_OD)

# Define the directory path for Lum data
folder_lum <- "data/glycogen/Lum"

# Create a list of full file names for all '_Lum.csv' files in the directory
file_list_lum <- list.files(folder_lum, pattern = "_Lum.csv", full.names = TRUE)

# Read each file into a list of dataframes
df_list_lum <- lapply(file_list_lum, read.csv)

# Combine all dataframes in the list into a single dataframe
df_lum <- bind_rows(df_list_lum)

# Convert the 'Species' column in all average data frames to a factor with specified levels
df_OD$Species <- factor(df_OD$Species, levels = c("Cpr", "Cps", "Cac", "Ctu", "Cgl"))
df_lum$Species <- factor(df_lum$Species, levels = c("Cpr", "Cps", "Cac", "Ctu", "Cgl"))

# Convert the 'TimePoint' column in Lum data frame to a factor with default level ordering
df_lum$TimePoint <- factor(df_lum$TimePoint)

# Define colors
species_colors <- c("#FF8C00", "#FF0000", "#6A329F", "#3B54A4", "grey50")