Extension of the SAEM algorithm to left-censored data in nonlinear mixed-effects model: Application to HIV dynamics model

A modelling framework to characterize the impact of antibiotics on the gut microbiota diversity

Metagenomic sequencing deepened our knowledge about the role of the intestinal microbiota in human health, and several studies with various methodologies explored its dynamics during antibiotic treatments. We compared the impact of four widely used antibiotics on the gut bacterial diversity. We used plasma and fecal samples collected during and after treatment from healthy volunteers assigned to a 5-day treatment either by ceftriaxone (1 g every 24 h through IV route), ceftazidime/avibactam (2 g/500 mg every 8 h through IV route), piperacillin/tazobactam (1 g/500 mg every 8 h through IV route) or moxifloxacin (400 mg every 24 h through oral route). Antibiotic concentrations were measured in plasma and feces, and bacterial diversity was assessed by the Shannon index from 16S rRNA gene profiling. The relationship between the evolutions of antibiotic fecal exposure and bacterial diversity was modeled using non-linear mixed effects models. We compared the impact of antibiotics on gut microbiota diversity by simulation, using various reconstructed pharmacodynamic indices. Piperacillin/tazobactam was characterized by the highest impact in terms of intensity of perturbation (maximal [IQR] loss of diversity of 27.3% [1.9; 40.0]), while moxifloxacin had the longest duration of perturbation, with a time to return to 95% of baseline value after the last administration of 13.2 d [8.3; 19.1]. Overall, moxifloxacin exhibited the highest global impact, followed by piperacillin/tazobactam, ceftazidime/avibactam and ceftriaxone. Their AUC between day 0 and day 42 of the change of diversity indices from day 0 were, respectively, -13.2 Shannon unit.day [-20.4; -7.9], -10.9 Shannon unit.day [-20.4; -0.6] and -10.1 Shannon unit.day [-18.3; -4.6]. We conclude that antibiotics alter the intestinal diversity to varying degrees, both within and between antibiotics families. Such studies are needed to help antibiotic stewardship in using the antibiotics with the lowest impact on the intestinal microbiota.

Mechanistic Learning for Predicting Survival Outcomes in Head and Neck Squamous Cell Carcinoma

We employed a mechanistic learning approach, integrating on-treatment tumor kinetics (TK) modeling with various machine learning (ML) models to address the challenge of predicting post-progression survival (PPS)-the duration from the time of documented disease progression to death-and overall survival (OS) in Head and Neck Squamous Cell Carcinoma (HNSCC). We compared the predictive power of model-derived TK parameters versus RECIST and assessed the efficacy of nine TK-OS ML models against conventional survival models. Data from 526 advanced HNSCC patients treated with chemotherapy and cetuximab in the TPExtreme trial were analyzed using a double-exponential model. TK parameters from the first line and maintenance (TKL1) or after four cycles (TK4) were used to predict PPS and post-cycle 4 OS (OS4), combined with 12 baseline parameters. While ML algorithms underperformed compared to the Cox model for PPS, a random survival forest was superior for OS prediction using TK4 and surpassed RECIST-based metrics. This model demonstrated unbiased OS4 prediction, suggesting its potential for improving HNSCC treatment evaluation. Trial Registration: ClinicalTrials.gov identifier: NCT02268695.

Fast standard error estimation for joint models of longitudinal and time-to-event data based on stochastic EM algorithms

Maximum likelihood inference can often become computationally intensive when performing joint modeling of longitudinal and time-to-event data, due to the intractable integrals in the joint likelihood function. The computational challenges escalate further when modeling HIV-1 viral load data, owing to the nonlinear trajectories and the presence of left-censored data resulting from the assay's lower limit of quantification. In this paper, for a joint model comprising a nonlinear mixed-effect model and a Cox Proportional Hazards model, we develop a computationally efficient Stochastic EM (StEM) algorithm for parameter estimation. Furthermore, we propose a novel technique for fast standard error estimation, which directly estimates standard errors from the results of StEM iterations and is broadly applicable to various joint modeling settings, such as those containing generalized linear mixed-effect models, parametric survival models, or joint models with more than two submodels. We evaluate the performance of the proposed methods through simulation studies and apply them to HIV-1 viral load data from six AIDS Clinical Trials Group studies to characterize viral rebound trajectories following the interruption of antiretroviral therapy (ART), accounting for the informative duration of off-ART periods.

Identifying optimal dosing strategies for meropenem in the paediatric intensive care unit through modelling and simulation

BACKGROUND

Meropenem, a β-lactam antibiotic commonly prescribed for severe infections, poses dosing challenges in critically ill patients due to highly variable pharmacokinetics.

OBJECTIVES

We sought to develop a population pharmacokinetic model of meropenem for critically ill paediatric and young adult patients.

PATIENTS AND METHODS

Paediatric intensive care unit patients receiving meropenem 20-40 mg/kg every 8 h as a 30 min infusion were prospectively followed for clinical data collection and scavenged opportunistic plasma sampling. Nonlinear mixed effects modelling was conducted using Monolix®. Monte Carlo simulations were performed to provide dosing recommendations against susceptible pathogens (MIC ≤ 2 mg/L).

RESULTS

Data from 48 patients, aged 1 month to 30 years, with 296 samples, were described using a two-compartment model with first-order elimination. Allometric body weight scaling accounted for body size differences. Creatinine clearance and percentage of fluid balance were identified as covariates on clearance and central volume of distribution, respectively. A maturation function for renal clearance was included. Monte Carlo simulations suggested that for a target of 40% fT > MIC, the most effective dosing regimen is 20 mg/kg every 8 h with a 3 h infusion. If higher PD targets are considered, only continuous infusion regimens ensure target attainment against susceptible pathogens, ranging from 60 mg/kg/day to 120 mg/kg/day.

CONCLUSIONS

We successfully developed a population pharmacokinetic model of meropenem using real-world data from critically ill paediatric and young adult patients with an opportunistic sampling strategy and provided dosing recommendations based on the patients' renal function and fluid status.

Tinzaparin, an alternative to subcutaneous unfractionated heparin, in patients with severe and end-stage renal impairment: a retrospective observational single-center study

BACKGROUND

Tinzaparin could be easier to manage than unfractionated heparin in patients with severe renal impairment. However, clinical and pharmacologic data regarding its use in such patients are lacking.

OBJECTIVES

The aims of this study were to determine, in patients with estimated glomerular filtration rate (eGFR) of <30 mL.min⁻1, tinzaparin pharmacokinetics (PK) parameters using a population PK approach and bleeding and thrombotic complications.

METHODS

We performed a retrospective observational single-center study, including in-patients with eGFR of <30 mL.min⁻1 receiving prophylactic (4500 IU.d⁻1) or therapeutic (175 IU.kg⁻1.d⁻1) tinzaparin. Measured anti-Xa levels were analyzed using a nonlinear mixed-effects modeling approach. Individual predicted tinzaparin exposure markers at steady state were calculated for each patient and dosing regimen. The PK was also evaluated through Monte Carlo simulations based on the final covariate model parameter estimates.

RESULTS

Over a 22-month period, 802 tinzaparin treatment periods in 623 patients were analyzed: two-thirds received a prophylactic dose, 66% had an eGFR of <20 mL.min⁻1, and 25% were on renal replacement therapy. In patients for whom anti-Xa measurements were performed (n = 199; 746 values), PK parameters, profiles, and maximum plasma concentrations were comparable with those in patients without renal impairment or in healthy volunteers. In the whole population, major bleeding occurred in 2.4% and 3.5% of patients receiving prophylactic and therapeutic doses over a median 9- and 7-day treatment period, respectively. No patients had thrombotic complications.

CONCLUSION

Tinzaparin PK parameters and profiles were not affected by renal impairment. This suggests that tinzaparin, at therapeutic or prophylactic dose, could be an alternative to unfractionated heparin in hospitalized patients with severe renal impairment.

A Population Pharmacokinetic Model of Tocilizumab in Kidney Transplant Patients Treated for Chronic Active Antibody-Mediated Rejection: Comparison of Plasma Exposure Between Intravenous and Subcutaneous Administration Schemes

BACKGROUND

Tocilizumab prevents the clinical worsening of chronic active antibody-mediated rejection (CAAMR) in kidney transplant recipients. Following a global shortage of the intravenous pharmaceutical form in 2022, patients were switched from monthly intravenous administration of 8 mg/kg to weekly subcutaneous injection of 162 mg, raising the question of bioequivalence between these schemes of administration.

AIMS

We aimed to compare the areas under the curve (AUC) of tocilizumab in virtual simulations of populations treated with the two administration schemes and to identify the covariates that could contribute to pharmacokinetic variability of tocilizumab in kidney transplant patients with CAAMR who received tocilizumab as salvage treatment.

METHODS

This retrospective monocentric study included 43 kidney transplant patients (202 tocilizumab concentrations) with CAAMR treated with intravenous or subcutaneous tocilizumab between December 2020 and January 2023. We developed a population pharmacokinetic model using nonlinear mixed effects modeling and identified the covariates that could contribute to tocilizumab AUC variability. Monte Carlo simulations were then performed to assess the subcutaneous and intravenous tocilizumab AUC for 0-28 days (M1), 56-84 days (M3), 140-168 days (M6), and 308-336 days (M12). Bioequivalence was defined by SC/IV AUC geometric mean ratios (GMRs) between 0.80 and 1.25.

RESULTS

A two-compartment model with parallel linear and nonlinear elimination best described the concentration-time data. Significant covariates for tocilizumab clearance were body weight, urinary albumin-to-creatinine ratio (ACR), and inflammation status [C-reactive protein (CRP) ≥ 5 mg/L]. The GMR values and their 90% confidence intervals at M3, M6, and M12 were within the 0.8-1.25 margin for equivalence. Conversely, the 90% prediction intervals of the GMR were much wider than the 90% confidence intervals and did not fall within 0.8 and 1.25.

CONCLUSIONS

From month 3 of treatment, the subcutaneous and intravenous tocilizumab administration schemes provided average bioequivalent pharmacokinetic exposure at the population level but not at the individual level. Body weight, inflammation, ACR, and administration scheme should be considered to personalize the dose of tocilizumab for patients with CAAMR. Further studies are required to determine the target of tocilizumab exposure in kidney transplant patients with CAAMR.

Assessment of body mass-related covariates for rifampicin pharmacokinetics in healthy Caucasian volunteers

PURPOSE

Currently, body weight-based dosing of rifampicin is recommended. But lately, fat-free mass (FFM) was reported to be superior to body weight (BW). The present evaluation aimed to assess the influence of body mass-related covariates on rifampicin's pharmacokinetics (PK) parameters in more detail using non-linear mixed effects modeling (NLMEM).

METHODS

Twenty-four healthy Caucasian volunteers were enrolled in a bioequivalence study, each receiving a test and a reference tablet of 600 mg of rifampicin separated by a wash-out period of at least 9 days. Monolix version 2023R1 was used for NLMEM. Monte Carlo simulations (MCS) were performed to visualize the relationship of body size descriptors to the exposure to rifampicin.

RESULTS

A one-compartment model with nonlinear (Michaelis-Menten) elimination and zero-order absorption kinetics with a lag time best described the data. The covariate model including fat-free mass (FFM) on volume of distribution (V/F) and on maximum elimination rate (Vmax/F) lowered the objective function value (OFV) by 56.4. The second-best covariate model of sex on V/F and Vmax/F and BW on V/F reduced the OFV by 51.2. The decrease in unexplained inter-individual variability on Vmax/F in both covariate models was similar. For a given dose, MCS showed lower exposure to rifampicin with higher FFM and accordingly in males compared to females with the same BW and body height.

CONCLUSION

Our results indicate that beyond BW, body composition as reflected by FFM could also be relevant for optimized dosing of rifampicin. This assumption needs to be studied further in patients treated with rifampicin.

Modeling COVID-19 vaccine booster-elicited antibody response and impact of infection history

The 3-dose COVID-19 vaccine (booster vaccination) has been offered worldwide. As booster vaccinations continue, it is important to understand the antibody dynamics elicited by booster vaccination in order to evaluate and develop vaccination needs and strategies. Here, we investigated longitudinal data by monitoring IgG antibodies against the receptor binding domain (RBD) in health care workers. We extended our previously developed mathematical model to booster vaccines and successfully fitted antibody titers over time in the absence and presence of past SARS-CoV-2 infection. Quantitative analysis using our mathematical model indicated that anti-RBD IgG titers increase to a comparable extent after booster vaccination, regardless of the presence or absence of infection, but infection history extends the duration of antibody response by 1.28 times. Such a mathematical modeling approach can be used to inform future vaccination strategies on the basis of an individual's immune history. Our simple quantitative approach can be extended to any kind of vaccination and therefore can form a basis for policy decisions regarding the distribution of booster vaccines to strengthen immunity in future pandemics.

Nonlinear Mixed Effects Modeling of Glucagon Kinetics in Healthy Subjects

OBJECTIVE

To date, the lack of a model of glucagon kinetics precluded the possibility of estimating and studying glucagon secretion in vivo, e.g., using deconvolution, as done for other hormones like insulin and C-peptide. Here, we used a nonlinear mixed effects technique to develop a robust population model of glucagon kinetics, able to describe both the typical population kinetics (TPK) and the between-subject variability (BSV), and relate this last to easily measurable subject characteristics.

METHODS

Thirty-four models of increasing complexity (variably including covariates and correlations among random effects) were identified on glucagon profiles obtained from 53 healthy subjects, who received a constant infusion of somatostatin to suppress endogenous glucagon production, followed by a continuous infusion of glucagon (65 ng/kg/min). Model selection was performed based on its ability to fit the data, provide precise parameter estimates, and parsimony criteria.

RESULTS

A two-compartment model was the most parsimonious. The model was able to accurately describe both the TPK and the BSV of model parameters as function of body mass and body surface area. Parameters were precisely estimated, with central volume of distribution V1 = 5.46 L and peripheral volume of distribution V2 = 5.51 L. The introduction of covariates resulted in a significant shrinkage of the unexplained BSV and considerably improved the model fit.

CONCLUSION

We developed a robust population model of glucagon kinetics.

SIGNIFICANCE

This model provides a deeper understanding of glucagon kinetics and is usable to estimate glucagon secretion in vivo by deconvolution of plasma glucagon concentration data.

Modeling the kinetics of the neutralizing antibody response against SARS-CoV-2 variants after several administrations of Bnt162b2

Because SARS-CoV-2 constantly mutates to escape from the immune response, there is a reduction of neutralizing capacity of antibodies initially targeting the historical strain against emerging Variants of Concern (VoC)s. That is why the measure of the protection conferred by vaccination cannot solely rely on the antibody levels, but also requires to measure their neutralization capacity. Here we used a mathematical model to follow the humoral response in 26 individuals that received up to three vaccination doses of Bnt162b2 vaccine, and for whom both anti-S IgG and neutralization capacity was measured longitudinally against all main VoCs. Our model could identify two independent mechanisms that led to a marked increase in measured humoral response over the successive vaccination doses. In addition to the already known increase in IgG levels after each dose, we identified that the neutralization capacity was significantly increased after the third vaccine administration against all VoCs, despite large inter-individual variability. Consequently, the model projects that the mean duration of detectable neutralizing capacity against non-Omicron VoC is between 348 days (Beta variant, 95% Prediction Intervals PI [307; 389]) and 587 days (Alpha variant, 95% PI [537; 636]). Despite the low neutralization levels after three doses, the mean duration of detectable neutralizing capacity against Omicron variants varies between 173 days (BA.5 variant, 95% PI [142; 200]) and 256 days (BA.1 variant, 95% PI [227; 286]). Our model shows the benefit of incorporating the neutralization capacity in the follow-up of patients to better inform on their level of protection against the different SARS-CoV-2 variants. Trial registration: This clinical trial is registered with ClinicalTrials.gov, Trial IDs NCT04750720 and NCT05315583.

Smoothed simulated pseudo-maximum likelihood estimation for nonlinear mixed effects models with censored responses

Nonlinear mixed effects models have been widely applied to analyses of data that arise from biological, agricultural, and environmental sciences. Estimation of and inference on parameters in nonlinear mixed effects models are often based on the specification of a likelihood function. Maximizing this likelihood function can be complicated by the specification of the random effects distribution, especially in the presence of multiple random effects. The implementation of nonlinear mixed effects models can be further complicated by left-censored responses, representing measurements from bioassays where the exact quantification below a certain threshold is not possible. Motivated by the need to characterize the nonlinear human immunodeficiency virus RNA viral load trajectories after the interruption of antiretroviral therapy, we propose a smoothed simulated pseudo-maximum likelihood estimation approach to fit nonlinear mixed effects models in the presence of left-censored observations. We establish the consistency and asymptotic normality of the resulting estimators. We develop testing procedures for the correlation among random effects and for testing the distributional assumptions on random effects against a specific alternative. In contrast to the existing variants of expectation-maximization approaches, the proposed methods offer flexibility in the specification of the random effects distribution and convenience in making inference about higher-order correlation parameters. We evaluate the finite-sample performance of the proposed methods through extensive simulation studies and illustrate them on a combined dataset from six AIDS Clinical Trials Group treatment interruption studies.

Evaluation of model-based bioequivalence approach for single sample pharmacokinetic studies

In a traditional pharmacokinetic (PK) bioequivalence (BE) study, a two-way crossover study is conducted, PK parameters (namely the area under the time-concentration curve [AUC] and the maximal concentration [C max ]) are obtained by noncompartmental analysis (NCA), and the BE analysis is performed using the two one-sided test (TOST) method. For ophthalmic drugs, however, only one sample of aqueous humor, in one eye, per eye can be obtained in each patient, which precludes the traditional BE analysis. To circumvent this issue, the U.S. Food and Drug Administration (FDA) has proposed an approach coupling NCA with either parametric or nonparametric bootstrap (NCA bootstrap). The model-based TOST (MB-TOST) has previously been proposed and evaluated successfully for various settings of sparse PK BE studies. In this paper, we evaluate, via simulations, MB-TOST in the specific setting of single sample PK BE study and compare its performance to NCA bootstrap. We performed BE study simulations using a published PK model and parameter values and evaluated multiple scenarios, including study design (parallel or crossover), sampling times (5 or 10 spread across the dosing interval), and geometric mean ratio (of 0.8, 0.9, 1, and 1.25). Using the simulated structural PK model, MB-TOST performed similarly to NCA bootstrap for AUC. ForC max , the latter tended to be conservative and less powerful. Our research suggests that MB-TOST may be considered as an alternative BE approach for single sample PK studies, provided that the PK model is correctly specified and the test drug has the same structural model as the reference drug.

A French Real-World Evidence Study Evaluating the Efficacy, Safety, and Pharmacokinetic Parameters of rVIII-SingleChain in Patients with Hemophilia A Receiving Prophylaxis

BACKGROUND

 rVIII-SingleChain is a recombinant factor VIII (FVIII) with increased binding affinity to von Willebrand factor compared with other FVIII products. rVIII-SingleChain is indicated for the treatment and prevention of bleeding episodes in patients with hemophilia A.

OBJECTIVES

 To collect real-world evidence data from patients treated with rVIII-SingleChain to confirm the efficacy and safety established in the clinical trial program and carry out a population pharmacokinetic (PK) analysis.

PATIENTS/METHODS

 This interim analysis includes data, collected between January 2018 - September 2021, from patients treated with rVIII-SingleChain prophylaxis at French Hemophilia Treatment centers. Data on annualized bleeding rates, dosing frequency, and consumption before and after switching to rVIII-SingleChain were recorded. A population PK analysis was also conducted to estimate PK parameters.

RESULTS

 Overall, 43 patients switched to prophylaxis with rVIII-SingleChain either from a previous prophylaxis regimen or from on-demand treatment. Following the switch to rVIII-SingleChain, patients maintained excellent bleed control. After switching to rVIII-SingleChain, most patients maintained or reduced their regimen. Interestingly, a majority of patients treated >2 ×/weekly with a standard half-life FVIII reduced both injection frequency and FVIII consumption with rVIII-SingleChain. A PK analysis revealed a lower clearance of rVIII-SingleChain (1.9 vs. 2.1 dL/h) and a longer half-life both in adolescents/adults (n = 28) and pediatric (n = 6) patients (15.5 and 11.9 hours, respectively vs. 14.5 and 10.3 hours) than previously reported.

CONCLUSIONS

 Patients who switched to rVIII-SingleChain prophylaxis demonstrated excellent bleed control and a reduction in infusion frequency. A population PK analysis revealed improved PK parameters compared with those reported in the clinical trial.

Population pharmacokinetics and pharmacodynamics of imipenem in neutropenic adult patients

OBJECTIVE

Imipenem is recommended in patients with chemotherapy-induced febrile neutropenia. Although alterations of antibiotic pharmacokinetic parameters have been reported in such patients, little data is available on imipenem.

METHODS

Prospective, single-center, non-interventional pharmacokinetic cohort study in adults with chemotherapy-induced febrile neutropenia. Critically ill patients were excluded. Imipenem was administered as a 30-min infusion of 1000 mg/8h. Total imipenem plasma concentrations were assayed by high-performance liquid chromatography during neutropenia and just after neutrophil recovery. We estimated population pharmacokinetic parameters of imipenem by non-linear mixed-effect modelling using the SAEM algorithm.

RESULTS

Sixteen patients were included in the study, including nine women (56.3%), median age 37 years (range, 18.3; 78.3). Eight patients had an hematological malignancy (50.0%) and seven had a solid tumor (43.8%). Imipenem pharmacokinetics were best described by a one-compartment model with first-order elimination. Mean values for imipenem were: clearance 14.3L/h and 10.9L/h and volume of distribution 20.7L and 14.5 L during neutropenia and after recovery, respectively. Imipenem plasma area under the curve at steady state was reduced by 23% during neutropenia. However, all patients achieved a pharmacodynamic target of %fT_>MIC ≥ 40% with a regimen of 1000 mg/8 h or 500 mg/6 h, for MICs up to 2 mg/L. The pharmacodynamics profile for a target of %fT > MIC = 100% was however less favorable with 500 mg/6 h or 1000 mg/8 h either during or after neutropenia.

CONCLUSION

Pharmacokinetic/pharmacodynamic goals for imipenem were similar in patients during and after neutropenia, despite reduced plasma exposure.

A simple model for viral decay dynamics and the distribution of infected cell life spans in SHIV-infected infant rhesus macaques

The dynamics of HIV viral load following the initiation of antiretroviral therapy is not well-described by simple, single-phase exponential decay. Several mathematical models have been proposed to describe its more complex behavior, the most popular of which is two-phase exponential decay. The underlying assumption in two-phase exponential decay is that there are two classes of infected cells with different lifespans. However, with the exception of CD4+ T cells, there is not a consensus on all of the cell types that can become productively infected, and the fit of the two-phase exponential decay to observed data from SHIV.C.CH505 infected infant rhesus macaques was relatively poor. Therefore, we propose a new model for viral decay, inspired by the Gompertz model where the decay rate itself is a dynamic variable. We modify the Gompertz model to include a linear term that modulates the decay rate. We show that this simple model performs as well as the two-phase exponential decay model on HIV and SIV data sets, and outperforms it for the infant rhesus macaque SHIV.C.CH505 infection data set. We also show that by using a stochastic differential equation formulation, the modified Gompertz model can be interpreted as being driven by a population of infected cells with a continuous distribution of cell lifespans, and estimate this distribution for the SHIV.C.CH505-infected infant rhesus macaques. Thus, we find that the dynamics of viral decay in this model of infant HIV infection and treatment may be explained by a distribution of cell lifespans, rather than two distinct cell types.

Estimation of timing of infection from longitudinal SARS-CoV-2 viral load data: mathematical modelling study

BACKGROUND

Multiple waves of the COVID-19 epidemic have hit most countries by the end of 2021. Most of those waves are caused by emergence and importation of new variants. To prevent importation of new variants, combination of border control and contact tracing is essential. However, the timing of infection inferred by interview is influenced by recall bias and hinders the contact tracing process.

METHODS

We propose a novel approach to infer the timing of infection, by employing a within-host model to capture viral load dynamics after the onset of symptoms. We applied this approach to ascertain secondary transmission which can trigger outbreaks. As a demonstration, the 12 initial reported cases in Singapore, which were considered as imported because of their recent travel history to Wuhan, were analyzed to assess whether they are truly imported.

RESULTS

Our approach suggested that 6 cases were infected prior to the arrival in Singapore, whereas other 6 cases might have been secondary local infection. Three among the 6 potential secondary transmission cases revealed that they had contact history to previously confirmed cases.

CONCLUSIONS

Contact trace combined with our approach using viral load data could be the key to mitigate the risk of importation of new variants by identifying cases as early as possible and inferring the timing of infection with high accuracy.

Modeling the bacterial dynamics in the gut microbiota following an antibiotic‐induced perturbation

Recent studies have highlighted the importance of ecological interactions in dysbiosis of gut microbiota, but few focused on their role in antibiotic‐induced perturbations. We used the data from the CEREMI trial in which 22 healthy volunteers received a 3‐day course of ceftriaxone or cefotaxime antibiotics. Fecal samples were analyzed by 16S rRNA gene profiling, and the total bacterial counts were determined in each sample by flux cytometry. As the gut exposure to antibiotics could not be experimentally measured despite a marked impact on the gut microbiota, it was reconstructed using the counts of susceptible Escherichia coli. The dynamics of absolute counts of bacterial families were analyzed using a generalized Lotka–Volterra equations and nonlinear mixed effect modeling. Bacterial interactions were studied using a stepwise approach. Two negative and three positive interactions were identified. Introducing bacterial interactions in the modeling approach better fitted the data, and provided different estimates of antibiotic effects on each bacterial family than a simple model without interaction. The time to return to 95% of the baseline counts was significantly longer in ceftriaxone‐treated individuals than in cefotaxime‐treated subjects for two bacterial families: Akkermansiaceae (median [range]: 11.3 days [0; 180.0] vs. 4.2 days [0; 25.6], p = 0.027) and Tannerellaceae (13.7 days [6.1; 180.0] vs. 6.2 days [5.4; 17.3], p = 0.003). Taking bacterial interaction as well as individual antibiotic exposure profile into account improves the analysis of antibiotic‐induced dysbiosis.

Analysis of the Equilibrium Phase in Immune-Controlled Tumors Provides Hints for Designing Better Strategies for Cancer Treatment

When it comes to improving cancer therapies, one challenge is to identify key biological parameters that prevent immune escape and maintain an equilibrium state characterized by a stable subclinical tumor mass, controlled by the immune cells. Based on a space and size structured partial differential equation model, we developed numerical methods that allow us to predict the shape of the equilibrium at low cost, without running simulations of the initial-boundary value problem. In turn, the computation of the equilibrium state allowed us to apply global sensitivity analysis methods that assess which and how parameters influence the residual tumor mass. This analysis reveals that the elimination rate of tumor cells by immune cells far exceeds the influence of the other parameters on the equilibrium size of the tumor. Moreover, combining parameters that sustain and strengthen the antitumor immune response also proves more efficient at maintaining the tumor in a long-lasting equilibrium state. Applied to the biological parameters that define each type of cancer, such numerical investigations can provide hints for the design and optimization of cancer treatments.

A workflow for the joint modeling of longitudinal and event data in the development of therapeutics: Tools, statistical methods, and diagnostics

Clinical trials investigate treatment endpoints that usually include measurements of pharmacodynamic and efficacy biomarkers in early-phase studies and patient-reported outcomes as well as event risks or rates in late-phase studies. In recent years, a systematic trend in clinical trial data analytics and modeling has been observed, where retrospective data are integrated into a quantitative framework to prospectively support analyses of interim data and design of ongoing and future studies of novel therapeutics. Joint modeling is an advanced statistical methodology that allows for the investigation of clinical trial outcomes by quantifying the association between baseline and/or longitudinal biomarkers and event risk. Using an exemplar data set from non-small cell lung cancer studies, we propose and test a workflow for joint modeling. It allows a modeling scientist to comprehensively explore the data, build survival models, investigate goodness-of-fit, and subsequently perform outcome predictions using interim biomarker data from an ongoing study. The workflow illustrates a full process, from data exploration to predictive simulations, for selected multivariate linear and nonlinear mixed-effects models and software tools in an integrative and exhaustive manner.

Complex decay dynamics of HIV virions, intact and defective proviruses, and 2LTR circles following initiation of antiretroviral therapy

In persons living with HIV-1 (PLWH) who start antiretroviral therapy (ART), plasma virus decays in a biphasic fashion to below the detection limit. The first phase reflects the short half-life (<1 d) of cells that produce most of the plasma virus. The second phase represents the slower turnover (t1/2 = 14 d) of another infected cell population, whose identity is unclear. Using the intact proviral DNA assay (IPDA) to distinguish intact and defective proviruses, we analyzed viral decay in 17 PLWH initiating ART. Circulating CD4+ T cells with intact proviruses include few of the rapidly decaying first-phase cells. Instead, this population initially decays more slowly (t1/2 = 12.9 d) in a process that largely represents death or exit from the circulation rather than transition to latency. This more protracted decay potentially allows for immune selection. After ∼3 mo, the decay slope changes, and CD4+ T cells with intact proviruses decay with a half-life of 19 mo, which is still shorter than that of the latently infected cells that persist on long-term ART. Two-long-terminal repeat (2LTR) circles decay with fast and slow phases paralleling intact proviruses, a finding that precludes their use as a simple marker of ongoing viral replication. Proviruses with defects at the 5' or 3' end of the genome show equivalent monophasic decay at rates that vary among individuals. Understanding these complex early decay processes is important for correct use of reservoir assays and may provide insights into properties of surviving cells that can constitute the stable latent reservoir.

Nonlinear mixed-effects models for HIV viral load trajectories before and after antiretroviral therapy interruption, incorporating left censoring

OBJECTIVES

Characterizing features of the viral rebound trajectories and identifying host, virological, and immunological factors that are predictive of the viral rebound trajectories are central to HIV cure research. We investigate if key features of HIV viral decay and CD4 trajectories during antiretroviral therapy (ART) are associated with characteristics of HIV viral rebound following ART interruption.

METHODS

Nonlinear mixed effect (NLME) models are used to model viral load trajectories before and following ART interruption, incorporating left censoring due to lower detection limits of viral load assays. A stochastic approximation EM (SAEM) algorithm is used for parameter estimation and inference. To circumvent the computational intensity associated with maximizing the joint likelihood, we propose an easy-to-implement three-step method.

RESULTS

We evaluate the performance of the proposed method through simulation studies and apply it to data from the Zurich Primary HIV Infection Study. We find that some key features of viral load during ART (e.g., viral decay rate) are significantly associated with important characteristics of viral rebound following ART interruption (e.g., viral set point).

CONCLUSIONS

The proposed three-step method works well. We have shown that key features of viral decay during ART may be associated with important features of viral rebound following ART interruption.

Mathematical Modeling of Within-Host, Untreated, Cytomegalovirus Infection Dynamics after Allogeneic Transplantation

Cytomegalovirus (CMV) causes significant morbidity and mortality in recipients of allogeneic hematopoietic cell transplantation (HCT). Whereas insights gained from mathematical modeling of other chronic viral infections such as HIV, hepatitis C, and herpes simplex virus-2 have aided in optimizing therapy, previous CMV modeling has been hindered by a lack of comprehensive quantitative PCR viral load data from untreated episodes of viremia in HCT recipients. We performed quantitative CMV DNA PCR on stored, frozen serum samples from the placebo group of participants in a historic randomized controlled trial of ganciclovir for the early treatment of CMV infection in bone marrow transplant recipients. We developed four main ordinary differential Equation mathematical models and used model selection theory to choose between 38 competing versions of these models. Models were fit using a population, nonlinear, mixed-effects approach. We found that CMV kinetics from untreated HCT recipients are highly variable. The models that recapitulated the observed patterns most parsimoniously included explicit, dynamic immune cell compartments and did not include dynamic target cell compartments, consistent with the large number of tissue and cell types that CMV infects. In addition, in our best-fitting models, viral clearance was extremely slow, suggesting severe impairment of the immune response after HCT. Parameters from our best model correlated well with participants’ clinical risk factors and outcomes from the trial, further validating our model. Our models suggest that CMV dynamics in HCT recipients are determined by host immune response rather than target cell limitation in the absence of antiviral treatment.

Incomplete antiviral treatment may induce longer durations of viral shedding during SARS-CoV-2 infection

The duration of viral shedding is determined by a balance between de novo infection and removal of infected cells. That is, if infection is completely blocked with antiviral drugs (100% inhibition), the duration of viral shedding is minimal and is determined by the length of virus production. However, some mathematical models predict that if infected individuals are treated with antiviral drugs with efficacy below 100%, viral shedding may last longer than without treatment because further de novo infections are driven by entry of the virus into partially protected, uninfected cells at a slower rate. Using a simple mathematical model, we quantified SARS-CoV-2 infection dynamics in non-human primates and characterized the kinetics of viral shedding. We counterintuitively found that treatments initiated early, such as 0.5 d after virus inoculation, with intermediate to relatively high efficacy (30-70% inhibition of virus replication) yield a prolonged duration of viral shedding (by about 6.0 d) compared with no treatment.

Revisiting the guidelines for ending isolation for COVID-19 patients

Since the start of the COVID-19 pandemic, two mainstream guidelines for defining when to end the isolation of SARS-CoV-2-infected individuals have been in use: the one-size-fits-all approach (i.e. patients are isolated for a fixed number of days) and the personalized approach (i.e. based on repeated testing of isolated patients). We use a mathematical framework to model within-host viral dynamics and test different criteria for ending isolation. By considering a fixed time of 10 days since symptom onset as the criterion for ending isolation, we estimated that the risk of releasing an individual who is still infectious is low (0-6.6%). However, this policy entails lengthy unnecessary isolations (4.8-8.3 days). In contrast, by using a personalized strategy, similar low risks can be reached with shorter prolonged isolations. The obtained findings provide a scientific rationale for policies on ending the isolation of SARS-CoV-2-infected individuals.

Pharmacokinetics of enoxaparin in COVID-19 critically ill patients

Background

In intensive-care unit (ICU) patients, pathophysiological changes may affect the pharmacokinetics of enoxaparin and result in underdosing.

Objectives

To develop a pharmacokinetic model of enoxaparin to predict the time-exposure profiles of various thromboprophylactic regimens in COVID-19 ICU-patients.

Methods

This was a retrospective study in ICUs of two French hospitals. Anti-Xa activities from consecutive patients with laboratory-confirmed SARS-CoV-2 infection treated with enoxaparin for the prevention or the treatment of venous thrombosis were used to develop a population pharmacokinetic model using non-linear mixed effects techniques. Monte Carlo simulations were then performed to predict enoxaparin exposure at steady-state after three days of administration.

Results

A total of 391 anti-Xa samples were measured in 95 patients. A one-compartment model with first-order kinetics best fitted the data. The covariate analysis showed that enoxaparin clearance (typical value 1.1 L.h-1) was related to renal function estimated by the CKD-EPI formula and volume of distribution (typical value 17.9 L) to actual body weight.

Simulation of anti-Xa activities with enoxaparin 40 mg qd indicated that 64% of the patients had peak levels within the range 0.2 to 0.5 IU.mL-1 and 75% had 12-hour levels above 0.1 IU.mL-1. Administration of a total daily dose of at least 60 mg per day improved the probability of target attainment.

Conclusion

In ICU COVID-19 patients, exposure to enoxaparin is reduced due to an increase in the volume of distribution and clearance. Consequently, enoxaparin 40 mg qd is suboptimal to attain thromboprophylactic anti-Xa levels.

Detection of significant antiviral drug effects on COVID-19 with reasonable sample sizes in randomized controlled trials: A modeling study

BACKGROUND

Development of an effective antiviral drug for Coronavirus Disease 2019 (COVID-19) is a global health priority. Although several candidate drugs have been identified through in vitro and in vivo models, consistent and compelling evidence from clinical studies is limited. The lack of evidence from clinical trials may stem in part from the imperfect design of the trials. We investigated how clinical trials for antivirals need to be designed, especially focusing on the sample size in randomized controlled trials.

METHODS AND FINDINGS

A modeling study was conducted to help understand the reasons behind inconsistent clinical trial findings and to design better clinical trials. We first analyzed longitudinal viral load data for Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) without antiviral treatment by use of a within-host virus dynamics model. The fitted viral load was categorized into 3 different groups by a clustering approach. Comparison of the estimated parameters showed that the 3 distinct groups were characterized by different virus decay rates (p-value < 0.001). The mean decay rates were 1.17 d-1 (95% CI: 1.06 to 1.27 d-1), 0.777 d-1 (0.716 to 0.838 d-1), and 0.450 d-1 (0.378 to 0.522 d-1) for the 3 groups, respectively. Such heterogeneity in virus dynamics could be a confounding variable if it is associated with treatment allocation in compassionate use programs (i.e., observational studies). Subsequently, we mimicked randomized controlled trials of antivirals by simulation. An antiviral effect causing a 95% to 99% reduction in viral replication was added to the model. To be realistic, we assumed that randomization and treatment are initiated with some time lag after symptom onset. Using the duration of virus shedding as an outcome, the sample size to detect a statistically significant mean difference between the treatment and placebo groups (1:1 allocation) was 13,603 and 11,670 (when the antiviral effect was 95% and 99%, respectively) per group if all patients are enrolled regardless of timing of randomization. The sample size was reduced to 584 and 458 (when the antiviral effect was 95% and 99%, respectively) if only patients who are treated within 1 day of symptom onset are enrolled. We confirmed the sample size was similarly reduced when using cumulative viral load in log scale as an outcome. We used a conventional virus dynamics model, which may not fully reflect the detailed mechanisms of viral dynamics of SARS-CoV-2. The model needs to be calibrated in terms of both parameter settings and model structure, which would yield more reliable sample size calculation.

CONCLUSIONS

In this study, we found that estimated association in observational studies can be biased due to large heterogeneity in viral dynamics among infected individuals, and statistically significant effect in randomized controlled trials may be difficult to be detected due to small sample size. The sample size can be dramatically reduced by recruiting patients immediately after developing symptoms. We believe this is the first study investigated the study design of clinical trials for antiviral treatment using the viral dynamics model.

Population pharmacokinetics and dosing simulations of amoxicillin in obese adults receiving co-amoxiclav

BACKGROUND

Pneumonia, skin and soft tissue infections are more frequent in obese patients and are most often treated by co-amoxiclav, using similar dosing regimens to those used for non-obese subjects. No data are available on amoxicillin pharmacokinetics among obese subjects receiving co-amoxiclav.

MATERIALS AND METHODS

Prospective, single-centre, open-label, non-randomized, crossover pharmacokinetic trial having enrolled obese otherwise healthy adult subjects. A first dose of co-amoxiclav (amoxicillin/clavulanate 1000/200 mg) was infused IV over 30 min, followed by a second dose (1000/125 mg) administered orally, separated by a washout period of ≥24 h. We assayed concentrations of amoxicillin by a validated ultra HPLC-tandem MS technique. We estimated population pharmacokinetic parameters of amoxicillin by non-linear mixed-effect modelling using the SAEM algorithm developed by Monolix.

RESULTS

Twenty-seven subjects were included in the IV study, with 24 included in the oral part of the study. Median body weight and BMI were 109.3 kg and 40.6 kg/m2, respectively. Amoxicillin pharmacokinetics were best described by a two-compartment model with first-order elimination. Mean values for clearance, central volume, intercompartmental clearance and peripheral volume were, respectively, 14.6 L/h, 9.0 L, 4.2 L/h and 6.4 L for amoxicillin. Oral bioavailability of amoxicillin was 79.7%. Amoxicillin Cmax after oral administration significantly reduced with weight (P = 0.013). Dosing simulations for amoxicillin predicted that most of the population will achieve the pharmacodynamic target of fT>MIC ≥40% with the regimen of co-amoxiclav 1000/200 mg (IV) or 1000/125 mg (oral) q8h for MICs titrated up to 0.5 mg/L (IV) and 1 mg/L (oral).

CONCLUSIONS

Pharmacokinetic/pharmacodynamic goals for amoxicillin can be obtained in obese subjects.

FDA-Industry Scientific Exchange on assessing quantitative systems pharmacology models in clinical drug development: a meeting report, summary of challenges/gaps, and future perspective

The pharmaceutical industry is actively applying quantitative systems pharmacology (QSP) to make internal decisions and guide drug development. To facilitate the eventual development of a common framework for assessing the credibility of QSP models for clinical drug development, scientists from US Food and Drug Administration and the pharmaceutical industry organized a full-day virtual Scientific Exchange on July 1, 2020. An assessment form was used to ensure consistency in the evaluation process. Among the cases presented, QSP was applied to various therapeutic areas. Applications mostly focused on phase 2 dose selection. Model transparency, including details on expert knowledge and data used for model development, was identified as a major factor for robust model assessment. The case studies demonstrated some commonalities in the workflow of QSP model development, calibration, and validation but differ in the size, scope, and complexity of QSP models, in the acceptance criteria for model calibration and validation, and in the algorithms/approaches used for creating virtual patient populations. Though efforts are being made to build the credibility of QSP models and the confidence is increasing in applying QSP for internal decisions at the clinical stages of drug development, there are still many challenges facing QSP application to late stage drug development. The QSP community needs a strategic plan that includes the ability and flexibility to Adapt, to establish Common expectations for model Credibility needed to inform drug Labeling and patient care, and to AIM to achieve the goal (ACCLAIM).

Time variation in the probability of failing to detect a case of polymerase chain reaction testing for SARS-CoV-2 as estimated from a viral dynamics model

Viral tests including polymerase chain reaction (PCR) tests are recommended to diagnose COVID-19 infection during the acute phase of infection. A test should have high sensitivity; however, the sensitivity of the PCR test is highly influenced by viral load, which changes over time. Because it is difficult to collect data before the onset of symptoms, the current literature on the sensitivity of the PCR test before symptom onset is limited. In this study, we used a viral dynamics model to track the probability of failing to detect a case of PCR testing over time, including the presymptomatic period. The model was parametrized by using longitudinal viral load data collected from 30 hospitalized patients. The probability of failing to detect a case decreased toward symptom onset, and the lowest probability was observed 2 days after symptom onset and increased afterwards. The probability on the day of symptom onset was 1.0% (95% CI: 0.5 to 1.9) and that 2 days before symptom onset was 60.2% (95% CI: 57.1 to 63.2). Our study suggests that the diagnosis of COVID-19 by PCR testing should be done carefully, especially when the test is performed before or way after symptom onset. Further study is needed of patient groups with potentially different viral dynamics, such as asymptomatic cases.

SARS-CoV-2 viral dynamics in non-human primates

Non-human primates infected with SARS-CoV-2 exhibit mild clinical signs. Here we used a mathematical model to characterize in detail the viral dynamics in 31 cynomolgus macaques for which nasopharyngeal and tracheal viral load were frequently assessed. We identified that infected cells had a large burst size (>104 virus) and a within-host reproductive basic number of approximately 6 and 4 in nasopharyngeal and tracheal compartment, respectively. After peak viral load, infected cells were rapidly lost with a half-life of 9 hours, with no significant association between cytokine elevation and clearance, leading to a median time to viral clearance of 10 days, consistent with observations in mild human infections. Given these parameter estimates, we predict that a prophylactic treatment blocking 90% of viral production or viral infection could prevent viral growth. In conclusion, our results provide estimates of SARS-CoV-2 viral kinetic parameters in an experimental model of mild infection and they provide means to assess the efficacy of future antiviral treatments.

A quantitative model used to compare within-host SARS-CoV-2, MERS-CoV, and SARS-CoV dynamics provides insights into the pathogenesis and treatment of SARS-CoV-2

The scientific community is focused on developing antiviral therapies to mitigate the impacts of the ongoing novel coronavirus disease 2019 (COVID-19) outbreak. This will be facilitated by improved understanding of viral dynamics within infected hosts. Here, using a mathematical model in combination with published viral load data, we compare within-host viral dynamics of SARS-CoV-2 with analogous dynamics of MERS-CoV and SARS-CoV. Our quantitative analyses using a mathematical model revealed that the within-host reproduction number at symptom onset of SARS-CoV-2 was statistically significantly larger than that of MERS-CoV and similar to that of SARS-CoV. In addition, the time from symptom onset to the viral load peak for SARS-CoV-2 infection was shorter than those of MERS-CoV and SARS-CoV. These findings suggest the difficulty of controlling SARS-CoV-2 infection by antivirals. We further used the viral dynamics model to predict the efficacy of potential antiviral drugs that have different modes of action. The efficacy was measured by the reduction in the viral load area under the curve (AUC). Our results indicate that therapies that block de novo infection or virus production are likely to be effective if and only if initiated before the viral load peak (which appears 2-3 days after symptom onset), but therapies that promote cytotoxicity of infected cells are likely to have effects with less sensitivity to the timing of treatment initiation. Furthermore, combining a therapy that promotes cytotoxicity and one that blocks de novo infection or virus production synergistically reduces the AUC with early treatment. Our unique modeling approach provides insights into the pathogenesis of SARS-CoV-2 and may be useful for development of antiviral therapies.

Population Pharmacokinetics of Intravenous Ganciclovir and Oral Valganciclovir in a Pediatric Population To Optimize Dosing Regimens

Ganciclovir is indicated for curative or preventive treatment of cytomegalovirus (CMV) infections. This study aimed to characterize ganciclovir pharmacokinetics, following intravenous ganciclovir and oral valganciclovir administration, to optimize dosing schemes. All children aged <18 years receiving ganciclovir or valganciclovir were included in this study. Pharmacokinetics were described using nonlinear mixed-effect modeling. Monte Carlo simulations were used to optimize the dosing regimen to maintain the area under the concentration-time curve (AUC) in the preventive or therapeutic target. Among the 105 children (374 concentration-time observations) included, 78 received intravenous (i.v.) ganciclovir, 19 received oral valganciclovir, and 6 received both drugs. A two-compartment model with first-order absorption for valganciclovir and first-order elimination best described the data. An allometric model was used to describe the bodyweight (BW) effect. Estimated glomerular filtration rate (eGFR) and medical status of critically ill children were significantly associated with ganciclovir elimination. Recommended doses were adapted for prophylactic treatment. To obtain a therapeutic exposure, doses should be increased to 40 mg/kg of body weight/day oral or 15 to 20 mg/kg/day i.v. in children with normal eGFR and to 56 mg/kg/day oral or 20 to 25 mg/kg/day i.v. in children with augmented eGFR. These doses should be prospectively confirmed, and therapeutic drug monitoring could be used to refine them individually. (This study has been registered at ClinicalTrials.gov under identifier NCT02539407.).

What drives the dynamics of HBV RNA during treatment?

Hepatitis B virus RNA (HBV RNA)-containing particles are encapsidated pre-genomic RNA (pgRNA) detectable in chronically infected patients in addition to virions (HBV DNA) that have been suggested as a marker of the treatment efficacy. This makes promising the use of core protein allosteric modulators, such as RG7907, which disrupt the nucleocapsid assembly and profoundly reduce HBV RNA. Here, we developed a multiscale model of HBV extending the standard viral dynamic models to analyse the kinetics of HBV DNA and HBV RNA in 35 patients treated with RG7907 for 28 days. We compare the predictions with those obtained in patients treated with the nucleotide analog tenofovir. RG7907 blocked 99.3% of pgRNA encapsidation (range: 92.1%-99.9%) which led to a decline of both HBV DNA and HBV RNA. As a consequence of its mode of action, the first phase of decline of HBV RNA was rapid, uncovering the clearance of viral particles with half-life of 45 min. In contrast, HBV DNA decline was predicted to be less rapid, due to the continuous secretion of already formed viral capsids (t_1/2  = 17 ± 6 h). After few days, both markers declined at the same rate, which was attributed to the loss of infected cells (t_1/2  ≅ 6 ± 0.8 days). By blocking efficiently RNA reverse transcription but not its encapsidation, nucleotide analog in contrast was predicted to lead to a transient accumulation of HBV RNA both intracellularly and extracellularly. The model brings a conceptual framework for understanding the differences between HBV DNA and HBV RNA dynamics. Integration of HBV RNA in viral dynamic models may be helpful to better quantify the treatment effect, especially in viral-suppressed patients where HBV DNA is no longer detectable.

Zika virus dynamics: Effects of inoculum dose, the innate immune response and viral interference

Experimental Zika virus infection in non-human primates results in acute viral load dynamics that can be well-described by mathematical models. The inoculum dose that would be received in a natural infection setting is likely lower than the experimental infections and how this difference affects the viral dynamics and immune response is unclear. Here we study a dataset of experimental infection of non-human primates with a range of doses of Zika virus. We develop new models of infection incorporating both an innate immune response and viral interference with that response. We find that such a model explains the data better than models with no interaction between virus and the immune response. We also find that larger inoculum doses lead to faster dynamics of infection, but approximately the same total amount of viral production.

The relationship between the infecting dose of a pathogen and the subsequent viral dynamics is unclear in many disease settings, and this relationship has implications for both the timing and the required efficacy of antiviral therapy. Since experimental challenge studies often employ higher doses of virus than would generally be present in natural infection assessment of this relationship is particularly important for translation of findings. In this study we used mathematical modelling of viral load data from a multi-dose study of Zika virus infection in a macaque model to describe the impact of varying the dose of Zika virus on model parameters, and developed a novel mathematical model incorporating viral interference with the innate immune response.

Lassa viral dynamics in non-human primates treated with favipiravir or ribavirin

Lassa fever is an haemorrhagic fever caused by Lassa virus (LASV). There is no vaccine approved against LASV and the only recommended antiviral treatment relies on ribavirin, despite limited evidence of efficacy. Recently, the nucleotide analogue favipiravir showed a high antiviral efficacy, with 100% survival obtained in an otherwise fully lethal non-human primate (NHP) model of Lassa fever. However the mechanism of action of the drug is not known and the absence of pharmacokinetic data limits the translation of these results to the human setting. Here we aimed to better understand the antiviral effect of favipiravir by developping the first mathematical model recapitulating Lassa viral dynamics and treatment. We analyzed the viral dynamics in 24 NHPs left untreated or treated with ribavirin or favipiravir, and we put the results in perspective with those obtained with the same drugs in the context of Ebola infection. Our model estimates favipiravir EC50 in vivo to 2.89 μg.mL-1, which is much lower than what was found against Ebola virus. The main mechanism of action of favipiravir was to decrease virus infectivity, with an efficacy of 91% at the highest dose. Based on our knowledge acquired on the drug pharmacokinetics in humans, our model predicts that favipiravir doses larger than 1200 mg twice a day should have the capability to strongly reduce the production infectious virus and provide a milestone towards a future use in humans.

On the curved exponential family in the Stochastic Approximation Expectation Maximization Algorithm

The Expectation-Maximization Algorithm (EM) is a widely used method allowing to estimate the maximum likelihood of models involving latent variables. When the Expectation step cannot be computed easily, one can use stochastic versions of the EM such as the Stochastic Approximation EM. This algorithm, however, has the drawback to require the joint likelihood to belong to the curved exponential family. To overcome this problem, [16] introduced a rewriting of the model which “exponentializes” it by considering the parameter as an additional latent variable following a Normal distribution centered on the newly defined parameters and with fixed variance. The likelihood of this new exponentialized model now belongs to the curved exponential family. Although often used, there is no guarantee that the estimated mean is close to the maximum likelihood estimate of the initial model. In this paper, we quantify the error done in this estimation while considering the exponentialized model instead of the initial one. By verifying those results on an example, we see that a trade-off must be made between the speed of convergence and the tolerated error. Finally, we propose a new algorithm allowing a better estimation of the parameter in a reasonable computation time to reduce the bias.

Population physiologically based modeling of pirlimycin milk concentrations in dairy cows

Predictions of drug residues in milk are critical in food protection and are a major consideration in the economics of treatment of mastitis in dairy cows. Nonlinear mixed-effects modeling (NLME) has been advocated as a suitable pharmaco-statistical method for the study of drug residues in milk. Recent developments in physiologically based pharmacokinetic (PBPK) modeling of intramammary drugs allow the combination of a mechanistic description of milk pharmacokinetics with NLME methods. The PBPK model was applied to NLME analysis of a data set consisting of milk drug concentrations from 78 healthy cows and 117 with clinical mastitis. Pirlimycin milk pharmacokinetics were adequately described by the model across the range of observed concentrations. Mastitis was characterized by increased variance in milk production volume. Udder residual volume was larger in cows with 1, or 2 or greater diseased mammary glands than in the healthy cows. Low-producing cows had a greater risk of prolonged milk residues. With the exclusion of the low-production cows, the model predicted that healthy cows required a milk discard time 12 h longer than that indicated by the label, and the diseased cows 36 h longer than indicated by the label. More pirlimycin was systemically absorbed in the gram-positive infected compared with the gram-negative infected or healthy cows, suggesting a greater risk of violative meat residues in gram-positive infected cows. Using NLME and PBPK models, we identified factors associated with changes in pirlimycin milk residues that may affect food safety. This model extends the verification of a simple physiologically based framework for the study of intramammary drugs.

Comparison of empirical and dynamic models for HIV viral load rebound after treatment interruption

OBJECTIVE

To compare empirical and mechanistic modeling approaches for describing HIV-1 RNA viral load trajectories after antiretroviral treatment interruption and for identifying factors that predict features of viral rebound process.

METHODS

We apply and compare two modeling approaches in analysis of data from 346 participants in six AIDS Clinical Trial Group studies. From each separate analysis, we identify predictors for viral set points and delay in rebound. Our empirical model postulates a parametric functional form whose parameters represent different features of the viral rebound process, such as rate of rise and viral load set point. The viral dynamics model augments standard HIV dynamics models-a class of mathematical models based on differential equations describing biological mechanisms-by including reactivation of latently infected cells and adaptive immune response. We use Monolix, which makes use of a Stochastic Approximation of the Expectation-Maximization algorithm, to fit non-linear mixed effects models incorporating observations that were below the assay limit of quantification.

RESULTS

Among the 346 participants, the median age at treatment interruption was 42. Ninety-three percent of participants were male and sixty-five percent, white non-Hispanic. Both models provided a reasonable fit to the data and can accommodate atypical viral load trajectories. The median set points obtained from two approaches were similar: 4.44 log10 copies/mL from the empirical model and 4.59 log10 copies/mL from the viral dynamics model. Both models revealed that higher nadir CD4 cell counts and ART initiation during acute/recent phase were associated with lower viral set points and identified receiving a non-nucleoside reverse transcriptase inhibitor (NNRTI)-based pre-ATI regimen as a predictor for a delay in rebound.

CONCLUSION

Although based on different sets of assumptions, both models lead to similar conclusions regarding features of viral rebound process.

Concentration-response model of rilpivirine in a cohort of HIV-1-infected naive and pre-treated patients

BACKGROUND

Rilpivirine is widely prescribed in people living with HIV. Although trough plasma concentrations have been associated with virological response, the drug pharmacodynamics remain incompletely characterized.

OBJECTIVES

To develop the first pharmacodynamic model of rilpivirine in order to establish the rilpivirine concentration-response relationship for future treatment optimization.

METHODS

A retrospective observational study was conducted in patients receiving the once-daily rilpivirine/tenofovir disoproxil fumarate/emtricitabine regimen. Individual rilpivirine trough plasma concentrations over time were predicted using a previous pharmacokinetic model. An established susceptible, infected, recovered model was used to describe HIV dynamics without assuming disease steady-state. Population analysis was performed with MONOLIX 2018 software. Simulations of the viral load evolution as a function of time and rilpivirine trough plasma concentration were performed.

RESULTS

Overall, 60 naive and 39 pre-treated patients were included with a follow-up ranging from 2 to 37 months. The final model adequately described the data and the pharmacodynamic parameters were estimated with a good precision. The population typical value of rilpivirine EC50 was estimated at 65 ng/mL. A higher infection rate constant of CD4 cells for HIV-1 was obtained in pre-treated patients. Consequently, the time to obtain virological suppression was longer in pre-treated than in naive patients.

CONCLUSIONS

The concentration-response relationship of rilpivirine was satisfactorily described for the first time using an original population pharmacodynamic model. Simulations performed using the final model showed that the currently used 50 ng/mL rilpivirine trough plasma concentration efficacy target might need revision upwards, particularly in pre-treated patients.

A flexible nonlinear mixed effects model for HIV viral load rebound after treatment interruption

Characterization of HIV viral rebound after the discontinuation of antiretroviral therapy is central to HIV cure research. We propose a parametric nonlinear mixed effects model for the viral rebound trajectory, which often has a rapid rise to a peak value followed by a decrease to a viral load set point. We choose a flexible functional form that captures the shapes of viral rebound trajectories and can also provide biological insights regarding the rebound process. Each parameter can incorporate a random effect to allow for variation in parameters across individuals. Key features of viral rebound trajectories such as viral set points are represented by the parameters in the model, which facilitates assessment of intervention effects and identification of important pretreatment interruption predictors for these features. We employ a stochastic expectation-maximization (StEM) algorithm to incorporate HIV-1 RNA values that are below the lower limit of assay quantification. We evaluate the performance of our model in simulation studies and apply the proposed model to longitudinal HIV-1 viral load data from five AIDS Clinical Trials Group treatment interruption studies.

Parameter estimation and treatment optimization in a stochastic model for immunotherapy of cancer

Adoptive Cell Transfer therapy of cancer is currently in full development and mathematical modeling is playing a critical role in this area. We study a stochastic model developed by Baar et al. (2015) for modeling immunotherapy against melanoma skin cancer. First, we estimate the parameters of the deterministic limit of the model based on biological data of tumor growth in mice. A Nonlinear Mixed Effects Model is estimated by the Stochastic Approximation Expectation Maximization algorithm. With the estimated parameters, we return to the stochastic model and calculate the probability of complete T cells exhaustion. We show that for some relevant parameter values, an early relapse is due to stochastic fluctuations (complete T cells exhaustion) with a non negligible probability. Then, focusing on the relapse related to the T cell exhaustion, we propose to optimize the treatment plan (treatment doses and restimulation times) by minimizing the T cell exhaustion probability in the parameter estimation ranges.

Population Pharmacokinetic/Pharmacodynamic Modeling of Glenzocimab (ACT017) a Glycoprotein VI Inhibitor of Collagen-Induced Platelet Aggregation

Glenzocimab (ACT017) is a humanized monoclonal antigen‐binding fragment (Fab) directed against the human platelet glycoprotein VI, a key receptor for collagen and fibrin that plays a major role in thrombus growth and stability. Glenzocimab is being developed as an antiplatelet agent to treat the acute phase of ischemic stroke. During a phase I study in healthy volunteers, the population pharmacokinetics (PK) and pharmacodynamics (PD) of glenzocimab were modeled using Monolix software. The PK/PD model thus described glenzocimab plasma concentrations and its effects on ex vivo collagen‐induced platelet aggregation. Glenzocimab was found to have dose‐proportional, 2‐compartmental PK with a central distribution volume of 4.1 L, and first and second half‐lives of 0.84 and 9.6 hours. Interindividual variability in clearance in healthy volunteers was mainly explained by its dependence on body weight. The glenzocimab effect was described using an immediate effect model with a dose‐dependent half maximal inhibitory concentration: Larger doses resulted in a stronger effect at the same glenzocimab plasma concentration. The mechanism of the overproportional concentration effect at higher doses remained unexplained. PK/PD simulations predicted that 1000‐mg glenzocimab given as a 6‐hour infusion reduced platelet aggregation to 20% in 100% of subjects at 6 hours and in 60% of subjects at 12 hours after dosing. Simulations revealed a limited impact of creatinine clearance on exposure, suggesting that no dose adjustments were required with respect to renal function. Future studies in patients with ischemic stroke are now needed to establish the relationship between ex vivo platelet aggregation and the clinical effect.

Pharmacokinetics for haemophilia treaters: Meaning of PK parameters, interpretation pitfalls, and use in the clinic

Replacement therapy with concentrates of factor VIII or IX remains the gold standard for severe haemophilia management. The recent development of clotting factor products with extended half-life, widely available on the market since 2 years, facilitates adherence, improves considerably the patients' quality of life, and simplifies the management of breakthrough bleedings or surgery. These molecules have also brought to the limelight the concepts of optimization and personalization of anti-haemophilic prophylaxis. Pharmacokinetics (PK) is one of the tools that can help haematologists to adapt in a more objective and precise manner the prophylaxis regimen to each individual patient's specific needs. For many years, clinicians at haemophilia centres have been using some simple PK parameters, such as recovery and residual level. However, recently, they have been confronted with an important number of new PK parameters they were not familiar with, but that can be used to improve patient management. Due to the accumulation of PK data and their relative complexity, it is now necessary to analyse the relevance of the different PK parameters relative to haemophilia specificities, and also to know their limits to better use them.

Real Life Population Pharmacokinetics Modelling of Eight Factors VIII in Patients with Severe Haemophilia A: Is It Always Relevant to Switch to an Extended Half-Life?

We retrospectively analysed the data files of 171 adults and 87 children/adolescents with severe haemophilia, except for 14 patients (moderate; minor) (1), to develop a global population pharmacokinetic (PK) model for eight factors VIII (FVIII) that could estimate individual PK parameters for targeting the desired level of FVIII activity (FVIII:C); and (2) to compare half-life (HL) in patients switching from a standard half-life (SHL) to an extended half-life (EHL) and evaluate the relevance of the switch. One-stage clotting assay for the measurement of FVIII activity (FVIII:C, IU/mL) was used for population PK modelling. The software, Monolix version 2019R1, was used for non-linear mixed-effects modelling. A linear two-compartment model best described FVIII:C. The estimated PK parameters (between-subject variability) were: 2640 mL (23.2%) for volume of central compartment (V1), 339 mL (46.8%) for volume of peripheral compartment (V2), 135 mL/h for Q (fixed random effect), and 204 mL/h (34.9%) for clearance (Cl). Weight, age, and categorical covariate EHL were found to influence Cl and only weight for V1. This model can be used for all of the FVIII cited in the study. Moreover, we demonstrated, in accordance with previous studies, that Elocta had longer half-life (EHL) than SHL (mean ratio: 1.48) as compared to Advate, Factane, Kogenate, Novoeight, and Refacto.

Efficient Pharmacokinetic Modeling Workflow With the MonolixSuite: A Case Study of Remifentanil

MonolixSuite is a software widely used for model‐based drug development. It contains interconnected applications for data visualization, noncompartmental analysis, nonlinear mixed effect modeling, and clinical trial simulations. Its main assets are ease of use via an interactive graphical interface, computation speed, and efficient parameter estimation even for complex models. This tutorial presents a step‐by‐step pharmacokinetic (PK) modeling workflow using MonolixSuite, including how to visualize the data, set up a population PK model, estimate parameters, and diagnose and improve the model incrementally.

Impact of Antibiotic Gut Exposure on the Temporal Changes in Microbiome Diversity

Although the global deleterious impact of antibiotics on the intestinal microbiota is well known, temporal changes in microbial diversity during and after an antibiotic treatment are still poorly characterized. We used plasma and fecal samples collected frequently during treatment and up to one month after from 22 healthy volunteers assigned to a 5-day treatment by moxifloxacin (n = 14) or no intervention (n = 8). Moxifloxacin concentrations were measured in both plasma and feces, and bacterial diversity was determined in feces by 16S rRNA gene profiling and quantified using the Shannon index and number of operational taxonomic units (OTUs). Nonlinear mixed effect models were used to relate drug pharmacokinetics and bacterial diversity over time. Moxifloxacin reduced bacterial diversity in a concentration-dependent manner, with a median maximal loss of 27.5% of the Shannon index (minimum [min], 17.5; maximum [max], 27.7) and 47.4% of the number of OTUs (min, 30.4; max, 48.3). As a consequence of both the long fecal half-life of moxifloxacin and the susceptibility of the gut microbiota to moxifloxacin, bacterial diversity indices did not return to their pretreatment levels until days 16 and 21, respectively. Finally, the model characterized the effect of moxifloxacin on bacterial diversity biomarkers and provides a novel framework for analyzing antibiotic effects on the intestinal microbiome.

Dose tailoring of human cell line-derived recombinant factor VIII simoctocog alfa: Using a limited sampling strategy in patients with severe haemophilia A

AIMS

The use of factor VIII (FVIII) prophylaxis in haemophilia A is considered the standard of care, particularly in children. Despite adjustment of doses for body weight and/or age, a large pharmacokinetic (PK) variability between patients has been observed. PK-tailored prophylaxis may help clinicians adjust coagulation factor FVIII activity (FVIII:C) to the desired level, which may differ in individual patients. The objective was to develop a population PK model for simoctocog alfa based on pooled clinical trial data and to develop a Bayesian estimator to allow PK parameters in individual patients to be estimated using a reduced number of blood samples.

METHODS

PK data from 86 adults and 29 children/adolescents with severe haemophilia A were analysed. The FVIII data measured using 2 different assays (chromogenic and the 1-stage clotting assay) were fit to separate develop population PK models using nonlinear mixed-effect models. A Bayesian estimator was then developed to estimate the time above the threshold of 1%.

RESULTS

The PK data for chromogenic and the 1-stage clotting assays were both best described by a 2-compartment models. Simulations demonstrated good predictive capacity. The limited sampling strategy using blood sample at 3 and 24 hours allowed an accurate estimation of the time above the threshold of 1% FVIII:C (mean bias 0.01 and 0.11, mean precision 0.18 and 0.45 for 2 assay methods).

CONCLUSION

In this study, we demonstrated that a Bayesian approach can help to reduce the number of samples required to estimate the time above the threshold of 1% FVIII:C with good accuracy.

Flexible longitudinal linear mixed models for multiple censored responses data

In biomedical studies and clinical trials, repeated measures are often subject to some upper and/or lower limits of detection. Hence, the responses are either left or right censored. A complication arises when more than one series of responses is repeatedly collected on each subject at irregular intervals over a period of time and the data exhibit tails heavier than the normal distribution. The multivariate censored linear mixed effect (MLMEC) model is a frequently used tool for a joint analysis of more than one series of longitudinal data. In this context, we develop a robust generalization of the MLMEC based on the scale mixtures of normal distributions. To take into account the autocorrelation existing among irregularly observed measures, a damped exponential correlation structure is considered. For this complex longitudinal structure, we propose an exact estimation procedure to obtain the maximum-likelihood estimates of the fixed effects and variance components using a stochastic approximation of the EM algorithm. This approach allows us to estimate the parameters of interest easily and quickly as well as to obtain the standard errors of the fixed effects, the predictions of unobservable values of the responses, and the log-likelihood function as a byproduct. The proposed method is applied to analyze a set of AIDS data and is examined via a simulation study.

Population pharmacokinetics of serelaxin in patients with acute or chronic heart failure, hepatic or renal impairment, or portal hypertension and in healthy subjects

AIMS

Serelaxin is a recombinant human relaxin-2 peptide being developed for the treatment of acute heart failure (AHF). The present analyses aimed to evaluate serelaxin pharmacokinetics following intravenous administration and to identify covariates that may explain pharmacokinetic variability in healthy subjects and patients.

METHODS

Serum concentration-time data for 613 subjects from nine phase I and II studies were analysed using a nonlinear mixed-effects model to estimate population pharmacokinetics and identify significant covariates. A quantile regression analysis was also conducted to assess the relationship between clearance and covariates by including sparse data from a phase III study.

RESULTS

A three-compartment disposition model was established to describe serelaxin pharmacokinetics. Three out of 23 covariates, including baseline body mass index (BMI) and estimated glomerular filtration rate (eGFR) and study A1201, were identified as significant covariates for clearance but with a moderate impact on steady-state concentration, reducing the intersubject variability from 44% in the base model to 41% in the final model with covariates. The steady-state volume of distribution (Vss) was higher in patients with AHF (544 ml kg^-1 ) or chronic heart failure (434 ml kg^-1 ), compared with typical nonheart failure subjects (347 ml kg^-1 ). Quantile regression analysis showed that a 20% increase in BMI or a 20% decrease in eGFR decreased serelaxin clearance by 9.2% or 5.2%, respectively.

CONCLUSIONS

Patients with HF showed higher Vss but similar clearance (and therefore steady-state exposure) vs. non nonheart failure subjects. BMI and eGFR were identified as the main covariates explaining intersubject variability in clearance; however, the impact of these covariates on steady-state concentration was moderate and therefore unlikely to be clinically relevant.