Population simulations

Population simulations can be easily performed in R by combining the simulation loaded from a *.pkml file with the population information created in PK-Sim and exported to CSV format (for details, please refer to OSPS online documentation) or created directly in R (see Creating populations).

Loading population file

The method loadPopulation creates an object of the Population class that can be passed to the runSimulations() method (see Running simulations and retrieving the results).

library(ospsuite)
# Load population information from csv
popFilePath <- system.file("extdata", "pop.csv", package = "ospsuite")
myPopulation <- loadPopulation(csvPopulationFile = popFilePath)
print(myPopulation)
#> <Population>
#>   • Number of Individuals: 10

Creating populations

Similar to creating individual parameter sets (see Creating individuals), a population is created from population characteristics created by calling the method createPopulationCharacteristics(). To see the list of available values for the arguments species and population (only for human), use the enums Species and HumanPopulation, respectively. The returned object of type PopulationCharacteristics is then passed to the function createPopulation to generate a set of parameter values. The algorithm behind is the same used in PK-Sim when creating an population. Molecule ontogenies can be added as described in the vignette Creating individuals.

# If no unit is specified, the default units are used. For "height" it is "dm", for "weight" it is "kg", for "age" it is "year(s)".
populationCharacteristics <- createPopulationCharacteristics(
  species = Species$Human,
  population = HumanPopulation$Asian_Tanaka_1996,
  numberOfIndividuals = 50,
  proportionOfFemales = 50,
  weightMin = 30,
  weightMax = 98,
  weightUnit = "kg",
  heightMin = NULL,
  heightMax = NULL,
  ageMin = 0,
  ageMax = 80,
  ageUnit = "year(s)"
)
print(populationCharacteristics)
#> <PopulationCharacteristics>
#>   • Species: Human
#>   • Population: Asian_Tanaka_1996
#>   • Number of individuals: 50
#>   • Proportion of females: 50
#>   • Age: [0.00 year(s)..80.00 year(s)]
#>   • Gestational age: ]-Inf..+Inf[
#>   • Weight: [30.00 kg..98.00 kg]
#>   • Height: ]-Inf..+Inf[
#>   • BMI: ]-Inf..+Inf[

# Create population from population characteristics
result <- createPopulation(
  populationCharacteristics = populationCharacteristics
)
myPopulation <- result$population
print(myPopulation)
#> <Population>
#>   • Number of Individuals: 50

Running population simulation

To run a population simulation, the Population object created by the createPopulation method must be passed to the runSimulation() method:

# Load simulation
simFilePath <- system.file("extdata", "Aciclovir.pkml", package = "ospsuite")
sim <- loadSimulation(simFilePath)

# Run population simulation
simulationResults <- runSimulations(
  simulations = sim,
  population = myPopulation
)[[1]]
print(simulationResults)
#> <SimulationResults>
#>   • Number of individuals: 50
#> For paths:
#>   • Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)
#>   • Organism|VenousBlood|Plasma|Aciclovir|Plasma Unbound

Population simulations are run in parallel on multi-core machines - one core simulates a subset of all individuals defined in the population. By default, the number of cores used equals the maximal number of logical cores available minus one.

The user can change the default behavior by providing custom SimulationRunOptions().

# Load simulation
simFilePath <- system.file("extdata", "Aciclovir.pkml", package = "ospsuite")
sim <- loadSimulation(simFilePath)

# Create a SimulationRunOptions object
simRunOptions <- SimulationRunOptions$new()
print(simRunOptions)
#> <SimulationRunOptions>
#>   • numberOfCores: 3
#>   • showProgress: FALSE

# Change the maximal number of cores to use and show a progress bar during simulation
simRunOptions$numberOfCores <- 3
simRunOptions$showProgress <- TRUE

# Run population simulation with custom options
populationResults <- runSimulations(
  simulations = sim,
  population = myPopulation,
  simulationRunOptions = simRunOptions
)[[1]]
print(populationResults)
#> <SimulationResults>
#>   • Number of individuals: 50
#> For paths:
#>   • Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)
#>   • Organism|VenousBlood|Plasma|Aciclovir|Plasma Unbound

Assigning a population to a simulation

Instead of passing the population to runSimulations() every time, you can assign a population directly to a simulation. Once a population is assigned, the simulation is a population simulation, and running it with runSimulations(simulation) (without the population argument) runs it for the whole population:

sim <- loadSimulation(simFilePath, loadFromCache = FALSE, addToCache = FALSE)

# Assign the population to the simulation
sim$population <- myPopulation

# `isPopulation` tells you how the simulation will be run
sim$isPopulation
#> [1] TRUE

# Run it as a population simulation - no `population` argument needed
simulationResults <- runSimulations(sim)[[1]]
print(simulationResults)
#> <SimulationResults>
#>   • Number of individuals: 50
#> For paths:
#>   • Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)
#>   • Organism|VenousBlood|Plasma|Aciclovir|Plasma Unbound

You can read the assigned population back from the simulation at any time with sim$population.

To turn the simulation back into an ordinary individual simulation, assign NULL to remove the population:

# Switch back to an individual simulation
sim$population <- NULL
sim$isPopulation
#> [1] FALSE

This is simply an alternative to the population argument shown above, with one important difference: aging data can only be supplied via the agingData argument of runSimulations(). There is no aging-data field on the simulation, and assigning a population (or NULL) to sim$population clears any aging data previously set on the simulation, because aging data is only meaningful together with the population it was generated for. For aging population simulations, use the population/agingData arguments of runSimulations(). Assigning the population to the simulation is useful when you want to run several population simulations in a single call: assign a population to each simulation and pass them all to runSimulations() together. (The population argument, by contrast, only works with a single simulation.)

# Load two separate simulations (here from the same file, for illustration)
sim1 <- loadSimulation(simFilePath, loadFromCache = FALSE, addToCache = FALSE)
sim2 <- loadSimulation(simFilePath, loadFromCache = FALSE, addToCache = FALSE)

# Assign a population to each simulation. The population is stored by
# reference: here both simulations share the same `myPopulation` object, so
# modifying it later (e.g. with the `Population` class method
# `myPopulation$setParameterValues()`) would affect both.
# Assign distinct `Population` objects if the simulations should differ.
sim1$population <- myPopulation
sim2$population <- myPopulation

# Both population simulations are run in one call. Naming the input list gives
# the result list readable names; with an unnamed list, results are keyed by
# the simulations' auto-generated ids instead.
results <- runSimulations(list(firstSim = sim1, secondSim = sim2))

# `results` is a named list with one `SimulationResults` object per simulation
results$firstSim
#> <SimulationResults>
#>   • Number of individuals: 50
#> For paths:
#>   • Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)
#>   • Organism|VenousBlood|Plasma|Aciclovir|Plasma Unbound

Note that the population argument of runSimulations() and a population assigned to the simulation can be combined: if you pass a population argument for a simulation that already has a population assigned, the argument takes precedence for that run only, and the simulation’s assigned population is left untouched afterwards.

Simulated time-value pairs for a specific output from the SimulationResults-object returned by the runSimulation method can be accessed with the method getOutputValues. The user can provide either the path(s) of the output (which can be a molecule, a parameter, or an observer), or the object(s) of the type Molecule, Parameter, or Quantity (for observers) with the argument quantitiesOrPaths. If no output is specified, all outputs available in the simulation results are returned.

The paths of all available outputs can be accessed via

populationResults$allQuantityPaths
#> [1] "Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)"
#> [2] "Organism|VenousBlood|Plasma|Aciclovir|Plasma Unbound"

getOutputValues() returns a list with two entries: data and metadata:

data is a dataframe with two predefined columns (IndividualId and Time) as well as one column for each requested output
- IndividualId
- Time a vector with simulated time values (in minutes, equal for all outputs)
- a vector with simulated entries for each output requested.

The values of IndividualId, Time, and the simulated outputs, are appended for each simulated individual. Note that this results in non-monotonously increasing column Time.

# Get simulated results by path
resultsPath <- populationResults$allQuantityPaths[[1]]
print(resultsPath)
#> [1] "Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)"

resultsData <- getOutputValues(
  populationResults,
  quantitiesOrPaths = resultsPath
)

resultsTime <- resultsData$data$Time
resultsValues <- resultsData$data$`Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)`

plot(resultsTime, resultsValues, type = "l")

To get the results for a specific individual or a set of individuals, the argument individualIds of the method getOutputValues() can be specified:

# Get simulated results by path
resultsPath <- populationResults$allQuantityPaths[[1]]
print(resultsPath)
#> [1] "Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)"

# Get only the results for individuals with IDs 1 and 2
resultsData <- getOutputValues(
  populationResults,
  quantitiesOrPaths = resultsPath,
  individualIds = c(1, 2)
)

resultsTime <- resultsData$data$Time
resultsValues <- resultsData$data$`Organism|PeripheralVenousBlood|Aciclovir|Plasma (Peripheral Venous Blood)`

plot(resultsTime, resultsValues, type = "l")

For more information about running simulations, please refer to Running simulations and retrieving the results.