🧪 Diving into PK/PD for the first time — simulating ceftriaxone with mrgsolve in R. Free drug levels were… surprisingly high? Even pushed it to q48h dosing out of curiosity and the results left me with more questions than answers 🤔📈 Motivations Learning pharmacokinetics (PK) and pharmocodynamics (PD) have always been an interest of mine. It’s always challenging to read through these population PK papers with all the numbers etc. What’s a better way of diving into the surface of these and see if we can at least know how to code a simulation to obtain the probability target attainment (PTA) of different minimal inhibitory concentration (mic) and learn the basics via code! Let’s dive on! Disclaimer: I am not a pharmacist and not an expert in PK/PD. This is a documentation for my own learning and for educational purposes only. Not a medical advice. If you noticed anything wrong here, please let me know! Objectives: What Is Population PK What Are The Parameters of Interest On a Paper? Let’s Code Different CrCl Low Albumin ?q48 Dosing Oppotunities For Improvement Lessons Learnt What Is Population PK Population pharmacokinetics (popPK) is a statistical approach that describes how medications behave in the body across groups of people, accounting for variability between individuals. Instead of studying one person intensively, popPK analyzes sparse data from many patients to understand typical medication behavior and why people differ in their medication exposure. We’ll use Garot D et al Population pharmacokinetics of ceftriaxone in critically ill septic patients: a reappraisal as an example for learning. What Are The Parameters of Interest On a Paper? From the paper, we can see that there are a lot of parameters and numbers. But what are the parameters of interest? We will focus on the following parameters: Looking at their table 3, we see these values: CL = $\theta_1 + \theta_2 . (CL_{cr}/4.26)$ $\theta_1$ : Non-renal (baseline) clearance component $\theta_2$ : Renal clearance scaling coefficient V1 : Volume of distribution of the central compartment. V2 : Volume of distribution of the peripheral compartment. Q : Inter-compartmental clearance. $\omega^2 (CL)$ : Between-subject variability of clearance. $\omega^2 (V1)$ : Between-subject variability of volume of distribution of the central compartment. $\omega^2 (V2)$ : Between-subject variability of volume of distribution of the peripheral compartment. These are the parameters we’ll use in our mrgsolve model. I’ve always wondered what these parameters represent and it was a bit difficult to conceptualize until we dove into the code and finally understood the rationale! It’s a mixed effect model where the estimates were modeled as a function of the fixed effect (theta) and the random effect (eta), you will see this in the code later. The fixed effect represents the typical value of the parameter in the population, while the random effect represents the variability between individuals. The random effect is assumed to be normally distributed with a mean of zero and a variance of omega squared. If we were to draw a flow chart of the above, it will look something like this: One starts with medication being administered into the central compartment, and then it goes into either peripheral compartment (tissue etc) and medication clearance. Notice that the Q is a bidirectional flow between central and peripheral, whereas all other directions are either into central or out from central to clearance. This is very helpful for me to get a surface understanding of the distribution. Let’s get on with the code! Let’s Code library(mrgsolve) library(tidyverse) mod as_tibble() ETA in the above means the random effect, which is assumed to be normally distributed with a mean of zero and a variance of omega squared. ETA is a greek letter (eh-ta). EPS here is Epsilon. For ev, amount is dosing in mg; rate is amount given per hour; ii is frequency; addl is number of additional doses; cmt is the compartment where the dose is given. In this case, we are giving 2000 mg of ceftriaxone as a 30 minute infusion every 24 hours for 3 doses (1 initial dose + 2 additional doses) into the central compartment. We will then need to set seed for reproducibility, pipe in your initial model with dosing, then nid is number of individuals you want to simulate, end is the end time of the simulation in hours, and delta is the time interval for the simulation output in hours. In this case, we are simulating 1000 individuals for 72 hours with a time interval of 0.25 hours (15 minutes). Next, we’ll calculate the probability of target attainment (PTA) for different minimal inhibitory concentration (MIC) values. The PTA is the probability that the free drug concentration exceeds the MIC for a certain percentage of the dosing interval. MIC filter(time >= 48) |> group_by(ID) |> summarise(fT = mean(Cp_free > MIC)) |> summarise(PTA = mean(fT >= 0.50)) |> pull())) ## [1] "Probability of Target Attainment: 0.996" sims |> ggplot(aes(x=time,y=Cp_free,group=ID)) + geom_line(alpha=0.01) + geom_hline(yintercept = MIC, color = "red") + theme_bw() In the above we choose mic of 1, filtered off time after 48 hours for steady state, then calculate the average free ceftriaxone that is above the mic, then assess the mean of times where free ceftriaxone is above 50% per simulated subject. We can see that the probability of target attainment is around 99.6%. We can also visualize the free ceftriaxone concentration over time with a red line indicating the mic of 1. Not too shabby! Now let’s assess when there is a difference in CrCl and albumin. I’ll spare you the code. Changes in CrCl code library(glue) crcl_vec as_tibble() MIC = 48) |> group_by(ID) |> summarise(fT = mean(Cp_free > MIC)) |> summarise(PTA = mean(fT >= 0.50)) |> pull(), " ,CrCl: ", crcl_vec_i[crcl_vec==crcl], "ml/min") plot ggplot(aes(x=time,y=Cp_free,group=ID)) + geom_line(alpha=0.01) + geom_hline(yintercept = MIC, color = "red") + theme_bw() + ggtitle(pta) plot(plot) } That’s interesting! That makes sense, increased CrCl will increase clearance of ceftriaxone, hence decrease in PTA. It’s still pretty good though! Though, what is considered acceptable? 90%? 70%? 50%? Also, the above PTA is based on 50% of a time free ceftriaxone is above MIC. What is the acceptable number for that then? 🤷‍♂️ What if, since ceftriaxone is albumin bound, if we model albumin into the model as well? Changes With Hypoalbuminemia Notice that our initial model had a fixed fraction of unbound (fu) of 0.1, which is the middle of the range reported in the package insert. However, in critically ill patients, hypoalbuminemia is common and can lead to an increase in the fraction of unbound drug, which can affect the pharmacokinetics and pharmacodynamics of ceftriaxone. Let’s see how we can model this in our mrgsolve code. From the paper in method section, they used the formula below to estimate free ceftriaxone from total ceftriaxone: We’ll add that to our model and adjust the np (total concentration of protein binding sites) according to estimate with lower albumin (np = 295), this number again was from the paper in the discussion portion where their median albumin were ~25g/L. code mod mrgsim(nid = 1000, end = 72, delta = 0.25) |> as_tibble() MIC = 48) |> group_by(ID) |> summarise(fT = mean(Cp_free > MIC)) |> summarise(PTA = mean(fT >= 0.50)) |> pull(),", Albumin: ~25g/L, CrCl: ~63 ml/min") plot ggplot(aes(x=time,y=Cp_free,group=ID)) + geom_line(alpha=0.01) + geom_hline(yintercept = MIC, color = "red") + # geom_text(aes(x=20,y=150,label=pta)) + theme_bw() + ggtitle(pta) plot(plot) Wow, that’s interesting! After we correctly fit in the free ceftriaxone estimation, it actually improved the PTA even when albumin is lower. What if we make albumin even lower to ~15g/L (np=~172), and increase our CrCl to 180 ml/min, and increase our fT >= 0.7 (more than 70% of the time free ceftriaxone is above mic), and see if we’ll be able to clear the medication faster? PTA is still 100% !?!?! wow, ceftriaxone 2g really is a beast! Hmmm.. The free ceftriaxone is REALLY high, around ~200-300, can we simulate a q48h dosing and see what the PTA is like, even for our worse case scenario, low albumin, high CrCl, and stil cover ft>mic >= 70%? ?q48 Dosing code dosing ev(dosing) |> mrgsim(nid = 1000, end = 144, delta = 0.25) |> as_tibble() MIC = 48) |> group_by(ID) |> summarise(fT = mean(Cp_free > MIC)) |> summarise(PTA = mean(fT >= 0.70)) |> pull(),", Albumin: ~15g/L, CrCl: ~180 ml/min, fT > mic >= 70%, q48h dosing") plot ggplot(aes(x=time,y=Cp_free,group=ID)) + geom_line(alpha=0.01) + geom_hline(yintercept = MIC, color = "red") + theme_bw() + ggtitle(pta) plot(plot) Seriously!? PTA is still so high !? What does this actually mean? Is there literature on this? Maybe my code is not right… 🤔🤷‍♂️ Let’s look at fT > mic >= 99%. ## [1] "PTA: 0.979, Albumin: ~15g/L, CrCl: ~180 ml/min, fT > mic >= 99%, q48h dosing" If you know anything about this, please let me know! this is for organism with mic

R-bloggers

A beginner's walkthrough of pharmacokinetic/pharmacodynamic (PK/PD) simulation using the mrgsolve package in R. The author implements a two-compartment population PK model for ceftriaxone based on a published paper, simulating 1000 virtual patients via Monte Carlo methods to calculate Probability of Target Attainment (PTA) across different MIC values. The tutorial covers model setup with fixed and random effects (theta/eta), dosing event specification, and explores scenarios including varying creatinine clearance, hypoalbuminemia with a nonlinear protein binding equation, and experimental q48h dosing. Results showed surprisingly high PTA even under worst-case conditions, raising questions the author couldn't fully explain.

Learning PK/PD Simulation: A Beginner’s Monte Carlo Analysis With mrgsolve in R

What Are The Parameters of Interest On a Paper?