| clc {tci} | R Documentation |
Simulate closed-loop control
Description
Simulate closed-loop control using Bayesian updates. Infusion rates are calculated using 'pkmod_prior' to reach 'target_vals' at 'target_tms'. Data values are simulated using 'pkmod_true' at 'obs_tms'. 'pkmod_prior' and 'pkmod_true' do not need to have the same structure. Model parameters are updated at each update time using all available simulated observations. Processing delays can be added through the 'delay' argument, such that observations aren't made available to the update mechanism until 'update_tms >= obs_tms + delay'.
Usage
clc(
pkmod_prior,
pkmod_true,
target_vals,
target_tms,
obs_tms,
update_tms,
type = c("effect", "plasma"),
custom_alg = NULL,
resp_bounds = NULL,
delay = 0,
seed = NULL
)
Arguments
pkmod_prior |
'pkmod' object describing a PK/PK-PD model that is used to calculate TCI infusion rates and is updated as data are simulated and incorporated. Must have an associated Omega matrix. |
pkmod_true |
‘pkmod' object describing the patient’s "true" response. This model will be used to simulate observations. |
target_vals |
A vector of numeric values indicating PK or PD targets for TCI algorithm. |
target_tms |
A vector of numeric values indicating times at which the TCI algorithm should begin targeting each value. |
obs_tms |
Times at which data values should be simulated from 'pkmod_true'. |
update_tms |
Times at which 'pkmod_prior' should be updated using all available simulated observations. |
type |
Type of TCI algorithm to be used. Options are "plasma" and "effect". Defaults to "effect". Will be overwritten if 'custom_alg' is non-null. |
custom_alg |
Custom TCI algorithm to overwrite default plasma- or effect-site targeting. |
resp_bounds |
Optional vector of two values indicating minimum and maximum values possible for the response. |
delay |
Optional numeric value indicating a temporal delay between when observations are simulated and when they should be made available for updating 'pkmod_prior'. For example, a delay should be set to account for a processing time delay in Bispectral Index measurements or the time required to measure drug concentrations from collected samples. |
seed |
An integer used to initialize the random number generator. |
Examples
prior_vcov <- matrix(diag(c(0.265,0.346,0.209,0.610,0.565,0.597,0.702,0.463)),
8,8, dimnames = list(NULL,c('cl','q2','q3','v','v2','v3','ke0','sigma_add')))
pkmod_prior <- pkmod(pars_pk = c(cl = 10, q2 = 2, q3 =20, v = 15, v2 = 30, v3 = 50, ke0 = 1.2),
sigma_add = 0.2, log_response = TRUE, Omega = prior_vcov)
pkmod_true <- pkmod(pars_pk = c(cl = 16, q2 = 4, q3 =10, v = 20, v2 = 20, v3 = 80, ke0 = 0.8),
sigma_add = 0.1, log_response = TRUE)
target_vals <- c(2,3,4,3,3)
target_tms <- c(0,5,10,36,60)
obs_tms <- c(1,2,4,8,12,16,24,36,48)
update_tms <- c(5,15,25,40,50)
sim <- clc(pkmod_prior, pkmod_true, target_vals, target_tms, obs_tms, update_tms, seed = 200)
len <- 500
tms <- seq(0,60,length.out = len)
# true, prior, posterior plasma concentrations
df <- data.frame(time = rep(tms,3),
value = c(predict(pkmod_true, sim$inf,tms)[,1],
predict(pkmod_prior, sim$inf,tms)[,1],
predict(sim$pkmod_post, sim$inf, tms)[,1]),
type = c(rep("true",len),rep("prior",len),rep("posterior",len)))
library(ggplot2)
ggplot(df, aes(x = time, y = value, color = type)) +
geom_line() +
geom_point(data = sim$obs, aes(x = time, y = obs), inherit.aes = FALSE) +
geom_step(data = data.frame(time = target_tms, value = target_vals),
aes(x = time, y = value), inherit.aes = FALSE)
# PK-PD example with observation delay (30 sec)
prior_vcov <- matrix(diag(c(0.265,0.346,0.209,0.702,0.242,0.230)),6,6,
dimnames = list(NULL,c('cl','q2','q3','ke0','c50','sigma_add')))
pkpdmod_prior <- update(pkmod_prior, pars_pd = c(c50 = 2.8, gamma = 1.47, e0 = 93, emx = 93),
pdfn = emax, pdinv = emax_inv, sigma_add = 4, log_response = FALSE, Omega = prior_vcov)
pkpdmod_true <- update(pkmod_true, pars_pd = c(c50 = 3.4, gamma = 1.47, e0 = 93, emx = 93),
pdfn = emax, pdinv = emax_inv, sigma_add = 3, log_response = FALSE)
target_vals <- c(75,60,50,50)
target_tms <- c(0,3,6,10)
obs_tms <- seq(1/6,10,1/6)
update_tms <- seq(1,10,0.5)
## Not run:
sim_pkpd <- clc(pkpdmod_prior, pkpdmod_true, target_vals, target_tms, obs_tms,
update_tms, seed = 201, delay = 0.5)
# plot results
tms <- seq(0,10,length.out = len)
df <- data.frame(time = rep(tms,3),
value = c(predict(pkpdmod_true, sim_pkpd$inf,tms)[,"pdresp"],
predict(pkpdmod_prior, sim_pkpd$inf,tms)[,"pdresp"],
predict(sim_pkpd$pkmod_post, sim_pkpd$inf, tms)[,"pdresp"]),
type = c(rep("true",len),rep("prior",len),rep("posterior",len)))
library(ggplot2)
ggplot(df, aes(x = time, y = value, color = type)) +
geom_line() +
geom_point(data = sim_pkpd$obs, aes(x = time, y = obs), inherit.aes = FALSE) +
labs(x = "Hours", y = "Bispectral Index") + theme_bw() +
geom_vline(xintercept = update_tms, linetype = "dotted", alpha = 0.6) +
geom_step(data = data.frame(time = target_tms, value = target_vals),
aes(x = time, y = value), inherit.aes = FALSE)
## End(Not run)