--- title: "Getting started" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Getting started} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", eval = FALSE ) ``` The *OSP Global Sensitivity* package facilitates the implementation of one-at-a-time (OAT) and global sensitivity analyses (GSA) of physiologically based pharmacokinetic (PBPK) models built in the [Open Systems Pharmacology (OSP) Suite](https://www.open-systems-pharmacology.org/). The package evaluates the sensitivity of user-selected pharmacokinetic (PK) parameters, such as the maximal concentration (`C_max`) and the area under the curve (`AUC`), for user-selected model output time profiles, with respect to user-selected model input parameters. This article describes how to install the package and walks through the building blocks of a sensitivity analysis workflow. It corresponds to *Supplementary Materials 3* of the accompanying publication: > Najjar A, Hamadeh A, Krause S, Schepky A, Edginton A. Global sensitivity > analysis of Open Systems Pharmacology Suite physiologically based > pharmacokinetic models. *CPT Pharmacometrics Syst Pharmacol.* > 2024;13:2052-2067. doi: > [10.1002/psp4.13256](https://doi.org/10.1002/psp4.13256) ## Installation The *OSP Global Sensitivity* package can be installed and then loaded using the commands: ```{r install} install.packages("path/to/ospsuite.globalsensitivity.zip") library(ospsuite.globalsensitivity) ``` ## Overview of the workflow The figure below provides an overview of the R functions and R6 classes of the *OSP Global Sensitivity* package to which the user has access and which are needed to run a sensitivity analysis workflow. The workflow can be decomposed into three steps: 1. **Workflow setup**: Here the user provides: a. the path to the PKML simulation file and, if performing sensitivity analyses of drug-drug interaction (DDI) ratios, the path to the PKML DDI simulation file, via the `loadSimulation()` function. b. An R list of the model parameters to be analyzed. These parameters are provided in the form of objects of the R6 `SAParameter` class. During the instantiation of a parameter object of this class, a probability distribution is provided in the form of an R6 object of either the `UniformDistribution`, `LogUniformDistribution`, `NormalDistribution`, or `LogNormalDistribution` class. c. An R list of the model outputs to be analyzed as `SAOutput` objects. Each such object includes the path of the model output in the PKML simulations as well as the PK parameters to be evaluated for each of the model outputs. 2. **Sensitivity calculation**: Here, the user may conduct a sensitivity analysis using any of four sensitivity calculation functions: a. `runSU()` conducts one-at-a-time local sensitivity and uncertainty analyses. b. `runMorris()` conducts a Morris screening. c. `runSobol()` and `runEFAST()` conduct variance-based sensitivity analyses using the Sobol and EFAST algorithms respectively. 3. **Result plotting**: Each of the sensitivity calculation methods has its dedicated function for plotting sensitivity analysis results: a. `generateTornadoPlot()` is used to plot the results of a one-at-a-time local sensitivity or uncertainty analysis returned by `runSU()`. b. `generateMorrisPlot()` is used to plot the results of a Morris screening returned by `runMorris()`. c. `generateSobolBarGraph()` and `generateEFASTBarGraph()` plot the results of variance-based sensitivity analyses returned by `runSobol()` and `runEFAST()` respectively. ```{r workflow-overview, echo = FALSE, eval = TRUE, out.width = "100%", fig.cap = "Overview of the R functions and R6 classes of the OSP Global Sensitivity package."} knitr::include_graphics("figures/workflow-overview.png") ``` ## Loading PKML models A model in the required PKML simulation format may be generated from a PK-Sim or MoBi simulation. The user enters the full path to this file to load it using the `loadSimulation()` function: ```{r load} simFilePath <- "path/to/simulation-file.pkml" simulation <- loadSimulation(simFilePath) # A simulation with a DDI can similarly be loaded if sensitivity analysis of # DDI ratios is required: DDIsimFilePath <- "path/to/DDI-simulation-file.pkml" DDIsimulation <- loadSimulation(DDIsimFilePath) ``` ## Specifying model parameters for analysis A convenient way to access the paths of parameters and outputs in the simulation is to use the `getSimulationTree()` function from the *ospsuite* package: ```{r tree} tree <- getSimulationTree(simulation) ``` The user next provides the R list of the model parameters as `SAParameter` objects. Each such object includes the path of every input parameter in the PKML simulations as well as the probability distribution of each input parameter. The inputs for the specification of an input parameter are: - `simulation` (required): the simulation object previously loaded, in which the parameter must be located. - `DDIsimulation` (optional): the previously loaded DDI simulation object, in which the parameter must be located at the exact same path as in the `simulation` object. - `path` (required): the path of the input parameter in the PKML simulation objects `simulation` and `DDIsimulation` (if provided). - `displayName` (optional): a convenient short name for the input parameter that can be used in plots and tables. - `parameterDistribution` (optional): a probability distribution object for the input parameter, selected from among the four R6 classes: `UniformDistribution`, `LogUniformDistribution`, `NormalDistribution`, and `LogNormalDistribution`. - `unit` (optional): for non-dimensionless parameters, this input allows the user to specify the units of the descriptors of the probability distribution of the parameter, such as the mean and variance of a normal distribution. If not provided, the units of the probability distribution are assumed to be the base units of the dimension of the parameter in the PKML simulation. A list of input parameters can be created as follows: ```{r parameters} parametersList <- list( SAParameter$new(simulation = simulation, path = tree$path$to$parameter1$path, displayName = "parameter1", unit = ospUnits$Length$µm, parameterDistribution = LogNormalDistribution$new(mean = 10, CV = 0.5) ), SAParameter$new(simulation = simulation, path = tree$path$to$parameter2$path, displayName = "parameter2", parameterDistribution = UniformDistribution$new(minimum = 0, maximum = 1) ), SAParameter$new(simulation = simulation, path = tree$path$to$parameter3$path, displayName = "parameter3" ), SAParameter$new(simulation = simulation, path = tree$path$to$parameter4$path, displayName = "parameter4" ) ) ``` Here, four input parameters are defined. The first has a path within the simulation encoded by the output of `tree$path$to$parameter1$path`. It has a dimension of `Length` and a `LogNormal` distribution that is best described in units of micrometers (µm) with mean 10 and CV (coefficient of variation) of 0.5. The unit (µm) can easily be specified using the `ospUnits` list provided in *ospsuite*. In contrast, the second parameter has a `Uniform` distribution. This second parameter is dimensionless, and therefore no units for its distribution are specified. Note that the descriptors for the two distributions are different. The `Uniform` and `LogUniform` distributions take parameters `minimum` and `maximum`. The `Normal` distribution takes a `mean` and standard deviation (`stdv`) input. If no distribution is specified for a parameter, as for `parameter3` and `parameter4`, a default `LogUniform` distribution is used that ranges within a factor of 10% of the nominal value of the parameter in the `simulation` object. Note that if the default value of the input parameter is zero, then an upward or downward scaling of the parameter by a factor of 10% will not yield any variation. ## Specifying model outputs for analysis To specify a model quantity for which sensitivity to parameters is to be evaluated, the user inputs a: - `path` (required): the path of the model output in the PKML simulation object, which must also exist in the `DDIsimulation` object if it is provided, - `displayName` (optional): a convenient short name for the model output that can be used in plots and tables, - The PK parameters to be evaluated for each of the model outputs (required). Two model outputs, `output1` and `output2`, and their PK parameters, are added in the following example: ```{r outputs} output1 <- SAOutput$new(simulation = simulation, DDIsimulation = DDIsimulation, path = tree$path$to$output1$path, displayName = "output1") output1$addPKParameter(standardPKParameter = "C_max") output1$addPKParameter(standardPKParameter = "AUC_tEnd") output2 <- SAOutput$new(simulation = simulation, DDIsimulation = DDIsimulation, path = tree$path$to$output2$path, displayName = "output2") output2$addPKParameter(standardPKParameter = "C_max") outputList <- list(output1, output2) ``` Here, the sensitivity of the PK parameters `C_max` and `AUC_tEnd` of the output `output1` will be evaluated with respect to the parameters defined in `parametersList` above. For output `output2`, only the PK parameter `C_max` will be evaluated. Note that `AUC_tEnd` denotes the AUC of the output up to the end of the simulation time, which is set in the simulation in MoBi or PK-Sim. A list of available PK parameters may be printed using the *ospsuite* command `allPKParameterNames()`. It is also possible to evaluate the PK parameters of each output over a specific time period during the simulation, where sensitivity is of particular interest. For example, the user may examine which model input parameters most impact the AUC of a drug's concentration time profile during the early phase after a drug is administered. In such a case, the time period over which the PK parameter is evaluated could be set to be between 0 and 100 minutes as follows: ```{r outputs-time} output1$addPKParameter(standardPKParameter = "AUC_tEnd", pkParameterDisplayName = "AUC_early", startTime = 0, endTime = 100) ``` ## Running one-at-a-time analyses (local and uncertainty analysis) The local and uncertainty analyses may be run together. In the following example, the `runSU()` function is used to compute a local sensitivity and uncertainty analysis for the outputs defined in `outputList` with respect to the parameters defined in `parametersList`. ```{r run-su} su <- runSU(simulation = simulation, DDIsimulation = DDIsimulation, customParameters = parametersList, outputs = outputList, evaluateForAllParameters = FALSE, # Sensitivity analysis parameters: variationRange = 0.2, numberOfSensitivityAnalysisSteps = 2, sensitivityThreshold = 0.1, # Uncertainty analysis parameters: runUncertaintyAnalysis = TRUE, runUncertaintlyOnlyForSensitiveParameters = TRUE, quantiles = c(0.25, 0.75), numberOfUncertaintyAnalysisSamples = 10, saveResults = TRUE, saveFileName = "sensitivityUncertaintyResults.xlsx", saveFolder = "path/to/folder/where/results/are/saved/") ``` Here: - the `simulation` argument takes the simulation object, while the optional `DDIsimulation` argument takes the DDI simulation object, - The `customParameters` argument optionally takes `parametersList` as input. If this argument is provided but `evaluateForAllParameters` is `FALSE`, then the local sensitivity analysis is only evaluated for the parameters in this list, using the probability distributions set for these parameters. If `customParameters` is not provided but `evaluateForAllParameters` is `TRUE`, then the analyses are performed for all constant parameters of the model, assuming a `LogUniform` distribution that scales the parameter upward or downward by up to `variationRange` with respect to its nominal value in the simulation. - The `outputs` argument takes the `outputList` list, which defines the output paths and PK parameters for which local sensitivity and uncertainty is to be analyzed. - When the input `evaluateForAllParameters` is set to `TRUE`, the local sensitivity and uncertainty analyses are evaluated for all constant parameters of the model. - The argument `variationRange` sets the fraction by which the parameters are perturbed in the local sensitivity analysis, and has a default value of 0.1. - The argument `numberOfSensitivityAnalysisSteps` sets the number of steps within the `variationRange` at which to evaluate the local sensitivity analysis. - The argument `sensitivityThreshold`, default value 0.1, can be used to set a local sensitivity value which is the minimum allowable for a parameter to be included among the sensitivity results. This argument can be used to exclude parameters that have little impact on the model outputs when perturbed. - When `runUncertaintyAnalysis` is set to `FALSE`, only the local sensitivity analysis is evaluated. - When `runUncertaintlyOnlyForSensitiveParameters` is set to `TRUE`, the uncertainty analysis is evaluated only for parameters for which the local sensitivity is above `sensitivityThreshold`. - The `quantiles` input is used to set the percentiles of the PK parameters to be evaluated in addition to the 50% and 95% percentiles, which are computed by default. This argument must be a vector of numbers between 0 and 1. - The `numberOfUncertaintyAnalysisSamples` argument is used to set the number of Monte Carlo samples to be drawn from the parameter distributions and at which the PK parameters are to be evaluated for the uncertainty analysis. - The `saveFolder` and `saveFileName` arguments set the folder and name of the Excel `.xlsx` file to which the local sensitivity and uncertainty analysis results are to be saved when `saveResults` is set to `TRUE`. The function `generateTornadoPlot()` may subsequently be used to generate a tornado plot of the local sensitivity and uncertainty analysis results separately, as follows: ```{r su-plots} # Generate a `ggplot` tornado plot for local sensitivity analysis: plt <- generateTornadoPlot(sensitivityDataFrame = su$Results, generateForUncertaintyAnalysis = FALSE) # Generate a `ggplot` tornado plot for uncertainty analysis: plt <- generateTornadoPlot(sensitivityDataFrame = su$Results, generateForUncertaintyAnalysis = TRUE) ``` ## Running the Morris (global) sensitivity method To run the Morris algorithm over 100 trajectories for the previously defined `simulation`, `DDIsimulation`, `parametersList` and `outputList`, the `runMorris()` function is used. The function `generateMorrisPlot()` may subsequently be used to generate a Morris plot: ```{r run-morris} # Run Morris sensitivity analysis: morrisResults <- runMorris(simulation = simulation, DDIsimulation = DDIsimulation, parameters = parametersList, outputs = outputList, numberOfSamples = 100) # Generate a `ggplot` Morris plot for the results: plt <- generateMorrisPlot(morrisResults$Results) ``` ## Running variance-based Global Sensitivity Analyses (Sobol and EFAST) To run the variance-based Sobol or EFAST global sensitivity analysis methods, the functions `runSobol()` and `runEFAST()` are used. The `generateSobolBarGraph()` and `generateEFASTBarGraph()` functions may then be used to generate `ggplot` bar graphs of the results of the two methods: ```{r run-variance} # Run Sobol sensitivity analysis: sobolResults <- runSobol(simulation = simulation, DDIsimulation = DDIsimulation, parameters = parametersList, outputs = outputList, numberOfSamples = 1000) # Run EFAST sensitivity analysis: EFASTresults <- runEFAST(simulation = simulation, DDIsimulation = DDIsimulation, parameters = parametersList, outputs = outputList, numberOfResamples = 1) # Generate a `ggplot` bar graph of Sobol sensitivity analysis results: pltSobol <- generateSobolBarGraph(sobolResults$Results) # Generate a `ggplot` bar graph of EFAST sensitivity analysis results: pltEFAST <- generateEFASTBarGraph(EFASTresults$Results) ``` ## Using the graphical user interface and R script generation An R Shiny graphical user interface (GUI) is provided with the *OSP Global Sensitivity* package to facilitate the setup and execution of the sensitivity algorithms. This app, which consists of four tabs, is launched by first loading the package using the `library(ospsuite.globalsensitivity)` command, followed by the `runGUI()` command: ```{r run-gui} library(ospsuite.globalsensitivity) runGUI() ``` The PKML simulation file, as well as the optional DDI simulation file, may be uploaded via the **Start** tab of the application. In the **Parameters** tab, a tree structure of the simulation paths is generated from the uploaded PKML. This tree enables the user to select the model input parameters to be analyzed. Once the parameters have been selected in the tree, clicking the **Specify distributions** button creates miniature forms through which the user may select the probability distribution of each parameter. In the **Outputs** tab, the user similarly selects the output quantities to be analyzed from a tree structure generated from the PKML simulation. Clicking the **Specify PK parameters** button then generates a miniature form through which the user can select the PK parameters to be analyzed for each model output. In the final **Run sensitivity analyses** tab, the user selects the sensitivity algorithm to run and specifies its run settings. Once run, a progress bar appears in the bottom right of the application. Upon completion, the user may download the sensitivity result in the form of Excel `.xlsx` files, or download visualizations of the sensitivity results. In addition, clicking the **Get code** button corresponding to any of the sensitivity algorithms generates an R script from which the analyses may be run from R. This R script is generated based on the user-specified input parameters and their distributions as well as the user-selected outputs and their PK parameters. ```{r gui-figures, echo = FALSE, eval = TRUE, out.width = "48%", fig.show = "hold", fig.cap = "The four tabs of the accompanying R Shiny app graphical user interface to the OSP Global Sensitivity package."} knitr::include_graphics(c( "figures/gui-start.png", "figures/gui-parameters.png", "figures/gui-outputs.png", "figures/gui-run.png" )) ``` ## Next steps Two fully worked case studies are available as further articles: - [UV filter formulation example](uv-filter-formulation.html) - [Midazolam DDI example](midazolam-ddi.html) A description of the mathematical foundations of the variance-based methods is available in the [Mathematical overview of variance-based methods](variance-based-methods.html) article.