NCA for simulated Cp-Time

library(tidyverse)
library(gastroPlusAPI)

if (rstudioapi::isAvailable()){
  current_working_directory <- dirname(rstudioapi::getSourceEditorContext()$path)
  setwd(current_working_directory)
}

gpx_service <- start_service(verbose=FALSE)

✔ Configured the GastroPlus Service

gpx_service$is_alive()

[1] TRUE

project_path = "../../ProjectFiles/GPX Library.gpproject"

simulation_of_interest <- "Piroxicam 20mg PO capsule"
compound_name <- "Piroxicam"

new_run_name <-       "simulation_run"
pkplus_input_name <-  "NCA_input"
pkplus_run_name <-    "NCA_run"
new_group_name <-     "DataforNCA"
new_series_name <- "CpTime"

adm_type <- PKPlusDoseRoute$Oral

dose_amount_mg <- 20
body_weight_kg <- 78
infusion_time_h <- 0

#Load the project
open_project(project_path)

#Checking user inputs if they exist
if(!(simulation_of_interest %in% get_simulations()$simulation_name)){
  cli::cli_abort("The assigned simulation name is not found in the project. Please check your input!")
  }


#Create run and execute simulation
create_run(run_name = new_run_name, run_type = RunType$BatchSimulation)

add_simulations_to_run(run_name = new_run_name, simulation_names = simulation_of_interest)

execute_run(run_name = new_run_name)

available_conc_series <- get_available_series_output(
  run_name = new_run_name,
  series_type = OutputSeriesType$Concentration,
  simulation_name = simulation_of_interest,
  simulation_iteration_key_name = NULL
)

#Check compartmental Cp-time descriptor first, if its not in the list, then look for PBPK Cp-time descriptor
plasma_concentration_series_descriptor <- paste(
  compound_name,
  CompartmentType$SystemicCirculation,
  StateType$ConcentrationPlasma,
  sep = " - "
)

venous_return_concentration_series_descriptor <- paste(
  compound_name,
  CompartmentType$VenousReturn,
  StateType$ConcentrationPlasma,
  sep = " - "
)

CpTime_descriptor <- ""

if(plasma_concentration_series_descriptor %in% available_conc_series$state_descriptor){
  CpTime_descriptor <- plasma_concentration_series_descriptor
}else if(venous_return_concentration_series_descriptor %in% available_conc_series$state_descriptor) {
   CpTime_descriptor <- venous_return_concentration_series_descriptor
}

series_data <- get_series_data(
  run_name = new_run_name,
  simulation_name = simulation_of_interest,
  series_descriptor = CpTime_descriptor
)

series_df <- series_data %>% select(independent, dependent) %>% filter(independent <= 1e-10 | independent >= 1e-2)  #Filter very small time points due to rounding issue in GPX observed data view, but keep time point of 0.

dependentUnit <- series_data$dependent_unit[1]
independentUnit <- series_data$independent_unit[1]

new_observed_series_data <- ObservedSeriesInformation$new(independent_unit = independentUnit,
                                                          dependent_unit = dependentUnit,
                                                          series = series_df,
                                                          series_name = new_series_name)

new_series_key <- SeriesKey$new(group_type = ObservedDataGroupType$ExposureData,
                                group_name = new_group_name,
                                series_type = ObservedDataSeriesType$UncertainConcentrationSeries,
                                series_name = new_series_name)


set_observed_series(new_series_key, new_observed_series_data)

#PKPlus Section
exposure_data_series_key <- get_experimental_data_inventory(group_type = ObservedDataGroupType$ExposureData) %>% filter(group_name == new_group_name)

dose_amount <- ScalarObservable$new(scalar_type=ScalarType$Mass,  scalar_value=ScalarValue$new(value = dose_amount_mg, unit = "milligram"))
body_weight <- ScalarObservable$new(scalar_type=ScalarType$Mass, scalar_value=ScalarValue$new(value = body_weight_kg, unit = "kilogram"))
infusion_time <- ScalarObservable$new(scalar_type=ScalarType$Time, scalar_value=ScalarValue$new(value = infusion_time_h, unit = "hour"))


series_key <- SeriesKey$new(group_type = exposure_data_series_key$group_type,
                              group_name = exposure_data_series_key$group_name,
                              series_type = exposure_data_series_key$series_type,
                              series_name = exposure_data_series_key$series_name)

#Assigns correct PKPlus function to create input based on administration type assigned by user"
input_config <- PKPlusInputConfiguration$new(dose_route = adm_type,
                                               exposure_series_key = series_key,
                                               dose_amount = dose_amount,
                                               body_weight = body_weight)


add_pkplus_input(pkplus_input_name, input_config)

run_config <- PKPlusRunConfiguration$new(kinetics_type = PKPlusEliminationModel$Linear,
                                         model_types = c(PKPlusAnalysisMethod$NonCompartmentalAnalysis,
                      PKPlusAnalysisMethod$OneCompartment,
                      PKPlusAnalysisMethod$TwoCompartment,
                      PKPlusAnalysisMethod$ThreeCompartment))


add_pkplus_run(pkplus_run_name, pkplus_input_name)
configure_pkplus_run(pkplus_run_name, run_config)

execute_pkplus_run(pkplus_run_name)

#Retrieve PKPlus output
pkplus_output <- get_pkplus_run_output(pkplus_run_name)

PKPlus_checkpoint <- tryCatch(
  {
    pkplus_output[["non_compartmental_analysis_output"]][[1]][[1]]
  },
  error = function(e) {
    warning("PKPlus is not executed properly!", conditionMessage(e))
    FALSE
  }
)

if (PKPlus_checkpoint == FALSE) {
  if(dose_amount_mg == 0){
    cli::cli_abort("Dose amount is 0, please check!")
  } else{}
  if(body_weight_kg == 0){
    cli::cli_abort("Body weight is 0, please check!")
  }else{}
}

#Data Manipulation
NCA_data <- pkplus_output[["non_compartmental_analysis_output"]][[1]][[2]]

NCA_data_longer <- as.data.frame(NCA_data) %>%
  mutate(across(everything(), as.character)) %>%
  pivot_longer(
    cols = everything(),
    names_to  = "column1", 
    values_to = "column2"
  ) %>% filter(column1 != "dose_route") #delete dose_route, unnecessary and violating data manipulation

NCA_df <- NCA_data_longer %>%
  separate(
    column1, 
    into = c("parameter", "attribute"), 
    sep = "\\.(?=[^.]+$)"   # split at the last dot
  ) %>%
  pivot_wider(
    names_from  = attribute,
    values_from = column2
  ) %>% select(parameter, value, unit)

write.table(
  NCA_df, 
  file      = "clipboard", 
  sep       = "\t", 
  row.names = FALSE
)

cli::cli_alert_info("NCA output is generated and copied to the clipboard. It can also be accessed through NCA_df in global environments!")

ℹ NCA output is generated and copied to the clipboard. It can also be accessed through NCA_df in global environments!

glimpse(NCA_df)

Rows: 13
Columns: 3
$ parameter <chr> "auc_first_moment_last", "auc_infinity", "auc_last", "body_m…
$ value     <chr> "12.348", "0.157897", "0.150247", "78", "126.665", "1.623910…
$ unit      <chr> "(mg/mL)*h*h", "(mg/mL)*h", "(mg/mL)*h", "kg", "mL/h", "(mL/…

gpx_service$kill()

[1] TRUE