ABSTRACT
Effective bacterial infection eradication requires not only potent antibacterial agents but also proper dosing strategies. Current practices generally utilize point estimates of the effects of therapeutic agents, even though the actual kinetics of exposure are much more complex and relevant. Here, we use a full time course of the observed in vitro effects to develop a semimechanistic pharmacokinetic-pharmacodynamic model for eravacycline against multiple Gram-negative bacterial pathogens. This model incorporates components such as pharmacokinetics, bacterial life cycle, and drug effects to quantitatively describe the time course of antibacterial killing and the emergence of resistance. Model discrimination was performed by comparing goodness of fit, convergence diagnostics, and objective function values. Models were validated by assessing their abilities to describe bacterial count time courses in visual predictive checks. The final model describes 576 bacterial counts (expressed in log10 CFU per milliliter) from 144 in vitro time-kill experiments with low residual error and high precision. We characterize antibacterial susceptibility as a function of the MIC and adaptive resistance. In doing so, we show that the MIC is proportional to initial susceptibility at 0 h and the development of resistance over the course of 16 h. Altogether, this model may be useful in supporting dose selection, since it incorporates in vitro pharmacodynamics and clinically observed individual drug susceptibilities.
TEXT
Current antibiotic drug development is a global issue as well as a challenge for both scientific and economic reasons (1). Regulatory mechanisms such as the Generating Antibiotic Incentives Now (GAIN) Act and the Limited Population Antibacterial Drug (LPAD) pathway have been put in place to maintain a strong pipeline of effective anti-infectives against current and future pathogens (2, 3). In addition, preclinical models with pharmacokinetic-pharmacodynamic (PK-PD) modeling are increasingly used in an attempt to make drug development more rational and efficient (4). However, PK-PD approaches differ in complexity regarding the method of translating information into dose and schedule selection for clinical study. Several approaches have been developed and are commonly used to inform dosing in antibiotic drug development.
The traditional approach to the selection of doses of antibacterial agents for clinical evaluation in patients relies on the pharmacokinetic-pharmacodynamic indices (PDI) derived from in vitro or in vivo infection systems (5–7). These PDI are composite ratios of a drug’s pharmacokinetic exposure parameter (e.g., the maximum concentration of the free, unbound fraction of the drug in serum [fCmax], the area under the concentration-time curve from 0 to 24 h for the free, unbound fraction of the drug [fAUC0-24], the time during which the free, unbound fraction of the drug exceeds the MIC [fT>MIC]) and the static in vitro microbiological endpoint, mainly the MIC. In general, the relationship between antimicrobial drug exposure and the MIC for an infectious microorganism is associated with microbiological and clinical responses (8). However, it has been well demonstrated that the use of PDI may not be optimal in all situations (9, 10). For instance, Regoes et al. have demonstrated that antimicrobial effectiveness may differ significantly among five different antibiotics against bacteria with the same MIC (10). Additionally, determination of the PDI is significantly influenced by experimental design issues such as the dosing frequency, the MIC assay, and the range of drug exposures studied. This method completely ignores bacterial dynamics, relying exclusively on the kinetic behavior of the drug in the experimental system (e.g., mouse, human) (11).
The use of a semimechanistic PK-PD model presents an alternative option for the translation of nonclinical PK-PD studies. While the PDI is conceptually simple and is commonly employed to determine antibacterial activity based on a 24-h experiment, it does not account for factors that may influence efficacy, such as drug kinetics, changes in bacterial load over time, and the emergence of resistance (8). In contrast, semimechanistic PK-PD models from time-kill experiments integrate these components and may more accurately depict drug-microbe interactions or predict synergistic drug combinations (12–14). Furthermore, these models can be used to forecast probable bacterial drug sensitivity (i.e., killing or growth inhibition) and the emergence of resistance simultaneously. The aim of the current work was to evaluate an in silico semimechanistic PK-PD model based on in vitro time-kill data for a bacteriostatic agent against various bacterial pathogens. This model may be used to evaluate the potential for drug tolerance following a single dose of drug.
RESULTS
Pharmacodynamic model.The final data set included 576 bacterial counts (expressed in log10 CFU per milliliter) from 144 in vitro time-kill experiments. Each in vitro time-kill experiment included 4 observations of bacterial counts over 24 h. There were 4 experiments per strain, each with a different concentration-to-MIC ratio, and 12 strains per pathogen.
A semimechanistic PD model was built to describe the data. The general scheme of the final model is shown in Fig. 1. Parameter estimates for the final model are shown in Table 1. The final structural model was an indirect PD model with inhibition of bacterial growth as the drug effect and with a dynamic susceptibility term (ZEC50). This choice of the structural model was based on the drug’s mechanism of action (15). Bacteriostatic drugs exert their effects by inhibiting bacterial growth instead of exhibiting a direct kill effect, as described in equations 1 and 3 (16). The inhibitory effect (Inh) is described using a standard Hill equation, with SIGM representing the Hill coefficient, C representing the concentration, and Emax representing the maximal effect of the drug. Bacterial growth (equations 2 and 3) is proportional to the product of the bacterial count or density (N) and the growth constant (Kg), with a saturation (Sat) function dependent on N and Nsat (the saturation coefficient). Bacterial death as described in equation 4 is a product of the death rate constant (Kd) and N. The change in bacterial density is then described as the difference between growth and death in equation 5 (where t stands for time). Maximal bacterial density (Nmax) is described in equation 6.
PK-PD model schematic. Squares represent the reservoirs of drug and bacteria. Pointed arrows represent stimulation; flat arrows represent inhibition.
Population parameter estimates for the final model
In the structural model, an adaptive resistance (AR) component was used to describe empirical trends in bacterial regrowth, as has been described previously (17). Instead of using a static susceptibility term, such as the half-maximal effective concentration (EC50), which is traditionally used in the Hill equation to relate concentration to pharmacodynamic effect, this modeling approach used ZEC50, a dynamic EC50 that can change with concentration and time. ZEC50 is defined as the product of initial susceptibility prior to drug exposure (the model-estimated EC50) and the development of resistance over time, as shown in equation 7. Adaptive resistance was defined in this model using equation 8.
During model development, the MIC was empirically found to be proportional to ZEC50 at 24 h, as shown in Fig. 2A (P < 0.05). This relationship was used to redefine EC50 as being proportional to the MIC divided by the development of adaptive resistance at 16 h (AR16), as shown in equations 9 and 10, where q represents a proportionality constant. Equation 8 was then used to redefine ZEC50 as shown in equation 11. The adaptive resistance time point was selected at 16 h postincubation instead of 24 h in order to reflect the timing of the MIC assay per CLSI recommendations (18). The use of 16 h as a time point instead of 24 h decreased the objective function in a non-statistically significant manner but was deemed to be more scientifically appropriate. In this manner, this usage of the adaptive resistance model incorporated the MIC into the dynamic susceptibility term. Normalizing the exponent in the adaptive resistance term by q was also found to decrease the objective function value (P < 0.05) and was included in the final model as shown in equation 12. The relationship between the MIC and ZEC50 in the final model is shown in Fig. 2B.
Relationships between ZEC50 (dynamic EC50) and the MIC in the initial (A) and final (B) model-building stages. These graphs show the relationship between ZEC50 at 24 h and the MIC. Blue lines represent linear regression, with the 95% confidence interval shown as the shaded area.
Furthermore, this model used species-specific covariates for various parameters, as displayed in Table 1. The pharmacodynamic terms represent the typical values for Escherichia coli, while the covariates represent the proportional changes for Acinetobacter baumannii and Klebsiella pneumoniae. No significant differences between species were found for Kg, Nsat, Kd, ARα, or ARβ. However, there were significant differences between species for SIGM, q, and Emax.
Model qualification.Goodness-of-fit (GOF) plots were generated for the model as shown in Fig. 3. These plots generally demonstrate the ability of the model to describe the data as noted by the agreement of the observed and predicted values for bacterial counts.
Goodness-of-fit plots. This figure shows the diagnostic GOF plots of the final model. The boldface dashed lines represent trends in the data. CWRES, conditional weighted residuals.
A prediction- and variability-corrected visual predictive check (VPC) for the full data set is shown in Fig. 4. A prediction-corrected VPC stratified by MICX (the ratio of the concentration of the drug to the MIC) is shown in Fig. 5. This figure demonstrates the ability of the model to describe the data in different scenarios:
• The panel with an MICX value of 0 demonstrates the observed values and model predictions in a situation where the bacteria are incubated without the test drug. In this situation, the bacteria grow to saturation.
• The panels with MICX values of 2 and 4 demonstrate the observed values and model predictions in situations where the bacteria are incubated with the drug at various concentrations. In these situations, there is a decrease in bacterial counts, with some bacterial regrowth for some isolates.
• The panel with an MICX value of 8 demonstrates the observed values and model predictions in a situation where the bacteria are incubated with a relatively high concentration of the drug. In this situation, there is a decrease in bacterial counts without regrowth.
Prediction- and variability-corrected visual predictive check. Circles represent observations of bacterial counts. The solid line represents the median of the data observed, while the dashed lines represent the 2.5th and 97.5th percentiles of the data observed. Shaded areas represent the 95% confidence interval of the median and the 2.5th and 97.5th percentiles of the model-predicted values.
Prediction-corrected visual predictive check stratified by the drug concentration/MIC ratio (MICX). Open circles represent observations of bacterial counts. Solid lines represent the medians of the observations, and dashed lines represent the 2.5th and 97.5th percentiles of the observations. Shaded areas represent the 95% confidence intervals of the median and the 2.5th and 97.5th percentiles of the model-predicted values.
Additionally, Fig. 6 displays a prediction-corrected VPC for the subset of MICX values equal to 0 or 4. These plots further demonstrate the ability of the modeling approach to describe individual observations in multiple scenarios of different MICX values and species. The top panel, with an MICX value of 0, shows each species of bacteria growing to saturation at approximately the same rate. The bottom panel, with an MICX value of 4, demonstrates different likelihoods of developing resistance among the three species of bacteria.
Prediction-corrected visual predictive check stratified by the drug concentration/MIC ratio (MICX) and the species identifier (SPECID), where 1, 2, and 3 refer to A. baumannii, E. coli, and K. pneumoniae, respectively. Open circles represent observations of bacterial counts. The solid line represents the median of the observations; dashed lines represent the 2.5th and 97.5th percentiles of the observations. Shaded areas represent the 95% confidence intervals of the median and the 2.5th and 97.5th percentiles of the model-predicted values.
DISCUSSION
In this analysis, we successfully built a PD model that adequately described bacterial growth and the pharmacodynamic effect of eravacycline on the target bacteria. Most of the model parameters are estimated with high precision (relative standard error, <25%). The GOF plots demonstrate that the model predicts the observed bacterial counts and variability reasonably well. The use of species-specific covariates allowed for better model fitting. The VPCs demonstrate that the model was able to describe bacterial growth, death, and regrowth over a wide range of exposures as signified by MICX and bacterial species.
The semimechanistic approach is an improvement on the traditional PDI approach, given the ability to describe the emergence of resistance in addition to the direct PD effect. A general assumption in antibacterial modeling is that the effect of the drug on the bacteria is the same in the test tube, the animal model, and the human patient. However, bacteria may behave differently in terms of growth and replication in different environments, which may then affect the pharmacodynamic properties. A semimechanistic approach, as outlined in this analysis, allows for the incorporation of multiple physiologically relevant effects beyond the pharmacodynamic effect, such as growth and saturation rates. Ultimately, a better description of the physiological and pharmacodynamic effects in the system leads to better understanding of the drug-bacteria interaction and may allow better bridging to other systems.
Using a simulation-based evaluation of the traditional PDI approach, Kristoffersson et al. observed that the relationships between PK-PD indices and the antibacterial activity (as measured by bacterial reduction, expressed in log CFU) of meropenem against Pseudomonas aeruginosa were sensitive to PK in the host system (e.g., mouse, human) and the exposures studied (dosages [i.e., size and frequency]) (19). The interplay between these factors influenced the identification of the PDI best correlated with antibacterial activity as well as the magnitude required for an explicit effect (e.g., a 2-log10 CFU reduction), impacting dosing choices. This suggests to us that the PDI concept may be applicable only in certain scenarios, such as when the PK in the experiment where PDI is quantified is similar to the PK in the target patient group. For example, critically ill patients may experience hyperdynamic states that will influence a drug’s PK significantly, and this PK is likely to be dissimilar to PK under the experimental conditions from which the applicable PDI was evaluated. These patients, then, are at risk for inadequate drug exposure, which may, in turn, decrease the effectiveness of bacterial killing and antibiotic resistance minimization.
The semimechanistic approach appeared to be sufficient in describing bacterial population dynamics. A more fully mechanistic approach may have included the proportions of resting and dividing subpopulations of bacteria so as to better reflect the bacterial growth cycle. While this approach may describe the underlying biology more extensively, the data requirement for such a model is steeper. For the data presented here, the adaptive resistance model described the trends in initial population decrease and regrowth observed in the experiments. Other reports have used the adaptive resistance approach to describe their data, although to our knowledge, this is the first time that this approach has been applied to a tetracycline-class antibacterial (17, 20). Studies with longer durations and incubation periods will be needed to further understand the exact mechanism for the reduced killing observed after a period of time following drug exposure.
During modeling, we empirically developed a novel relationship between susceptibility and the MIC. Our analysis showed that the MIC is proportional to the product of initial susceptibility (i.e., model-estimated EC50) and the adaptive resistance at 16 h (AR16). This relationship suggests that the MIC is not an instantaneous measure of susceptibility prior to treatment. Because bacteria are incubated with the drug for 16 to 20 h during the MIC assay, resistance likely developed and was propagated, so that the MIC value reflects a combination of susceptibility prior to treatment and adaptive drug resistance.
Additionally, we used the relationship between dynamic antibacterial sensitivity (i.e., ZEC50) and the MIC (equation 11) to update the calculation of adaptive resistance (equation 12). The term q reflects the susceptibility of the organism beyond the susceptibility measured by the MIC. Given that the MIC is measured only in exponents of 2 (e.g., 0.25, 0.5, 1, 2), q may account for the difference between the measured value of the MIC and the precise value of the MIC given an MIC assay with theoretically perfect precision. This assumption also agrees with the variability term associated with q, i.e., 0.336, which appears sufficient to describe a distribution between successive exponents of 2. Furthermore, q may also account for differences between conditions in the MIC assay and the in vitro time-kill experiment.
The inclusion of the MIC in the model allowed us to pool data across three different pathogens (12 isolates of each) while still capturing the pharmacodynamics of each isolate by its MIC level and by adding relevant covariates. This approach minimized overparameterization of individual fits and improved global parameter confidence. Additionally, this approach may allow the simulation of alternative dosing strategies and the resulting bacterial densities, which may provide an additional way to identify dose regimens that better prevent bacterial drug adaptation.
One potential limitation of our analysis is that our model assumes that regrowth of the bacteria is due to an adaptive mechanism. At the same time, we cannot exclude technical factors, such as drug degradation. However, negligible loss of the drug is assumed based on the in vitro static time-kill experimental design. The stability data for the test drug preclude significant chemical breakdown of the drug over 24 h. Additionally, unlike certain classes of antibacterials, such as beta-lactams, bacteriostatic agents are more likely to develop resistance by efflux pumps or binding site mutations, neither of which destroy the drug (21).
Moreover, while this model was able to describe the pharmacodynamics of eravacycline in the current experiment on the basis of in vitro time-kill data, further evaluation will be needed. This modeling strategy will need to be evaluated in the setting of dynamic concentration experiments, which will be more relevant for the expected clinical use of antibacterials. In a dynamic setting, the descriptive potential of the semimechanistic modeling approach could be compared to that of the traditional PDI approach. This modeling strategy could also be extended to other bacteriostatic antibacterials and to antibacterials from other classes.
Overall, our model allows us to capitalize on the strengths of the MIC assay while minimizing its weaknesses. Since the MIC assay is widely used in patient treatment to identify suitable antibacterial agents, a modeling strategy using the MIC assay has the potential to be more relevant to health care practitioners. By separating the effects of initial susceptibility and adaptive resistance as outlined in equation 9, we are better able to use the information from the MIC to inform the model. Using this model as a platform with necessary adjustments for the clinical environment, drug developers may be better able to select doses for the desired antibacterial activity and resistance suppression by simulating and evaluating changes in bacterial profiles under different dose regimens of the drug in silico. At the same time, it should be noted that this particular modeling strategy may not be optimal in all situations, and further investigation will be needed to identify how this model may be used most effectively. Altogether, this semimechanistic modeling approach to describing antibacterial pharmacodynamics presents an attractive alternative to the traditional PDI method.
MATERIALS AND METHODS
Data source.Eravacycline (Xerava) was approved for the treatment of complicated intra-abdominal infections in August 2016. Data sets were acquired from those submitted to the FDA by the eravacycline drug development program. One study consisted of in vitro time-kill data from E. coli, K. pneumoniae, and A. baumannii (12 isolates per species), three bacterial species that are common in complicated intra-abdominal infections. A mixture of clinical and laboratory isolates with MICs ranging from 0.016 to 2 μg/ml were included in the analysis. The bacteria exhibited various resistance mechanisms, which were confirmed by PCR.
Bacterial culturing, antibiotics and MIC determination, and time-kill curve experiments.Bacterial cultures and antibiotic MIC determination were performed in accordance with the guidelines of the Clinical and Laboratory Standards Institute (CLSI) (22). Briefly, bacteria from a 6-h logarithmic-growth-phase culture were diluted to obtain a start inoculum of 1.0 × 106 CFU/ml. Eravacycline was added to each inoculum to obtain a series of working concentrations corresponding to 0, 2×, 4×, or 8× MIC. Cultures were incubated at 37°C, and bacterial growth was determined at 0.25, 2, 7, and 24 h.
Pharmacodynamic analysis.Data were prepared using R, v4.3.3. Model fitting was performed using NONMEM, v7.3. The pharmacokinetic compartment was described by a static concentration of the drug, assuming negligible drug loss over the course of the experiment. Model development proceeded by assessing various model structures to link the concentration of eravacycline to bacterial counts, measured in log10 CFU, using initial structural models as reported in the literature (15). From there, alternate bacterial growth saturation and resistance pathways were considered for incorporation into the model (17, 23). Trends in the data were used to suggest other model structures. Models were selected by comparing goodness-of-fit (GOF) plots, convergence diagnostics, and the likelihood ratio test as described previously (24).
Additive, proportional, and combined error models were assessed to describe residual variability. The first-order conditional estimation (FOCE) method with or without interaction was used to estimate model parameters. Interaction was used with the FOCE method for estimating with either a proportional or a combined error model. An exponential error model was used to describe interexperiment variability.
Species, strain, MIC, and drug concentration-to-MIC ratio were considered as potential covariates on each pharmacodynamic parameter. The relationships between the covariates and pharmacodynamic parameters were described using linear equations. These covariates were then selected based on decreases in the objective function, relative standard error estimates, changes in GOF plots, and physiological plausibility.
Internal validation techniques were used to qualify the final model. Visual predictive checks (VPCs) assessed the ability of the model to describe changes in bacterial counts and the associated variability (25, 26). The visual predictive checks were performed by using the model to create 1,000 simulated data sets, which were then compared to the original data.
ACKNOWLEDGMENTS
We gratefully acknowledge the contributions of John Lazor and Kellie Reynolds in offering valuable comments on the manuscript.
No funding was received for this work.
The opinions or assertions contained herein are the private views of the authors and are not to be construed as official, or as reflecting true views of the Department of the Army, the Department of Defense, or the Department of Health and Human Services.
FOOTNOTES
- Received 23 June 2020.
- Accepted 24 June 2020.
- Accepted manuscript posted online 29 June 2020.
- This is a work of the U.S. Government and is not subject to copyright protection in the United States. Foreign copyrights may apply.
