2.1 Patient Population and Study Design

This monocentric prospective observational study was conducted at the department of nephrology of the Ghent University Hospital. The study was approved by the independent Ethics Committee of the Ghent University Hospital (2016/0973-B670201629362) and was registered on ClinicalTrials.gov (NCT03909698). Written informed consent was obtained from all patients prior to enrolment. Study completion was on 31 March 2022.

This observational study was part of a larger study that evaluated the pharmacokinetics and pharmacodynamics of several antibiotic treatments in patients with end-stage kidney disease treated with trice weekly high-flux intermittent haemodialysis. Up to 30 patients were planned for enrolment in the piperacillin/tazobactam arm. Additionally, patients were only eligible when there was a need for receiving antibiotic treatment for documented or presumed infection. The indication for starting piperacillin/tazobactam and the administered dosing regimen were decided upon by the treating physician. Patients who were below the age of 18 years at the time of enrolment, had known hypersensitivity reactions to β-lactam/β-lactamase inhibitors, or who were pregnant were excluded from participation.

Vascular access for intermittent haemodialysis consisted of a double lumen tunnelled central venous dialysis catheter or an arteriovenous fistula used in either single-needle or double-needle setup. Three different dialysis modes were used: haemodialysis (HD), which is mainly based on diffusion, and pre-dilution hemodiafiltration (HDFpre) and post-dilution hemodiafiltration (HDFpost), in which there is also an important contribution of convection. This information along with other haemodialysis characteristics and patient characteristics were extracted from the patient’s medical records, i.e., age, sex, total body weight (TBW), residual diuresis, dialysis duration, time that dialysis started and ended, dialysis mode, type of vascular access, blood flow rate through the dialyser and dialyser priming volume. The dose, duration of infusion and time of administration were documented. Moreover, the actual dose of piperacillin and tazobactam administered was corrected by determining the amount of piperacillin and tazobactam remaining in the used infusion lines. Piperacillin and tazobactam were administered intravenously in an 8:1 ratio as used in common practice.

2.2 Sampling and Bioanalytical Method

Blood samples were collected over a period of maximum 6 days and at expected time of maximum piperacillin/tazobactam concentration, during the distribution phase and at the estimated trough in between two haemodialysis sessions. These blood samples were taken from a central venous catheter, from an intravenous infusion line, or by direct venipuncture. During the dialysis sessions after study start, blood samples were also taken from the arterial blood line entering the dialyser (dialyser inlet) and venous blood line exiting the dialyser (dialyser outlet) simultaneously at different time points. The exact time of sample collection was recorded throughout the study. All blood samples were collected in lithium–heparin tubes and were centrifuged at 3000 rpm/2095g for 10 min at 4 °C (Beckman Coulter Allegra X-15R) and divided in Sarstedt microtubes 2 mL polypropylene (PP). Analytes were detected using a Q-Exactive Hybrid quadrupole orbitrap mass spectrometer (Thermo Scientific, Waltham, Massachusetts, USA). The lower limit of quantification was 0.5 mg/L for piperacillin and 0.25 mg/L for tazobactam. Within-run precision was < 3.7% for piperacillin and < 6.1% for tazobactam, and between-run precision was < 10.5% for piperacillin and < 15.0% for tazobactam.

2.3 Description of the Previously Developed General-Purpose Model

The development of the general-purpose model, description of the model structure and overall model performance has been reported previously [8]. In short, the piperacillin/tazobactam population pharmacokinetic model consisted of four compartments (a central and peripheral compartment for each piperacillin and tazobactam) with first-order elimination from the central compartments.

For model parameters of the central compartments, bodyweight, scaled to the bodyweight of a reference individual of 70 kg, was implemented following allometric theory and using the conventional fixed West–Brown–Enquist exponents of 0.75 for clearance and 1.00 for volume of distribution. This approach was also used for the peripheral volumes of distribution, but for the intercompartmental clearance (Q2), compartmental allometry was assumed, where Q2 is scaled by the individual estimated size of the peripheral volume of distribution (V2) with a fixed West–Brown–Enquist exponent of 0.75. Age was implemented in the model to account for maturation and decline in clearance pathways during early life and with increasing age using sigmoidal functions. To be aligned with the maturation and decline functions in the general-purpose model, postmenstrual age of the included patient population was calculated by adding an assumed gestational age of 40 weeks to the age of the patients from this study.

Serum creatinine (SCR), as measure of kidney function, was implemented on the piperacillin and tazobactam clearance parameters. The serum creatinine concentrations were standardised to a typical serum creatinine concentration (SCRstd) corresponding to the postmenstrual age of a patient. Between-subject variability was included on all model parameters and could be shared between piperacillin and tazobactam for the central and peripheral volumes of distribution and intercompartmental clearance, without worsening of the model fit to the data, indicating that these parameters are correlated between compounds (e.g. a patient with a high volume of distribution for piperacillin indicates that this patient would also have a high volume of distribution for tazobactam).

2.4 Expanding the General-Purpose Model

Expanding the general-purpose model of piperacillin/tazobactam to include patients with end-stage kidney disease undergoing intermittent haemodialysis using high-flux dialysers required a priori modifications to the model. The modifications were necessary to reflect anticipated differences in endogenous clearance between the general population and the new population and addition of a dialyser clearance to the model to describe the loss of piperacillin/tazobactam from the systemic circulation during intermittent haemodialysis.

Serum creatinine is effectively removed from the systemic circulation during intermittent haemodialysis [9]. Therefore, the influence of kidney function, as implemented in the general-purpose model by an effect of serum creatinine on the endogenous clearance parameters of piperacillin and tazobactam, was modified. For patients with end-stage kidney disease, the serum creatinine-based correction of clearance (FSCR) was fixed to 1, and a new correction term (FESKD) was introduced to the model to describe the residual endogenous clearance of piperacillin and tazobactam in these patients. These modifications are summarised in Equations 1–4.

$$}_}}} = e^}\left( }} \right)/100} \right)} \right)}}}\left( }} \right)/100} }}} \right]}}$$

(1)

$$F_}}} = \left\c} }}} \times \left( }\left( }^ } \right) - }_}}} } \right)}} , } } \\ } } \\ \end } \right.$$

(2)

$$}_}}} = e^}\left( }} \right)/100} \right)} \right)}}}\left( }} \right)/100} }}} \right]}}$$

(3)

$$} = }_}}} \times F_}}} \times F_}}} \times F_}}} \times F_}}} \times F_}}} \times e^}}} }} .$$

(4)

In these equations, serum creatinine (SCR) is standardised (SCRstd) for postmenstrual age (PMA). The influence of serum creatinine is based on the model parameter θSCR fixed to the value estimated in the development of the general-purpose model. The influence of end-stage kidney disease on the endogenous clearance is estimated by θESKD based on data from the end-stage kidney disease patient population. CL is the individual clearance of a patient, and CLnormal is the typical endogenous clearance (CLnormal) of a 35-year-old, 70 kg patient, with a standardised serum creatinine concentration. The influence of maturation (FMAT), decline (FDEC), bodyweight (FSIZE) and between-subject variability (ηCL), which describe a deviation in endogenous clearance from the typical patient, are parameterised as reported previously for the general-purpose model (model code provided in supplementary materials).

The dialyser clearance (CLdia) was implemented by additional separate clearance terms from the central compartment for both piperacillin and tazobactam that was activated by a time-varying indicator variable during the dialysis session. As a result, the total clearance is the sum of the endogenous clearance and the dialyser clearance. The dialyser clearance terms were derived by multiplying the recorded blood flow rate through the blood compartment of the dialyser by an estimated model parameter representing the extraction ratio (ER) as displayed in Equation 5 and 6. Given that the extraction ratio ranges between 0 and 1, these model parameters (θlgtER) were estimated using a logit transform. For dialysis sessions where a single-needle access strategy was used, a time-averaged blood flow rate was used to correct for the differences in blood flow rate as compared with the other types of access strategies.

$$}_}}} = } \times } \, \left( \right)$$

(5)

$$} = \frac}}} }} }}}}} }} }}$$

(6)

As a first step of the development of the modified general-purpose model, θESKD and θlgtER were estimated based on the data of ESKD patients undergoing haemodialysis using high-flux dialysers. In addition, model parameters that describe the residual error were also estimated. Original model parameters of the general-purpose model were kept fixed throughout the analysis to the values reported previously [8].

Subsequently, considering that both concentrations were determined in the same plasma sample, the L2 data item, implemented in NONMEM, was used to evaluate the correlation between model parameters that described the residual error in the piperacillin and tazobactam concentration. Finally, a stepwise covariate search was performed, in which residual diuresis (present/absent) was evaluated as a covariate influencing θESKD, and dialysis characteristics (dialysis mode, type of vascular access, blood flow rate and dialyser priming volume) were evaluated as covariates influencing θlgtER. Continuous covariates were centred on the median and evaluated additively in the logit domain (Equation 7). Categorical covariates were evaluated by estimating a typical parameter for each category for interpretation purposes.

$$\theta_}}} = }\theta_}}} + \theta_}}} \times \left( }}} } \right)$$

(7)

$$\theta_}}} = \left\c} , \, } \, COV = x_ } \\ , \, } \, COV = x_ } \\ \\ , \, } \, COV = x_ } \\ \end } \right.$$

(8)

In Eq. 7, θparameter is the model parameter corrected for covariate effects, TVθparameter is the model parameter for a typical patient, θCOV is the slope of the covariate effect in the logit domain, COV is the observed covariate value and COVmedian is the median of covariate value. In Eq. 8, x1 to xn represent individual categories, and θxi is the model parameter when the covariate value is xi. Additionally, differences in pharmacokinetic parameters between patients with end-stage kidney disease and patients included in the general-purpose model were also evaluated.

Modifications to the model were accepted on the basis of numerical improvement of the overall model fit to the data as derived from the objective function value (OFV) as well as improvement in visual trends observed in the standard goodness-of-fit (GOF) plots and prediction- and variance-corrected visual predictive check (pvcVPC) plots [10]. During structural model development and forward-inclusion step of the stepwise covariate search, modifications to the model were accepted on the basis of a p-value lower than 0.05, which corresponds to an improvement in OFV of at least 3.84 points per additional model parameter. A more conservative backward elimination step was conducted where covariates were removed on the basis of a p-value of 0.001 corresponding to an increase of less than 10.83 points for one additional model parameter.

2.5 Evaluation of Current Dosing Recommendations

Monte Carlo simulations were conducted with the modified general-purpose model to evaluate the influence of impaired kidney function and significant covariates on the pharmacokinetic profile of piperacillin and tazobactam and the probability of target attainment. For this, a total of 1000 template virtual patients were generated. For categorical covariates (e.g. sex), 1000 samples were drawn from the observed distribution in this study (e.g. 35% female to 65% male ratio). For continuous covariates (e.g. bodyweight), 1000 samples were drawn from a normal distribution using the observed mean and standard deviation, but samples were only accepted within the minimum and maximum values of this study. For simulations of patients with a normal kidney function, 1000 samples were drawn from a normal distribution using the observed mean and standard deviation, but samples were only accepted when the estimated glomerular filtration rate (Cockroft–Gault) was higher than 40 mL/min or the maximum observed value in the general-purpose model.

The simulations were conducted using current dosing recommendations (i.e. infusion rate and dose) as approved by the European Medicines Agency (EMA) and the US Food and Drug Administration (FDA). The dosing recommendations were obtained from the Summary of Product Characteristics (SmPC; and the FDA label [1, 2], which are displayed in Tables S1–2 of the Electronic Supplementary Material (ESM), respectively. Typical dialysis characteristics (e.g. dialysis schedule, blood flow rate in the dialyser) were derived from this study. A trice weekly dialysis schedule was assumed, with the duration of each dialysis session being 4 h. Piperacillin and tazobactam concentrations were simulated in approximate steady state (third week after start of treatment).

The bactericidal activity of piperacillin was assumed to be determined by the time fraction for unbound plasma concentration higher than the minimum inhibitory concentration (fT> MIC), and a fT> MIC of minimally 50% was required for clinical efficacy, which was a similar target as used in the simulations of the general-purpose model [8]. Similarly, higher efficacy with respect to preventing bacterial regrowth was assumed to be obtained with a fT> MIC of 100% compared with an fT> MIC of 50% [8]. Therefore, a 50% and 100% fT> MIC were both used as targets of piperacillin to evaluate the probability of target attainment with current dosing recommendations over a 1-week interval with three haemodialysis sessions and interdialytic periods. The piperacillin fraction unbound (fUNB_PIP) was fixed at 64.5% in the probability of target attainment simulations, which was based on the estimate derived in the general-purpose model analysis [8]. For tazobactam, the mean concentration (Cm_TAZ) was calculated according to Equation 9.

$$C_}\_} \left( }^ } \right) = \frac}_}\ }}}}.$$

(9)

In this equation, AUCss,24 h is the area under the total tazobactam concentration versus the time curve for 24 h in steady state.

2.6 Software

Model parameters were estimated using the first-order conditional estimation with interaction algorithm as implemented in NONMEM® (Version 7.5; Icon PLC, Dublin, Ireland), and the uncertainty of the model parameters was based on the covariance step in NONMEM®. Model development was managed with Pirana (Version 3.0.0; Princeton, New Jersey, USA). All pre- and post-processing and model simulations were performed in R® (R Foundation for Statistical Computing, Vienna, Austria), using the R packages “tidyverse” (Version 1.3.2; Wickham H. 2017) and “RxODE” (Version 1.1.5; Fidler M, Hallow M, Wang W. 2015).