knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) source(system.file("extdata", "vignette-helpers.R", package = "ospsuite.reportingengine"))
require(ospsuite.reportingengine)
This vignette focuses on time profiles and residuals plots (also called goodness of fit plots) generated using either mean model or population workflows.
The time profiles and residuals task first aims at generating time profile plots of specified simulation paths. The task can account for observed data, in which case residuals are calculated and plotted. For population workflows, aggregation functions are used, by default mean, median, 5th and 95th percentiles are computed and plotted in the time profiles.
Results obtained from the simulate
task and stored as csv files within the subdirectory SimulationResults from the workflowFolder
directory are required in order to perform the goodness of fit (plotTimeProfilesAndResiduals
) task.
As a consequence, if the workflow output folder does not include the appropriate simulations, the task simulate
needs to be active.
The objects SimulationSet
(or PopulationSimulationSet
for population workflows) and Output
define which and how the simulated and, if available, observed data will be plotted.
The time profiles and residuals task generates for each simulation set plots of:
For simulations with multiple applications, up to 6 time profile plots can be generated corresponding to simulated time profiles in linear and logarithmic scale for:
The input applicationRanges
from SimulationSet
objects can be used to select which ranges to include in the evaluation.
If observed data are available in the range selected in the output definition and selected in dataSelection
,
they will be added to the time profile plots.
Furthermore, residuals will be calculated according to the residual scale defined by the Output
objects included within the simulation sets
and 6 additional plots will be created figures:
The residual scale corresponds to one element from the enum ResidualScales
:
ResidualScales
The default uses "Logarithmic" and calculate residuals as: $Residuals=log(Observed)-log(Simulated)$
If the option "Linear" is used instead, the residuals are calculated as: $Residuals=Observed-Simulated$
A histogram and a qq-plot of residuals across all the simulation sets will also be added at the end of the analysis.
For population workflows, residuals are calculated from the aggregation of the Simulated output. By default, $median(Simulated)$ is used instead of $Simulated$. It is possible to modify this aggregation function by updating the fields of the object AggregationConfiguration.
The following sections will introduce template scripts including goodness of fit tasks as well as their associated reports. Table 1 shows the features tested in each template script.
examplesTable <- data.frame( Example = seq(1, 7), Workflow = c(rep("Mean model", 5), rep("Population model", 2)), `Simulation Sets` = c(1, 1, 1, 1, 2, 1, 1), Outputs = c(1, 1, 1, 2, 1, 1, 1), `Observed Data` = c("No", "Yes", "Yes", "Yes", "Yes", "No", "Yes"), check.names = FALSE ) knitr::kable(examplesTable, caption = "Table 1: Features tested by each template script")
Example 1 shows one of the simplest workflows including a goodness of fit. In this example, only a simulation time profile is performed since there is no observed data. Furthermore, only one simulation set and one output are requested.
Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Concentration of A", displayUnit = "nmol/ml" ) setA <- SimulationSet$new( simulationSetName = "A", simulationFile = simulationFile, outputs = outputA ) # Create the workflow instance workflow1 <- MeanModelWorkflow$new( simulationSets = setA, workflowFolder = "Example-1" ) # Set the workflow tasks to be run workflow1$activateTasks(c("simulate", "plotTimeProfilesAndResiduals")) # Run the workflow workflow1$runWorkflow()
The output report for Example 1 is shown below. There are only 2 figures in the report (linear and log scales) showing the time profile of the output path 'Organism|A|Concentration in container' displayed as 'Concentration of A'.
cat(includeReportFromWorkflow(workflow1))
Example 2 shows a workflow that used the simulation results performed in Example 1. This example shows 3 new features:
simulate
task Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) # Get the observation file and its dictionary observationFile <- system.file( "extdata", "SimpleData.nmdat", package = "ospsuite.reportingengine" ) dictionaryFile <- system.file( "extdata", "tpDictionary.csv", package = "ospsuite.reportingengine" ) dataSource <- DataSource$new( dataFile = observationFile, metaDataFile = dictionaryFile ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Simulated concentration of A", displayUnit = "µmol/ml", dataSelection = "MDV==0", dataDisplayName = "Observed concentration of A" ) setA <- SimulationSet$new( simulationSetName = "A", simulationFile = simulationFile, outputs = outputA, dataSource = dataSource ) # Create the workflow instance workflow2 <- MeanModelWorkflow$new( simulationSets = setA, workflowFolder = "Example-1" ) # Set the workflow tasks to be run workflow2$activateTasks("plotTimeProfilesAndResiduals") workflow2$inactivateTasks("simulate") # Run the workflow workflow2$runWorkflow()
It can be seen from the run results that the task simulate
was not run.
In this example, the content of the observed data file and its dictionary are respectively described Table 2 and 3.
obsTableExample2 <- readObservedDataFile(observationFile) knitr::kable(obsTableExample2, caption = "Table 2: Observed data content")
dictTableExample2 <- readObservedDataFile(dictionaryFile) knitr::kable(dictTableExample2, caption = "Table 3: Dictionary content")
The output report for Example 2 is shown below.
There are now many more figures in the report showing not only time profiles in µmol/ml
instead of nmol/ml
but also residuals.
cat(includeReportFromWorkflow(workflow2))
In example 3, the dictionary is changed and indicates that the lower limit of quantitation (LLOQ) is included in the observed data.
Like the previous example, there is no need to re-run the simulate
task.
Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) # Get the observation file and its dictionary dataSource <- DataSource$new( dataFile = observationFile, metaDataFile = system.file( "extdata", "tpDictionaryLoq.csv", package = "ospsuite.reportingengine" ) ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Simulated concentration of A", displayUnit = "µmol/ml", dataSelection = "MDV==0", dataDisplayName = "Observed concentration of A" ) setA <- SimulationSet$new( simulationSetName = "A", simulationFile = simulationFile, outputs = outputA, dataSource = dataSource ) # Create the workflow instance workflow3 <- MeanModelWorkflow$new( simulationSets = setA, workflowFolder = "Example-1" ) # Set the workflow tasks to be run workflow3$activateTasks("plotTimeProfilesAndResiduals") workflow3$inactivateTasks("simulate") # Run the workflow workflow3$runWorkflow()
The new dictionary content is described in Table 4.
dictTableExample2 <- readObservedDataFile(dictionaryFile) knitr::kable(dictTableExample2, caption = "Table 4: Dictionary content")
The output report for Example 3 is shown below. The time profiles now indicate the LLOQ.
cat(includeReportFromWorkflow(workflow3))
# Remove the workflow folders unlink(workflow3$workflowFolder, recursive = TRUE)
In example 4, the simulation set now includes 2 output paths.
Output path Organism|A|Concentration in container has observed data to be compared to whereas output path Organism|B|Concentration in container does not.
Since output path Organism|B|Concentration in container was not requested by the examples 1, 2 and 3, it is not included in the simulation results.
As a consequence, the simulate
task needs to be re-run.
Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) # Get the observation file and its dictionary dataSource <- DataSource$new( dataFile = observationFile, metaDataFile = dictionaryFile ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Simulated concentration of A", displayUnit = "µmol/ml", dataSelection = "MDV==0", dataDisplayName = "Observed concentration of A" ) outputB <- Output$new( path = "Organism|B|Concentration in container", displayName = "Simulated concentration of B", displayUnit = "µmol/ml" ) setAB <- SimulationSet$new( simulationSetName = "A and B", simulationFile = simulationFile, outputs = c(outputA, outputB), dataSource = dataSource ) # Create the workflow instance workflow4 <- MeanModelWorkflow$new( simulationSets = setAB, workflowFolder = "Example-4" ) # Set the workflow tasks to be run workflow4$activateTasks(c("simulate", "plotTimeProfilesAndResiduals")) # Run the workflow workflow4$runWorkflow()
The output report for Example 4 is shown below.
cat(includeReportFromWorkflow(workflow4))
# Remove the workflow folders unlink(workflow4$workflowFolder, recursive = TRUE)
In example 5, there are 2 simulation sets, each one includes one output path.
Simulation set A includes output path Organism|A|Concentration in container and has observed data.
Simulation set B includes output path Organism|B|Concentration in container and has no observed data.
The simulate
task needs to be re-run also in this case since the simulation results are saved in files that use the simulation set names.
Code
# Get the pkml simulation file: "MiniModel1.pkml" simulationFileA <- system.file("extdata", "MiniModel1.pkml", package = "ospsuite.reportingengine" ) # Get the pkml simulation file: "MiniModel2.pkml" simulationFileB <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) # Get the observation file and its dictionary dataSource <- DataSource$new( dataFile = observationFile, metaDataFile = dictionaryFile ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Simulated concentration of A", displayUnit = "µmol/ml", dataSelection = "MDV==0", dataDisplayName = "Observed concentration of A" ) outputB <- Output$new( path = "Organism|B|Concentration in container", displayName = "Simulated concentration of B", displayUnit = "µmol/ml" ) setA <- SimulationSet$new( simulationSetName = "A", simulationFile = simulationFileA, outputs = outputA, dataSource = dataSource ) setB <- SimulationSet$new( simulationSetName = "B", simulationFile = simulationFileB, outputs = outputB ) # Create the workflow instance workflow5 <- MeanModelWorkflow$new( simulationSets = c(setA, setB), workflowFolder = "Example-5" ) # Set the workflow tasks to be run workflow5$activateTasks(c("simulate", "plotTimeProfilesAndResiduals")) # Run the workflow workflow5$runWorkflow()
The output report for Example 5 is shown below.
cat(includeReportFromWorkflow(workflow5))
# Remove the workflow folders unlink(workflow5$workflowFolder, recursive = TRUE)
Example 6 shows a goodness of fit for a population workflow. In this example, there is no observed data, only one simulation set and one output are requested.
Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) populationFile <- system.file("extdata", "Pop500_p1p2p3.csv", package = "ospsuite.reportingengine" ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Concentration of A", displayUnit = "µmol/ml" ) setA <- PopulationSimulationSet$new( simulationSetName = "A", simulationFile = simulationFile, populationFile = populationFile, outputs = outputA ) # Create the workflow instance workflow6 <- PopulationWorkflow$new( workflowType = PopulationWorkflowTypes$parallelComparison, simulationSets = setA, workflowFolder = "Example-6" ) # Set the workflow tasks to be run workflow6$activateTasks(c("simulate", "plotTimeProfilesAndResiduals")) # Run the workflow workflow6$runWorkflow()
The output report for Example 6 is shown below. The 2 figures of the report (linear and log scales) include median, mean 5th and 95th percentiles of the population.
cat(includeReportFromWorkflow(workflow6))
Example 7 add observed data to the previous workflow.
Code
# Get the pkml simulation file: "MiniModel2.pkml" simulationFile <- system.file("extdata", "MiniModel2.pkml", package = "ospsuite.reportingengine" ) populationFile <- system.file("extdata", "Pop500_p1p2p3.csv", package = "ospsuite.reportingengine" ) # Get the observation file and its dictionary dataSource <- DataSource$new( dataFile = observationFile, metaDataFile = dictionaryFile ) # Define the input parameters outputA <- Output$new( path = "Organism|A|Concentration in container", displayName = "Simulated concentration of A", displayUnit = "µmol/ml", dataSelection = "MDV==0", dataDisplayName = "Observed concentration of A" ) setA <- PopulationSimulationSet$new( simulationSetName = "A", simulationFile = simulationFile, populationFile = populationFile, outputs = outputA, dataSource = dataSource ) # Create the workflow instance workflow7 <- PopulationWorkflow$new( workflowType = PopulationWorkflowTypes$parallelComparison, simulationSets = setA, workflowFolder = "Example-6" ) # Set the workflow tasks to be run workflow7$activateTasks(StandardPlotTasks$plotTimeProfilesAndResiduals) workflow7$inactivateTasks(StandardSimulationTasks$simulate) # Run the workflow workflow7$runWorkflow()
The output report for Example 7 is shown below. The residuals are calculated on the median simulated values.
cat(includeReportFromWorkflow(workflow7))
# Remove the workflow folders unlink(workflow7$workflowFolder, recursive = TRUE)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.