Model 1A GA optmilasation.Rmd

```{r}
library(tidyverse)
library(GillespieSSA2)
library(GA)
library(doParallel)

# Define initial states
ini_state <- c(
  DD = 10, 
  p53 = 10, 
  p53P = 17
)

# Define the objective function
objective_function <- function(params) {
  param_names <- c("CP_DD", "CD_DD", "CP_p53", "CD_p53", "phos", "dephos", "CD_p53P")
  names(params) <- param_names

  reactions1 <- list(
    reaction(~params["CP_DD"], c(DD = +1), name = "Prod_DD"),
    reaction(~params["CD_DD"] * p53P, c(DD = -1), name = "Decay_DD"),
    reaction(~params["CP_p53"] * DD, c(p53 = +1), name = "Prod_p53"),
    reaction(~params["CD_p53"] * p53, c(p53 = -1), name = "Decay_p53"),
    reaction(~params["phos"] * p53 * DD, c(p53 = -1, p53P = +1), name = "Phos_p53"),
    reaction(~params["dephos"] * p53P, c(p53P = -1, p53 = +1), name = "Dephos_p53P"),
    reaction(~params["CD_p53P"] * p53P, c(p53P = -1), name = "Decay_p53P")
  )

  out <- tryCatch({
    ssa(
      initial_state = ini_state,
      reactions = reactions1,
      params = params,
      method = ssa_exact(),
      final_time = 50,
      census_interval = .01,
      verbose = FALSE
    )
  }, error = function(e) {
    return(NULL)
  })
  
  if (is.null(out)) return(Inf)
  
  p53P_data <- out %>% filter(variable == "p53P")
  oscillation_metric <- sd(p53P_data$value)
  
  return(-oscillation_metric)
}
```

```{r}
# Define parameter bounds
lower_bounds <- c(0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001)
upper_bounds <- c(10, 10, 10, 10, 10, 10, 10, 10)

# Run the Genetic Algorithm optimization
ga_result <- ga(
  type = "real-valued",
  fitness = objective_function,
  lower = lower_bounds,
  upper = upper_bounds,
  popSize = 100,
  maxiter = 100,
  run = 50,
  keepBest = TRUE,
  parallel = TRUE, # Parallel for efficiency,
  seed = 123  # Set a seed for reproducibility
)

# Extract the best parameter set
best_params <- ga_result@solution
names(best_params) <- c("CP_DD", "CD_DD", "CP_p53", "CD_p53", "phos", "dephos", "CD_p53P")

# Print the best parameters
print(best_params)
```