Dynamic individual doses

This example calculates dynamic individual doses using the lixoftConnectors to automate Simulx to simulate the model for a specific interval, check outcomes according to a specified threshold, then adjust the dose for the next interval. The code and Simulx project were created with MonolixSuite version 2024.
In order to run the code:
Download and unzip the starting Simulx project: neutropenia_dynamic_dosing.zip
Put the following R script in the same folder
Set the working directory to the R script location
R
#------------------------------------------------------------
# example of dynamic individual dosing
#------------------------------------------------------------

library(lixoftConnectors)
initializeLixoftConnectors(software = "simulx")

# load neutrophil project
loadProject("neutropenia_dynamic_dosing.smlx")
project_name = "neutropenia_dynamic_dosing_simulations.smlx"
saveProject(project_name)

#------------------------------------------------------------
# scenario setup
#------------------------------------------------------------

# how many individuals to simulate 
nb_indivs = 20

# dose information
treatment_days = 300
starting_dose = 0.001
dose_increment = 0.0002
dose_ANC_threshold = 0.5

# time step for updating dose
interval_start = 0
interval_duration = 14

# define treatment as an external file that we can easily update per id
treatment_df = data.frame(id = rep(1:nb_indivs, each = treatment_days),
                          time = rep(1:treatment_days, times = nb_indivs),
                          amount = starting_dose)
treatment_file = "treatment.csv" 
write.csv(treatment_df, file = treatment_file, quote = FALSE, row.names = FALSE)

#------------------------------------------------------------
# simulx setup
#------------------------------------------------------------

# rename group to dynamic_dose and update accordingly
renameGroup("standard_dose", "dynamic_dose")

# remove full output element and add an output element to be updated to loop by the specified interval
names(getOutputElements()) # regularANC
defineOutputElement("interval_step", 
                    element = list(
                      data = data.frame(start = interval_start, interval = 1, final = interval_start + interval_duration - 1), 
                      output = "ANC"))
removeGroupElement("dynamic_dose", "regularANC")
setGroupElement("dynamic_dose", "interval_step")

# define treatment as external so easily adaptable
defineTreatmentElement(name = "Chemo_dynamic", element = list(data = treatment_file))

# set treatment
removeGroupElement("dynamic_dose", "Chemo_daily_dose")
setGroupElement("dynamic_dose", "Chemo_dynamic")

setGroupSize("dynamic_dose", nb_indivs)

# verify groups
getGroups()

# set outcome element to test if ANC < threshold during interval
defineOutcome("below_threshold", 
              element = list(
                output = "interval_step",
                processing = list(operator = "durationBelow", type = "percentTime", value = dose_ANC_threshold)
              ))
defineEndpoint("mean_below_threshold", element = list(outcome = "below_threshold"))

setScenario(c(simulation = TRUE, endpoints = TRUE))

#------------------------------------------------------------
# run dynamic dosing
#------------------------------------------------------------

# loop over the intervals, updating dose as necessary
while (interval_start < treatment_days)
{
  # run simulation for interval
  runScenario()

  # check conditions and update dose if needed
  outcome = getEndpointsResults()$outcomes$below_threshold
  
  for (id in outcome$id)
  {
    # if below threshold during interval
    if (outcome$below_threshold[outcome$id == id] > 0)
    {
      # update future doses
      current_dose = treatment_df$amount[treatment_df$id == id & treatment_df$time == interval_start]
      new_dose = max(0, current_dose - dose_increment)
      treatment_df$amount[treatment_df$id == id & treatment_df$time >= interval_start + interval_duration] = new_dose
    }
  }
  
  # set treatment and interval for next step
  interval_start = interval_start + interval_duration
  
  write.csv(treatment_df, file = treatment_file, quote = FALSE, row.names = FALSE)
  defineTreatmentElement(name = "Chemo_dynamic", element = list(data = treatment_file))

  defineOutputElement("interval_step", 
                      element = list(
                        data = data.frame(start = interval_start, interval = 1, final = interval_start + interval_duration - 1), 
                        output = "ANC"))
  
}

# save project if you want to look at it in Simulx
saveProject(project_name)


#------------------------------------------------------------
# rerun over the entire time comparing to a non-dynamic dose
#------------------------------------------------------------

removeGroupElement("dynamic_dose", "interval_step")
setGroupElement("dynamic_dose", "regularANC")
defineTreatmentElement(name = "Chemo_dynamic", element = list(data = treatment_file))

addGroup("standard_dose")
removeGroupElement("standard_dose", "Chemo_dynamic")
setGroupElement("standard_dose", "Chemo_daily_dose")

getGroups()

setSameIndividualsAmongGroups(TRUE)

saveProject(project_name)
runScenario()
results = getSimulationResults()$res$ANC

# verify individual parameters
getSimulationResults()$individualParameters

#------------------------------------------------------------
# plots
#------------------------------------------------------------

library(ggplot2)
results$original_id = factor((results$id %% nb_indivs) + 1)

ggplot(data = results, mapping = aes(x = time, y = ANC, color = original_id)) +
  geom_hline(yintercept = dose_ANC_threshold, color = "red", linewidth = 1) +
  geom_vline(xintercept = treatment_days, linetype = 2) +
  geom_line(aes(group = id)) +
  facet_wrap(vars(group)) +
  theme_classic()

ggplot(data = treatment_df,  mapping = aes(x = time, y = amount, color = factor(id))) +
  ylim(0, 1.1 * starting_dose) +
  geom_line(aes(group = id), show.legend = FALSE) +
  facet_wrap(vars(id)) +
  theme_classic()