Exposure–response modelling approaches for determining optimal dosing rules in children

Within paediatric populations, there may be distinct age groups characterised by different exposure–response relationships. Several regulatory guidance documents have suggested general age groupings. However, it is not clear whether these categorisations will be suitable for all new medicines and in all disease areas. We consider two model-based approaches to quantify how exposure–response model parameters vary over a continuum of ages: Bayesian penalised B-splines and model-based recursive partitioning. We propose an approach for deriving an optimal dosing rule given an estimate of how exposure–response model parameters vary with age. Methods are initially developed for a linear exposure–response model. We perform a simulation study to systematically evaluate how well the various approaches estimate linear exposure–response model parameters and the accuracy of recommended dosing rules. Simulation scenarios are motivated by an application to epilepsy drug development. Results suggest that both bootstrapped model-based recursive partitioning and Bayesian penalised B-splines can estimate underlying changes in linear exposure–response model parameters as well as (and in many scenarios, better than) a comparator linear model adjusting for a categorical age covariate with levels following International Conference on Harmonisation E11 groupings. Furthermore, the Bayesian penalised B-splines approach consistently estimates the intercept and slope more accurately than the bootstrapped model-based recursive partitioning. Finally, approaches are extended to estimate Emax exposure–response models and are illustrated with an example motivated by an in vitro study of cyclosporine.


Supplementary Appendix A: Worked example of methods
This supplementary appendix aims to give an illustration of the output that would be seen from fitting each of the methods in Section 4 to a single set of simulated data. For each of the approaches we estimate the relationship between intercept or slope and age. For i = 1, . . . , 100 subjects, we simulate the response as: where the C i exposure values are simulated as in Section 6 of the paper following Wadsworth et al. [1] and the i are random errors simulated from a normal distribution with mean 0 and variance 0.02. The A i age values are simulated from four Uniform distributions such that there are 25 subjects in each of four age groups: 0 to 4 years; 4 to 10 years; 10 to 14 years; and 14 to 18 years.
First, the linear model with categorical covariates as shown in Section 4.1 is fitted to the simulated example data. Using the true age groups to define A 1i , . . . , A Hi (the age groups used for the categorical covariates), the following intercepts and slopes are estimated in turn for each of the four age groups: intercept estimates are 5.13, 4.80, 4.42 and 3.84; and slope estimates are -0.010, -0.031, -0.080 and -0.120. Now, we fit a single PALM tree model as described in Section 4.2. Figure S1 shows the results of a PALM tree fitted to simulated data; this is standard output from the 'partykit' package [2,3]. Four nodes (here, age groups) have been found in the following age groupings: 0 to 3.89 years; 3.89 to 9.94 years; 9.94 to 14.00 years; and 14.00 to 18 years. Other than 0 and 18 (fixed based on the paediatric population), the age group limits are observed age values from the data; were there an age data point less than 4, but closer to 4 than 3.89, this age boundary could be even closer to the truth. Regardless, these age groups are very Figure S1: Example to demonstrate the structure of a single PALM tree fitted to simulated data, produced from the 'partykit' package [2,3].
close to the true age groups and estimate the underlying PK-PD parameters well also. For each age group in turn, the intercepts are 5.13, 4.80, 4.42 and 3.84 and the slopes are -0.010, -0.031, -0.080 and -0.120. For this data, this model gives identical estimates to the linear model with categorical covariates, to six decimal places.
We then extend to the bootstrapped PALM trees also described in Section 4.2. Figure  S2 presents plots of the bootstrapped PALM fits of intercept and slope parameters over age constructed by following the approach given in Section 4.2. We plot the relationship between intercept or slope against age using the bootstrap averaged median, 2.5th and 97.5th quantiles at a continuum of ages from 0 to 18 years, also highlighting the true underlying intercept/slope values by green dashed lines. The 2.5th and 97.5th quantile lines are asymmetric as the the distribution (over the bootstrap samples) of the intercept / slope values is asymmetric at many age values from 0 to 18.
Next, we apply the B-splines approach described in Section 4.3. Figure S3 presents plots of the median, 2.5th and 97.5th quantiles of the posterior distributions of the intercept/slope parameters at each A i from the MCMC output of the B-spline model fit, also highlighting the true underlying intercept/slope values by green dashed lines.

Supplementary Appendix B: Inclusion of additional covariate
For all scenarios given in the paper, the response has been modelled as in equation (1) without additional covariates x 1i , . . . , x P i . In this supplementary appendix, the data are generated such that there is a relationship between response and an additional covariate.
Assume we have data on body weight, x w , which is modelled as a linear function of age; the linear relationship we use, x w = 3A + 7, is based on the use of weight estimation in paediatrics (1 to 13 years, inclusive) suggested by Luscombe et al. [4], though other suggested weight/age relationships exist. We assume that, like age, this covariate has an effect on response. However, as body weight here is largely explained by age it feels more natural to regress body weight against age and to consider the effect of the fitted residuals, r w , in the model, essentially, what effect of body weight on response remains after adjusting for age. We simulate the response in this example by having the body weight residuals, r w (which have a standard normal distribution), included in the simulation model as follows: When fitting the model, we consider two approaches: one approach where body weight has been observed and is included in the model; and a variation where the effect of the body weight residuals still exists, though body weight is now an unobserved covariate and not included in the model. The directed acyclic graph shown in Figure S5 aims to visualise this causal relationship. For this second modelling approach, we seek to identify how the methods cope when there is an effect that we are unable to observe and control for. We will therefore fit the models in Section 4 to this scenario, modelling the PD response in two ways.
Comparing the panels in Figures S6 representing the supplementary scenario with and without r w , it is clear that when the additional covariate is included in the simulation model, but not included in the analysis model, all approaches do not perform as well at estimating the underlying relationship between age and the exposure-response model slope or intercept parameters. However, the B-splines approach again seems to perform better than the other approaches in terms of accuracy and precision. Figure S7 shows that when the additional covariate is not included in the analysis model, the accuracy of the K-group optimal dosing rule is lower and the true expected response (derived from the simulation model when children are dosed according to the estimated optimal dosing rule) is further from the target response, for both the bootstrapped PALM trees and Bayesian penalised B-splines approaches. However, the Bayesian penalised Bsplines approach does seem to provide better accuracy than the bootstrapped PALM trees. Figure S8 shows that when the additional covariate is included in the analysis model, the majority of simulated datasets would lead to the investigator selecting a global optimum dosing rule with K * = 4, especially for the Bayesian penalised B-splines approach.            Figure S9: True underlying E-R relationships for scenarios 1 to 5, 7 and 8, as described in Appendix A of the paper.