Assessing the utility of an anti-malarial pharmacokinetic-pharmacodynamic model for aiding drug clinical development

Abstract

Background

Mechanistic within-host models relating blood anti-malarial drug concentrations with the parasite-time profile help in assessing dosing schedules and partner drugs for new anti-malarial treatments. A comprehensive simulation study to assess the utility of a stage-specific pharmacokinetic-pharmacodynamic (PK-PD) model for predicting within-host parasite response was performed.

Methods

Three anti-malarial combination therapies were selected: artesunate-mefloquine, dihydroartemisinin-piperaquine, and artemether-lumefantrine. The PK-PD model included parameters to represent the concentration-time profiles of both drugs, the initial parasite burden and distribution across the parasite life cycle, and the parasite multiplication factor due to asexual reproduction. The model also included the maximal killing rate of each drug, and the blood drug concentration associated with half of that killing effect (in vivo EC50), derived from the in vitro IC50, the extent of binding to 0.5% Albumax present in the in vitro testing media, and the drugs plasma protein binding and whole blood to plasma partitioning ratio. All stochastic simulations were performed using a Latin-Hypercube-Sampling approach.

Results

The simulations demonstrated that the proportion of patients cured was highly sensitive to the in vivo EC50 and the maximal killing rate of the partner drug co-administered with the artemisinin derivative. The in vivo EC50 values that corresponded to on average 95% of patients cured were much higher than the adjusted values derived from the in vitro IC50. The proportion clinically cured was not strongly influenced by changes in the parameters defining the age distribution of the initial parasite burden (mean age of 4 to 16 hours) and the parasite multiplication factor every life cycle (ranging from 8 to 12 fold/cycle). The median parasite clearance times, however, lengthened as the standard deviation of the initial parasite burden increased (i.e. the infection became more asynchronous).

Conclusions

This simulation study demonstrates that the PD effect predicted from in vitro growth inhibition assays does not accord well with the PD effect of the anti-malarials observed within the patient. This simulation-based PK-PD modelling approach should not be considered as a replacement to conducting clinical trials but instead as a decision tool to improve the design of a clinical trial during drug development.

Background

Despite significant progress in control of malaria over the last decade, it remains a major global health problem. Almost 40% of the world’s population live in malaria endemic areas, with each year about a quarter of a billion people experiencing clinical malaria and an estimated 655,000 malaria-related deaths[1]. With no vaccine currently available, malaria control relies on preventative measures (i.e. insecticidal bed nets and indoor residual spraying) and effective treatment with artemisinin-based combination therapy (ACT). ACT involves treatment with two or more anti-malarials; a fast acting but short lived artemisinin derivative and a less effective, but of longer duration, partner drug. ACT is recommended by WHO as the first-line treatment of uncomplicated falciparum malaria[2], but recent reports from western Cambodia raise concerns that Plasmodium falciparum has developed reduced susceptibility to oral artesunate[3, 4]. In the context of emerging resistance to artesunate (the most widely used artemisinin derivative), it is critical that new anti-malarial treatments are developed and assessed.

Mechanistic within-host pharmacokinetic-pharmacodynamic (PK-PD) models that relate blood anti-malarial drug concentrations to the parasite-time profile have potential to aid anti-malarial drug development. Simulated parasite-time profiles for hypothetical patients can be generated from the mechanistic PK-PD model and incorporate between-patient variability in the drug concentration profiles. Comparisons of parasitological outcomes (e.g. distribution of parasite clearance times and proportion of patients cured) derived from these hypothetical individuals can then be used as a decision tool for assessing dosing schemes and potential partner drugs for new anti-malarial drugs. This simulation-based approach has been adopted previously using a within-host continuous-time PK-PD model for comparing dosing schemes of mefloquine[5, 6], artesunate[7], chloroquine[8], and the ACT, mefloquine and artesunate[9]. More recently a parasite stage-specific discrete-time within-host PK-PD model has been developed, and a stochastic simulation-based approach implemented to compare dosing schemes for artesunate[10].

In this paper, the above stage-specific model was extended to account for the action of two or more anti-malarial treatments, and the anti-malarial pharmacodynamic parameters were determined by extrapolating from in vitro data. Using a Latin-Hypercube-Sampling approach[11], the sensitivity of the PK-PD model to particular parameters was assessed, by comparing, across different sets of parameter values, the parasitological outcomes derived from simulated parasite-time profiles of hypothetical patients.

Methods

Within-host pharmacokinetic-pharmacodynamic model

The within-host PK-PD model is based on that described in Saralamba et al.[10], which is a discrete-time model that incorporates the age distribution of the parasite population within the malaria patient. This model determines how the distribution of the parasite age changes post treatment as a consequence of the concentration of the anti-malarial drug. The general form of the discrete-time model (see Saralamba et al.[10] supplemental information for more detail) and PK-PD parameter values used to simulate individual parasitaemia-time profiles in the presence of anti-malarial combination therapies is described below.

Prior to drug administration, the initial parasite load of each patient (P_0i) is distributed among the 48 hourly age intervals of the P. falciparum life cycle according to a Gaussian distribution with a mean age of μ hours and a standard deviation of σ hours (see Additional file1).

The expected number of parasites in patient i ($N_i\left(a,t\right)$) aged ‘a’ hours (where 1 ≤ a ≤ 48) at hourly time point t after drug administration, is expressed by the following difference equations,

$$N_i\left(a,t\right)=N_i\left(a-1,t-1\right)e^{-k_i\left(a-1,t-1\right)}\text{for\hspace{0.17em}1}<a\le 48$$

(1)

and

$$N_i\left(1,t\right)=PMF\times N_i\left(48,t-1\right)e^{-k_i\left(48,t-1\right)}\text{for}a=1\text{.}$$

(2)

In (1) and (2),$k_i\left(a,t\right)$ is the kill rate constant for the combination therapy of the parasites aged a hours at hourly time point t, and in (2), PMF is the parasite multiplication factor that represents the number of merozoites released per schizont at the end of the 48 hour life cycle that successfully reinvade red blood cells.

The kill rate constant for a particular combination therapy was calculated as follows:

$$k_i\left(a,t\right)=k_{dru{g}_1,i}\left(a,t\right)+k_{dru{g}_2,i}\left(a,t\right)$$

(3)

where$dru{g}_{\text{j}}$, for j = 1, 2, denotes one of the drugs comprising a combination therapy. The kill rate constant in (3) assumes that the effects of$dru{g}_1$ and$dru{g}_2$ on the parasite population in vivo are independent of one another.

The relationship between the kill rate constant for each drug and the drug concentration is given by:

$$k_i\left(a,t\right)=k_{\text{max}}\left(a\right)\times \left[c_i{\left(t\right)}^{\gamma }/\left(c_i{\left(t\right)}^{\gamma }+E{C}_{50}^{\gamma }\right)\right]$$

(4)

where$c_i\left(t\right)$ is the plasma drug concentration at time t for patient i, γ is the slope of the concentration-effect relationship,$E{C}_{50}$ is the blood drug concentration in vivo that gives 50% parasite killing and k_max is the maximal killing rate of each drug (which was assumed to be constant across the parasite ages ‘a’ where the drug is known to have an effect).

The drug concentration-time curve$c_i\left(t\right)$ for each drug assumes the form of a structural PK model, for example a one- or two- compartment PK model with first-order absorption and elimination from the central body compartment.

Simulation study

The model was implemented using R[12]. Parasite counts at several different time points and in the presence of anti-malarial combination therapies were simulated from this discrete-time model for hypothetical malaria patients. The summary measures that were derived from each simulated parasite count-time curve were: (i) the hypothetical patient was clinically cured (defined as the parasite count falling below 2.5 × 10⁸ parasites (i.e. 50 parasites per μL) and not reappearing above 2.5 × 10⁸ parasites by day 63 of follow-up); and (ii) the parasite clearance time (PCT) (hours), defined as the time for the circulating parasite count (parasites aged approximately 1 to 26 hours) to decrease below 2.5 × 10⁸ parasites. Circulating parasite counts were calculated from the total parasite counts following Saralamba et al.[10], where sequestration was estimated to start at a parasite age of 11 hours and the number of parasites older than 11 hours circulating in the blood was assumed to decrease exponentially.

The anti-malarial combination therapies selected for this simulation study were three artemisinin-based combination therapies: artesunate-mefloquine, dihydroartemisinin-piperaquine, and artemether-lumefantrine. Studies of in vitro interactions between the pharmacodynamic effects of the drugs have shown no interaction between dihydroartemisinin and piperaquine[13], and a small amount of synergy between artesunate and mefloquine[14], and between artemether and lumefantrine[15]. Thus, the assumption of independent pharmacodynamic effects of each drug in the combination therapies selected for this simulation study seems reasonable, especially considering the short amount of time (approximately six hours for dihydroartemisinin and twelve hours for artemether) that the drug concentrations of the artemisinin derivatives are present in the blood.

Pharmacokinetic-pharmacodynamic parameters

In order to assess the utility of the PK-PD model for determining patient outcomes (described above) the sensitivity of the model to parameter values was explored. This was implemented using Latin-Hypercube-Sampling (LHS). LHS is a method which is used to randomly sample over large parameter spaces in an evenly distributed manner[11].

Before carrying out the LHS sampling, pharmacokinetic profiles of each anti-malarial drug for the three artemisinin-based combination therapies were simulated for 100 hypothetical patients. The dosing regimen used was the regimen recommended by the WHO ([2]; see Table 1) and the PK profiles were simulated using parameter values and between-subject variability obtained from the literature (see Table 2 and Additional file2)[16–19].

Table 1 Dosing regimen for each artemisinin-based combination therapy (ACT)

Table 2 Population pharmacokinetic parameter values (BSV%) ^† for each drug

The PD parameters were varied across the LHS experiment to capture both biological and empirical uncertainty. The definition of each PD parameter and the statistical distribution selected for each PD parameter are given in Table 3 and Table 4, respectively. For the parameters that are drug independent, the distributions for the parameters which determine how the parasites are distributed across the 48 hours of the parasite life cycle (μ and σ) were sourced from PK-PD modelling of uncomplicated falciparum malaria patients[10] and the distribution for the parasite multiplication factor (PMF) was obtained from modelling of data collected from syphilis patients treated with an induced malaria infection[20, 21].

Table 3 Parameter definitions for the within-host pharmacokinetic-pharmacodynamic model

Table 4 Statistical distributions selected for each pharmacodynamic parameter

The maximal killing rate (k_max) of the drug was assumed to follow a triangular distribution. The middle values were taken from published clinical data[3, 22, 23]. Piperaquine was the only drug with no clinical studies of it administered as a monotherapy and the mode was set to a value equal to that derived from in vitro experiments[24]. The minimum and maximum values of the triangular distribution for each anti-malarial correspond to a 50-fold decrease and 50-fold increase in the number of parasites killed every 48 hours. Artesunate, dihydroartemisnin and artemether were assumed to kill parasites aged 6 to 44 hours; mefloquine 18 to 40 hours; piperaquine (assumed similar to chloroquine) 12 to 36 hours; and lumefantrine (assumed similar to mefloquine) 18 to 40 hours[25].

The slope of the concentration-effect curve (i.e. in vivo γ) was assumed to have the same statistical distribution (i.e. Log-normal) and parameter values (mean and standard deviation on log_e scale) as the in vitro γ which was derived from modelling of in vitro concentration-effect curves measured from a large number of parasite isolates.

The distribution of the EC50 for each anti-malarial drug is unknown, therefore, the conservative continuous-uniform distribution was chosen with the minimum value set to the adjusted in vitro IC50 and the maximum value equal to half of the maximum concentration of the population mean PK profile.

Five thousand parameter sets were selected from the above statistical distributions using LHS. For each parameter set, simulated parasite count versus time profiles for the 100 hypothetical patients (with the PK profiles determined above) were derived for each artemisinin-based combination therapy. The initial parasite burden for each of the 100 hypothetical patients was randomly selected from a log-normal distribution with a geometric mean of 1.14 × 10¹¹ (i.e. parasitaemia of 22746 parasites per μL) and a standard deviation on the log-scale of 1.13. The initial parasite burdens for the 100 patients did not vary with LHS parameter set or with the artemisinin combination therapy used in the simulation.

Determination of the in vitro IC50 and slope (γ) of concentration-effect relationship

Estimates of the in vitro IC50 (not corrected for binding) and γ for artesunate, dihydroartemisinin, mefloquine, piperaquine, and lumefantrine (refer to Table 5), were determined from statistical modelling of individual isolate effect versus drug concentration curves. The fresh P. falciparum parasite isolates were obtained from blood samples of 487 patients attending outpatient clinics in Papua, Indonesia between 2004 and 2010[26]. The in vitro drug susceptibility was determined using the World Health Organization guidelines for schizont maturation tests. The in vitro data for artemether were measured at the Swiss Tropical and Public Health Institute (Switzerland; Basel) against asynchronous intraerythrocytic forms of the P. falciparum strain NF54 (obtained from MR4) using the [3H]-hypoxanthine incorporation assay[27].

Table 5 Parameter values derived from in vitro experiments

The in vitro free drug IC50 was calculated from the measured in vitro IC50 (uncorrected for binding) by multiplying by the unbound fraction in the in vitro testing media. This value was then converted to an adjusted in vitro IC50 in whole blood (to represent concentrations comparable to an in vivo situation) by first dividing by the free fraction in plasma and then multiplying by the whole blood to plasma ratio (see Table 5).

Estimates of drug binding to the in vitro testing media (i.e. 0.5% Albumax in RPMI) and human plasma were determined using ultracentrifugation. Briefly, blank Albumax media and plasma were each spiked with compound and divided into six aliquots; three aliquots were subjected to ultracentrifugation (Beckman Optima XL-100K Ultracentrifuge, Rotor type 42.2 Ti; 223,000 x g) for 4.2 hours at 37°C to pellet the proteins whereas the remaining three aliquots served as controls and were incubated at 37°C for the same time period but without centrifugation. Controls were also stored at −20°C to confirm sample stability over the centrifugation period. Aliquots of the supernatants from the ultracentrifuged samples were first diluted 1:1 in acetonitrile, assayed by LC-MS and the responses compared to a calibration curve prepared in 50% aqueous acetonitrile to determine the unbound (i.e. free) concentration. Control concentrations in each matrix were determined using LC-MS by first precipitating the proteins with acetonitrile (3:1 acetonitrile: matrix) and then comparing the responses to a calibration curve from a blank matrix prepared using the same protein precipitation procedure. The free fraction in each matrix was then determined from the ratio of the average unbound (e.g. free) concentration to the average total concentration in each matrix.

The whole blood to plasma partitioning ratio (B/P) was determined by spiking aliquots of whole blood or plasma maintained at 37°C with compound, incubating for 2 min, and then centrifuging the whole blood sample to obtain the plasma fraction. Both the plasma fraction of whole blood and the plasma control were assayed for compound by LC-MS as described above. The blood to plasma ratio was calculated from the ratio of the concentration in the plasma control (used as a surrogate for the total whole blood concentration since whole blood assays were not available for each compound) to that in the plasma fraction of whole blood. The two-minute time point was used to avoid confounding issues due to potential blood instability; rapid equilibration between plasma and erythrocytes was assumed.

Results

Sensitivity of pharmacokinetic-pharmacodynamic (PK-PD) model to parameter values

Examined first was the model’s sensitivity to the three parameters that describe the function of the parasite: the parasite multiplication factor (PMF), the mean for the age distribution of the initial parasite burden (μ), and the standard deviation of the age distribution of parasites within each host (σ). For each drug combination therapy and model output (proportion clinically cured and PCT), Figure 1 shows tornado plots of the partial rank correlation coefficients (PRCCs). The magnitude of the PRCC indicates the importance of the uncertainty in estimating the parameters governing the age distribution of the initial parasite burden (μ and σ) and the PMF in contributing to the variability in the proportion clinically cured and the PCTs. A negative PRCC indicates that the model output tends to decrease as the parameter increases and a positive PRCC indicates that the model output tends to increase as the parameter increases.

Figure 1

Tornado plots of partial rank correlation coefficients, indicating the importance of each drug independent parameter’s (mean and standard deviation of the age distribution of the initial parasite burden, i.e. μ and σ , and the PMF) uncertainty in contributing to the variability in the proportion cured (left) and parasite clearance time (PCT) (right) for each artemisinin combination therapy.

Figure 1 (left panel) shows that the proportion clinically cured (for ease of exposition referred to subsequently as proportion cured) after treatment with artemether-lumefantrine was sensitive to changes in the standard deviation of the age distribution of the initial parasite burden, and that it tended to decrease as the standard deviation of the age distribution (i.e.σ) increases (i.e. as the infection becomes more asynchronous). The proportion cured after treatment with artemisinin-lumefantrine was not strongly influenced by changes in the mean of the age distribution of the initial parasite burden (i.e. μ) and the PMF. The proportion cured for the remaining artemisinin combination therapies was not very sensitive to changes in the parameters defining the age distribution of the initial parasite burden (i.e. μ and σ) and the PMF.

The PCTs for all three artemisinin combination therapies were sensitive to changes in the mean and standard deviation of the age distribution of the initial parasite burden (see Figure 1, right panel). The PRCCs plotted in the right panel of Figure 1 show that the PCTs for all three artemisinin combination therapies tended to lengthen as the standard deviation of the age distribution increased (or the infection became more asynchronous) and tended to shorten as the mean of the age distribution increased.

Tornado plots of the PRCCs between the model outputs and the drug dependent parameters (EC50, k_max and γ) for each artemisinin combination therapy are also provided in Additional file3.

Artesunate-mefloquine

Figure 2 depicts the distribution of the proportion cured for the combination therapy artesunate-mefloquine across deciles of the 5,000 parameter values for EC50 (drug concentration in vivo that corresponds to 50% parasite killing) and k_max (maximum killing rate constant) for both artesunate and mefloquine. The proportion of patients cured from the simulated parasite-time profiles was highly correlated with the in vivo EC50 and k_max for mefloquine, whereas the influence of the artesunate parameter values was marginal. Lower killing rates of mefloquine predicted a reduction in the percentage of patients being cured, the median value ranging from 10% to 95% for k_max<0.187 (i.e. Parasite Reduction Ratio (PRR) <10^1.33) up to k_max>0.301 (i.e. PRR>10^2.42). Only when the EC50 concentrations of mefloquine were between 447 to 554 ng/ml, much higher than the adjusted in vitro IC50 value (21 ng/ml), was the predicted median proportion cured similar to that observed in clinical studies[3, 28, 29]. With EC50 values between 661–874 ng/ml approximately 75% (on average) of the hypothetical patients were predicted to be cured. Figure 3 illustrates that PCTs are not vulnerable to the EC50 and k_max values of mefloquine and artesunate. Furthermore the values of the slope of either the artesunate or mefloquine concentration-effect curves had no association with any of the model outputs (see Additional files4A-B).

Figure 2

Distribution of proportion cured within the EC50 and k _max deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, artesunate (ARS) and mefloquine (MQ). Top panels are for artesunate (EC50 – left hand side, k_max right hand side) and bottom panels are for mefloquine (EC50 – left hand side, k_max right hand side). Each individual box (with whiskers) represents the distribution of the proportion cured simulated for 5000 parameter values of EC50 (or k_max) within the range of cut-off values of that decile grouping, all other parameters were varied.

Figure 3

Distribution of parasite clearance times (hours) within the EC50 and k _max deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, artesunate (ARS) and mefloquine (MQ). Top panels are for artesunate (EC50 – left hand side, k_max right hand side) and bottom panels are for mefloquine (EC50 – left hand side, k_max right hand side).

Dihydroartemisinin-piperaquine

Similarly for the combination of dihydroartemisinin-piperaquine, the proportion of patients cured correlated closely with the in vivo EC50 for the long acting partner drug, piperaquine (Figure 4), but was not influenced by values of EC50 or k_max of dihydroartemisinin. However in contrast to the mefloquine containing combination, values of piperaquine k_max were only weakly associated with the proportion cured. The median proportion cured changed dramatically from 100% to 40% across the deciles of EC50 values for piperaquine, with EC50 values of 45 to 53 ng/ml corresponding to a median of 95%. Figure 5 highlights that PCTs are not sensitive to k_max values of dihydroartemisinin and piperaquine but are marginally sensitive to the EC50 values of piperaquine. The values of the slope of either the dihydroartemisinin or piperaquine concentration-effect curves had no associations with the model outputs (see Additional files5A-B).

Figure 4

Distribution of proportion cured within the EC50 and k ₁ deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, dihydroartemisinin (DHA) and piperaquine (PQ). Top panels are for dihydroartemisinin (EC50 – left hand side, k_max right hand side) and bottom panels are for piperaquine (EC50 – left hand side, k_max right hand side).

Figure 5

Distribution of parasite clearance times (hours) within the EC50 and k ₁ deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, dihydroartemisinin (DHA) and piperaquine (PQ). Top panels are for dihydroartemisinin (EC50 – left hand side, k_max right hand side) and bottom panels are for piperaquine (EC50 – left hand side, k_max right hand side).

Artemether-lumefantrine

For artemether-lumefantrine, the proportion of patients cured was correlated with both parameters k_max and EC50 of lumefantrine and marginally sensitive to both parameter values of artemether (Figure 6). For lumefantrine, in vivo EC50 values of less than 235 ng/ml corresponded to 100% cured and values of 1632 ng/ml and above resulted in approximately 10-20% cured. For this combination therapy, artemether is given at 0, 8, 24, 36, 48 and 60 hours, whereas for the above two combination therapies the artemisinin derivative is given only at 0, 24 and 48 hours. The additional doses of artemether may explain why changes in the EC50 and k_max values of artemether are associated with the proportion cured for this artemisinin derivative. This was supported by the results of a simulation run on the first 500 of the 5,000 LHS parameter sets with artemether given at 0, 24 and 48 hours and lumefantrine given at 0, 8, 24, 36, 48 and 60 hours (see Additional file6). For PCTs, as with the other partner drugs, no associations were observed for the EC50 and k_max of lumefantrine (Figure 7). A slight gradient was observed such that longer PCTs were observed for those with higher artemether EC50 values. The values of the slope of either the artemether or lumefantrine concentration-effect curves had no associations with the model outputs (see Additional files7A-B).

Figure 6

Distribution of proportion cured within the EC50 and k ₁ deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, artemether (ART) and lumefantrine (LM). Top panels are for artemether (EC50 – left hand side, k_max right hand side) and bottom panels are for lumefantrine (EC50 – left hand side, k_max right hand side).

Figure 7

Distribution of parasite clearance times (hours) within the EC50 and k ₁ deciles derived from the 5000 parameter sets for the anti-malarial combination therapy, artemether (ART) and lumefantrine (LM). Top panels are for artemether (EC50 – left hand side, k_max right hand side) and bottom panels are for lumefantrine (EC50 – left hand side, k_max right hand side).

Figure 8

A-B Proportion of 100 hypothetical patients cured. Panel A: EC50 values of mefloquine (287, 554 and 821 ng/ml) for two different dosing schemes: 8.3 mg/kg at 0, 24 and 48 hours (WHO standard; circles); and 15 mg/kg at 48 hours (alternative scheme; crosses). Panel B: EC50 values of piperaquine (32, 53 and 74 ng/ml) for two different dosing schemes: 18 mg/kg at 0, 24 and 48 hours (WHO standard; circles); and 36 mg/kg at 0, 24 and 48 hours (alternative scheme; crosses).

Comparing alternative dosing regimens

The utility of the mechanistic within host PK-PD model as a decision tool for comparing dosing regimens was investigated for different dosing schemes for each of the partner drugs, mefloquine (WHO recommended dosing of 8.3 mg/kg at 0, 24 and 48 hours versus dosing of 15 mg/kg at 48 hours) and piperaquine (WHO recommended dosing of 18 mg/kg at 0, 24 and 48 hours versus 36 mg/kg at 0, 24 and 48 hours), administered with the artemisinin derivatives, artesunate and dihydroartemisinin respectively. Since the proportion of patients cured was highly dependent on the EC50 value of the partner drug (a parameter for which we don’t know the value), different dosing schemes were compared for the 25^th, 50^th and 75^th percentile values of the 5000 parameter values selected using LHS. All other parameters were fixed at their mean/mode value given in Table 4.

The proportion cured for 100 hypothetical patients with varying PK profiles for the two dosing schemes of mefloquine was 100% when the EC50 value of mefloquine was 287 ng/ml, however, under-dosing patients with only 15 mg/kg of mefloquine resulted in only 50% of the patients being cured if the EC50 is 554 ng/ml and 20% if the EC50 was 821 ng/ml. For the standard recommended mefloquine dose of 8.3 mg/kg at 0, 24 and 48 hours, approximately 95% and 80% would be expected to be cured at the respective EC50 values of 554 and 821 ng/ml (Figure8A).

Doubling the recommended WHO dose of piperaquine given each day over 3 days results in a higher proportion of patients cured at the EC50 values of 53 and 74 ng/ml for piperaquine (Figure8B).

Discussion

The parasitological outcomes simulated in this paper were proportion of patients cured and parasite clearance times, for three different artemisinin-based combination therapies currently recommended by the WHO as the first line treatment for uncomplicated falciparum malaria. This simulation study was comprehensive, randomly drawing from each distribution of the six key pharmacodynamic parameters using Latin-Hypercube-Sampling (LHS), and included between-patient variability in the pharmacokinetic profiles of each anti-malarial drug. The proportion of hypothetical patients cured was observed to be highly correlated to the in vivo EC50 and the killing rate (k_max) of the partner drug co-administered with the artemisinin derivative. However, in vivo EC50 values that corresponded to on average 95% of patients cured (a value observed in most clinical efficacy studies of these regimens) were much higher than the values we derived from in vitro data (i.e. adjusted in vitro IC50), even though the difference in protein binding in vitro and in vivo were taken into account in the model. In vitro experiments typically assess the pharmacodynamic effect of an anti-malarial drug by measuring inhibition of parasite growth in rising concentration of free drug for the length of one parasite life cycle (i.e. 48 hours). The duration of the in vitro assay is usually 48 to 72 hours and although this permits a reproducible estimate of parasite drug susceptibility, it may be too short for this estimate to reflect accurately clinical correlates. This may explain in part why in vitro measures accord so poorly with the observed pharmacodynamic effect in vivo.

There was evidence that the proportion cured decreased if asynchronous infections were treated with artemether-lumefantrine, but the synchronicity of the infection did not strongly influence the proportion cured following treatment with either artesunate-mefloquine or dihydroartemisinin-piperaquine. The PCT for all three artemisinin combination therapies tended to lengthen as the infection became more asynchronous and to shorten as the mean age of the initial parasite burden increased. The finding that asynchronous infections take longer to clear is plausible because it is more likely that there will be parasites outside the killing zone (i.e. early rings and schizonts) during the times when the patient is only exposed to the partner drug (e.g. approximately 7 to 24 hours for mefloquine or piperaquine). The distributions of PCTs in Additional files4B, Additional file5B and Additional file6B are bimodal because the hypothetical patients do not receive a single dose of the artemisinin derivative but are given multiple doses at 24 and 48 hrs (and also 8, 36 and 60 hours for artemether).

The proportion cured following treatment with artesunate-mefloquine was highly correlated to the in vivo EC50 value of mefloquine. This is not surprising given the brief time the parasite is exposed to dihydroartemisinin concentrations (approximately 6 hours following each dose of artesunate) compared to mefloquine which remains in the body, on average, for 40 days. This finding concurs with observations from deterministic simulated individual patient parasite versus time profiles using a continuous-time PK-PD model[9]. The association between the proportion cured and the in vivo EC50 values of piperaquine and lumefantrine was even stronger than that observed for mefloquine. Both piperaquine and lumefantrine have an enormous volume of distribution and an elimination profile that comprises a steep short distribution phase followed by a slow elimination phase from day 5–7 onwards[19, 30]. However, lumefantrine has a much lower volume of distribution compared with piperaquine and this explains why higher values of in vivo EC50 for lumefantrine are required before the simulated observations predict, on average, 10-20% cured. PCTs did not increase when the maximal killing rate (k_max) of artesunate/dihydroartemisinin decreased, although this was observed for the ring stage parasites in the discrete-time PK-PD model reported by Saralamba et al.[10]. These conflicting findings are likely to arise for a number of reasons. First, in this paper it was assumed that the maximal killing rate of the artemisinin derivatives was constant across the killing zone (i.e. the age range of parasites for which the drug kills) since stage-specific killing rates for each anti-malarial was not known. Second, k_max values for the artemisinin derivatives were randomly selected from a Parasite Reduction Ratio at 48 hours (PRR₄₈) ranging from 5 × 10^4.28 to 5 × 10^6.28 with a mode value of 10^5.28[3] parasites reduced every 48 hours. In the observations by Saralamba et al.[10], the maximal killing rate of artesunate observed for the ring stages with delayed PCTs was a mean of 62% /cycle corresponding to a much lower PRR (~10^0.42) than our minimal value. Third, in this study the partner drugs were administered at the same times as the doses of the artemisinin derivatives and therefore contributed to the parasite clearance times whereas in Saralamba et al. only artesunate was administered in the first 48 hours of treatment.

This simulation study has a number of strengths which includes: the method of LHS for randomly selecting 5000 sets of the pharmacodynamic parameter values combined with simulations of 100 pharmacokinetic profiles for each anti-malarial to capture between-patient variability in drug exposure. Moreover our comparison of alternative dosing regimens highlight the utility of PK-PD models to compare dosing schemes and have the capacity to examine the association between a range of PK parameters (e.g. time above therapeutic concentration, maximum concentration or area under the concentration-time profile) and parasitological outcome. The limitations of this study were: the within-host PK-PD model assumes that the background immunity of the hypothetical patients was low or absent; no pharmacodynamic synergism between the two anti-malarials of each combination therapy evaluated was assumed; and the a priori assumption that the maximal killing rate of the artemisinin derivatives and the partner drugs remained constant across the different ages of the parasite within the defined killing zone (e.g. 6 to 44 hours for artesunate). Furthermore, the median PCTs were approximately 24 hours whereas in many clinical studies approximately 48 hours is often observed. This difference may be due to one or a combination of factors including: an assumption that there was no synergy between the drugs; the maximal killing rate of the partner drug, piperaquine, was taken from an in vitro experiment and thus may be higher than observed in vivo; and the number of circulating parasites was calculated in time steps of one hour post initial treatment for determining PCT whereas in clinical efficacy studies this is often determined from blood smears collected only every 24 hours.

In conclusion, this simulation study demonstrates the utility of using mechanistic within-patient PK-PD models for comparing parasitological outcomes of different dosing schemes of anti-malarial treatments and different anti-malarial combination therapies. The findings of this study suggest that the parasitological outcomes be compared for a number of scenarios of the pharmacodynamic parameter values, especially the unknown in vivo EC50 value. These simulation studies should not be used as a replacement to conducting the clinical efficacy trials but instead used to assist in determining the best dosing schemes and potential partner drugs to be considered for new anti-malarial treatments. This simulation-based approach has the potential to reduce the number of clinical efficacy trials carried out in the Phase II and Phase III stages of drug development, which will reduce the cost of drug development, speed up the process of drug registration, and could help identify non-ethical trials of malaria patients.