Compute the area under the curve (AUC), the maximum and minimum values of a function of time over a given interval, or at steady state.

The AUC over interval \((a,b)\) is computed using the trapezoidal formula for integration \[ I(a,b) = \sum_{j=1}^{n-1} (t_{j+1} - t_j)\left(\frac{f(t_{j+1})+f(t_j)}{2}\right).\] where \(t_1=a < t_2 < \ldots < t_n=b\).

Cmax and Cmin are defined as \[C_{\rm max} = \max_{1\leq j\leq n} f(t_j) \quad ; \quad C_{\rm min} = \min_{1\leq j\leq n} f(t_j)\] Then, Tmax and Tmin are defined by \[C_{\rm max} = f(T_{\rm max}) \quad ; \quad C_{\rm min} = f(T_{\rm min})\]

The times \((t_j , 1 \leq j \leq n)\) are either provided as an input to exposure, or automatically determined if exposure is computed at steady state.


exposure(model, output, group = NULL, treatment = NULL, parameter = NULL, data = NULL, project = NULL, settings = NULL, regressor = NULL, varlevel = NULL)


Input arguments are the input arguments of Simulx.

Specific input arguments can be also used for computing the exposure at steady state, i.e. after the administration of an “infinite” number of doses:

a list with fields
a list with fields


A list with several elements:


Example 1: AUC for a one compartmental model

R script: doc_exposure1.R

Computing the exposure on a given interval

We first compute the exposure over 24 hours for a one-compartment model after a single oral administration.

The concentration is computed between 0 and 24 h every 0.2 h.

myModel <- inlineModel("
input = {ka, V, Cl}
Cc = pkmodel(ka, V, Cl)

p <- c(ka=0.5, V=10, Cl=1)

adm1 <- list(time=0, amount=100)
out1 <- list(name="Cc", time=seq(0, 24, by=0.2))

res1 <- exposure(model=myModel, parameter=p, output=out1, treatment=adm1)

res1 is a list with two elements:

## [1] "Cc"     "output"

res1$output is a list which contains only one data frame Cc with the values of the concentration computed between 0 and 24 h every 0.2 h.

ggplot(data=res1$output$Cc) + geom_line(aes(x=time,y=Cc))

res1$Cc is a data frame which contains several informations, such as the limits of the time interval, the time step, the area under the curve (AUC), the minimum and maximum concentration values over the interval.

##   t1 t2 step      auc tmax    cmax tmin cmin
## 1  0 24  0.2 88.64337    4 6.68731    0    0

We can also compute the exposure after multiple administrations, given every 8 hours in this example.

adm2 <- list(time=seq(0,160,by=8), amount=100)
out2 <- list(name="Cc", time=seq(160, 168, by=0.1))

res2 <- exposure(model=myModel, parameter=p, output=out2, treatment=adm2)
##    t1  t2 step      auc  tmax     cmax tmin    cmin
## 1 160 168  0.1 99.99583 162.6 14.03234  160 9.96636

Computing the exposure at steady state

When the number of doses increases, the system reaches its steady state: the PK profile tends to be the same between successive administrations. We can then compute the exposure between two doses by setting time=‘steady.state’ in the definition of the output.

adm3 <- list(tfd=0, ii=8, amount=100)
out3a <- list(name="Cc", time='steady.state')

res3a <- exposure(model=myModel, parameter=p, output=out3a, treatment=adm3)

ggplot(data=res3a$output$Cc) + geom_line(aes(x=time,y=Cc))

##   t1 t2       step      auc     tmax     cmax tmin     cmin
## 1 48 56 0.08080808 99.53506 50.58586 13.96767   48 9.882421

Steady state is approximated by a multiple dose administration. The number of dose is determined according to the desired precision of the approximation. By default, the tolerance, defined here as the relative error of the AUC, is 1%. Better approximation is obtained with a smaller value of tol.

out3b <- list(name="Cc", time='steady.state', tol=0.001)
res3b <- exposure(model=myModel, parameter=p, output=out3b, treatment=adm3)
##   t1 t2       step      auc     tmax    cmax tmin     cmin
## 1 72 80 0.08080808 99.95537 74.58586 14.0266   72 9.958746

Convergence to steady state is approximately geometric: the relative error decreases as \(e_j \approx \alpha^j\) where \(0 < \alpha < 1\). This relationship is exact in the case of repeated iv bolus administrations, with \(\alpha = 2^{-\delta/t_{1/2}}\), where \(\delta\) is the inter-dose time interval and \(t_{1/2}\) the half-life.

For a given tolerance tol, the number of doses to use is approximately \(\log(\text{tol})/\log(\alpha)\).

ngc is the number of doses used for estimating the convergence rate to steady state, i.e for estimating \(\alpha\). By default, ngc=5.

Using exposure with several individual and groups

We compare in this example the exposure for two different sets of PK parameters. The input argument parameter is a data frame with two rows.

p3 <- inlineDataFrame("
 id   ka   V   Cl
  1  0.5  10    1
  2  0.5   8  0.8
res3c <- exposure(model=myModel, parameter=p3, output=out3a, treatment=adm3)
##   id t1 t2       step       auc     tmax     cmax tmin      cmin
## 1  1 48 56 0.08080808  99.53506 50.58586 13.96767   48  9.882421
## 2  2 48 56 0.08080808 124.41883 50.58586 17.45959   48 12.353026
ggplot(data=res3c$output$Cc) + geom_line(aes(x=time,y=Cc, color=id))

We can also compare the exposure for different dosage regimens. One group receives 100 mg every 8 hours while the other group receives 50 mg every 4 hours.

Function exposure automatically determines a common time interval for both groups (an interval of 8 hours in this example).

g1 <- list(treatment=list(ii=8, amount=100))
g2 <- list(treatment=list(ii=4, amount=50))
out4 <- list(name="Cc", time='steady.state', tol=0.001)

res4 <- exposure(model=myModel, parameter=p, output=out4, group=list(g1,g2))
##   id group t1 t2       step      auc     tmax     cmax tmin      cmin
## 1  1     1 72 80 0.08080808 99.95537 74.58586 14.02660   72  9.958746
## 2  2     2 72 80 0.08080808 99.94708 77.57575 12.90107   72 11.720057
ggplot(data=res4$output$Cc) + geom_line(aes(x=time,y=Cc, color=group))

Example 2: AUC for multiple administrations

R script: doc_exposure2.R

When, for the same individual, several administrations are combined with different times of first dose and different inter-dose time intervals, exposure automatically determines a common time interval. The length of this interval is the least common multiple of the time intervals associated to each administration.

In this example, the inter-dose intervals are, respectively 8h and 12h. The PK profile obtained by combining these two administrations is therefore periodic at steady state and the period is 24h.

myModel <- inlineModel("
input = {ka, V, Cl}

depot(type=1, target=Ad)
depot(type=2, target=Ac)

ddt_Ad = -ka*Ad
ddt_Ac =  ka*Ad - (Cl/V)*Ac
Cc = Ac/V

p <- c(ka=0.5, V=10, Cl=1)

adm1 <- list(tfd=0, ii=8, amount=100, type=1)
adm2 <- list(tfd=3, ii=12, amount=50, type=2, tinf=1)
adm <- list(adm1, adm2)
out <- list(name="Cc", time='steady.state', ntp=200)

res <- exposure(model=myModel, parameter=p, output=out, treatment=adm)

##   t1 t2     step      auc    tmax     cmax tmin     cmin
## 1 48 72 0.120603 398.9503 51.9799 20.15223   48 12.91594
print(ggplot(data=res$output$Cc) + geom_line(aes(x=time,y=Cc)))

Example 3: AUC for multiple outputs and with inter individual variability

R script: doc_exposure3.R

We consider in this example a PKPD model, with inter individual variability on the PKPD parameters.


myModel <- inlineModel("

Cc  = pkmodel(ka, V, k)
Ec = Imax*Cc/(Cc+IC50)
PCA_0 = S0 
ddt_PCA = kout*((1-Ec)*S0- PCA)  


V_pred = V_pop*(w/w_pop)^beta_V
k_pred = k_pop*(w/w_pop)^beta_k

ka  ={distribution=lognormal, prediction=ka_pop,   sd=omega_ka}
V   ={distribution=lognormal, prediction=V_pred,   sd=omega_V}
k   ={distribution=lognormal, prediction=k_pred,   sd=omega_k}
S0  ={distribution=lognormal, prediction=S0_pop,   sd=omega_S0}
IC50={distribution=lognormal, prediction=IC50_pop, sd=omega_IC50}
kout={distribution=lognormal, prediction=kout_pop, sd=omega_kout}
Imax={distribution=logitnormal, prediction=Imax_pop, sd=omega_Imax}

input={w_pop, omega_w}

w={distribution=normal, mean=w_pop, sd=omega_w}

We provide some values for the population PKPD parameters:

p <- c(w_pop=70,     omega_w=10, 
       ka_pop=0.5,   omega_ka=0.2, 
       V_pop=10,     omega_V=0.1,    beta_V= 1,
       k_pop=0.1,    omega_k=0.15,   beta_k=-0.25,
       Imax_pop=0.8, omega_Imax=0.4,
       S0_pop=100,   omega_S0=0, 
       IC50_pop=1,   omega_IC50=0.2,
       kout_pop=0.1, omega_kout=0.1)

We want to compare the exposure for two different dosage regimens, taking into account the inter individual variability of the parameters. In order to investigate the distributions of the AUCs and Cmax, we simulate two arms with a large number of individuals \((N=1000)\).

Because we are interested in the PK and the PD, both outputs are defined.

N <- 1000
adm1 <- list(ii=24, amount=100)
adm2 <- list(ii=12, amount=50)
out <- list(name=c("Cc","PCA"), time='steady.state')

g1 <- list(treatment=adm1, size=N, level='covariate')
g2 <- list(treatment=adm2, size=N, level='covariate')

res <- exposure(model     = myModel, 
                 parameter = p, 
                 output    = out,
                 group     = list(g1,g2))

We can display the empirical distribution of the AUC and the Cmax using boxplots:

b11 <- ggplot(data=res$Cc)+geom_boxplot(aes(x=group,y=auc))
b12 <- ggplot(data=res$Cc)+geom_boxplot(aes(x=group,y=cmax))

We can also display the empirical distribution of the maximum and minimum PCA values over a 24 hours period:

b21 <- ggplot(data=res$PCA)+geom_boxplot(aes(x=group,y=cmin))
b22 <- ggplot(data=res$PCA)+geom_boxplot(aes(x=group,y=cmax))

Here, res$outputis a list with two elements: Cc and PCA computed over a 24 hours period. We can use prctilemlx to display the empirical distribution of these two functions of time (the 10 deciles are displayed for each distribution)

pl1 <- prctilemlx(res$output$Cc, labels=c("group 1", "group 2"))  +  theme(legend.position="none") 
pl2 <- prctilemlx(res$output$PCA, labels=c("group 1", "group 2"))  +  theme(legend.position="none")