Remimazolam (CNS 7056) is a new benzodiazepine that was designed to deliver a fast onset, a short, predictable duration of sedative action, and more rapid recovery than currently available drugs such as midazolam.^{1} This pharmacodynamic profile is achieved by the rapid hydrolysis of the drug's ester group by nonspecific tissue esterases to its pharmacologically inactive metabolite, CNS 7054. The clinical results of a single ascending dose study in healthy volunteers, with midazolam as a comparator, described in our accompanying report,^{2} provide strong evidence that remimazolam is likely to meet its target profile when administered to patients undergoing procedures such as upper and lower gastrointestinal endoscopies.

To use the data captured in the study to simulate suitable doses and dosing regimens for remimazolam in clinical use, it was necessary to develop pharmacokinetic and pharmacodynamic models and use these to predict the effects of varied courses of therapy. We report the results of the modeling and include simulations of both Bispectral Index (BIS) and Modified Observer's Assessment of Alertness/Sedation (MOAA/S) scores after single and multiple doses of remimazolam.

METHODS
Study Conduct, Subjects, and Procedures
The clinical conduct of the study is described in our accompanying report.^{2} This was a single-center, double-blind, randomized, single ascending dose study of remimazolam administered as a 1-minute IV infusion, compared with midazolam and placebo. Eighty-one healthy male and female adults were treated; 9 subjects received placebo, 18 midazolam (0.075 mg/kg), and 9 groups of 6, remimazolam (0.01–0.3 mg/kg).^{2} This study was conducted at PAREXEL International, Early Phase Clinical Unit in Baltimore, Maryland, in accordance with the Declaration of Helsinki, in compliance with Good Clinical Practice, local, and Food and Drug Administration regulatory requirements. The clinical study protocol and informed consent forms were reviewed and approved by an IRB (Chesapeake Research Review, Columbia, MD) before study start. Written informed consent was obtained from all subjects before the start of the study.

During the study, sedation was measured by BIS monitoring^{3} and MOAA/S score assessments at prespecified times up to 2 hours postdose.^{2} Arterial pharmacokinetic samples were collected for plasma level measurements at various times between 1 minute and 4 hours postdose and venous samples were collected between 2 and 12 hours postdose^{2} ; there was a period of overlap at the 2-, 3-, and 4-hour sample collection time points when both arterial and venous samples were taken to allow for comparison of their respective plasma concentrations. Plasma concentrations of midazolam and remimazolam were measured using high-performance liquid chromatography with tandem mass spectrometric detection as described in the accompanying report.^{2}

Population Pharmacokinetic Modeling
A population approach was used to model the pharmacokinetic/pharmacodynamic behavior of remimazolam and midazolam. A nonlinear, mixed-effect modeling program (NONMEM, version 6, Nonlinear Mixed-Effects Modeling; Globomax, Ellicott City, MD^{4} ) was applied to the pharmacokinetic data using the first-order conditional estimation with interaction method (FOCE-I). Compartmental models (coded with a “C,” e.g., Rem-C01 for remimazolam) and physiologically based recirculation models (coded with an “R,” e.g., Mdz-R01 for midazolam) were developed using NONMEM's ADVAN6 subprogram with tolerance (TOL) set to 5 to solve the differential equations; these were expressed as clearances and volumes. Example control files are included in online Appendix 1 (see Supplemental Digital Content 1, http://links.lww.com/AA/A339 ). Intersubject variability of the various pharmacokinetic parameters was modeled assuming log-normal distributions. For example with elimination clearance:

where CL_{i} is the elimination clearance of subject(i), TVCL is the typical value of clearance (i.e., the population mean) and η_{i} represents the difference between the clearance of subject(i) and the population mean. The individual pharmacokinetic parameters were derived by NONMEM from the population mean and individual “η” values using equations such as the above, and, hence, the individual-predicted concentrations at each time point were calculated. A combined additive and exponential error model was used for the residual random effects. Standard errors of the parameters were obtained from the covariance step of each NONMEM run and used to calculate percentage relative standard errors, which are a measure of the precision of the estimates.

The choice of the final structural model was governed by the objective function (OFV). Although no formal criterion for the choice of model or number of compartments was set, reductions in objective function (ΔOFV) of at least 20 (i.e., ΔOFV less than −20) were considered necessary for a more complex model to be considered. Body weight, sex, and heart rate were examined in turn as possible predictors of each of the 8 pharmacokinetic parameters (Table 1 ) of the final remimazolam model. The OFV was the criterion used to determine whether the introduction of a particular covariate was statistically significant; reductions of 3.841, 6.635, and 10.828 were statistically significant at the 0.05, 0.01, and 0.001 levels, respectively.

Table 1: Summary of Some Population Pharmacokinetic Models of Remimazolam and Midazolam (Fig. 2)

Systemic exposure (AUC[0–∞]) was calculated as dose/elimination clearance. Systemic exposure (AUC[0–t]) was also estimated, up until the last time point, by means of the trapezoidal rule using the individual-predicted concentrations and assuming a linear change in concentration between adjacent data points. The terminal rate constant (k ) was calculated as the slope of the last 2 individual-predicted concentrations on a semilog plot of concentration against time and used to extrapolate AUC to infinity according to the equation:

Population Pharmacodynamic Modeling
A “link” model was used to relate plasma concentrations of remimazolam and midazolam to levels at the tissue “effect” site and a sigmoid inhibitory effect model described the influence of these levels on the BIS and MOAA/S scores:

where k _{e0} is the equilibration rate constant between the arterial plasma and effect compartments, “conc” is the concentration of drug in the relevant compartment, E_{0} is the baseline score, I_{MAX} is the maximum possible reduction in score, IC_{50} is the concentration of drug that causes the scores to decrease to halfway between E_{0} and (E_{0} − I_{MAX} ), and γ (the Hill coefficient) is a coefficient that describes the shape of the sigmoidal curve. The population pharmacodynamic analysis of the BIS data was performed in a similar way to that of the pharmacokinetics using NONMEM, ADVAN6, and an additive residual error. Because the MOAA/S scale is categorical, rather than continuous, ordered categorical models were fitted to the data, again using NONMEM, but with a logistic function and the conditional Laplacian method of estimation. The cumulative logit for each category was calculated and hence the cumulative probabilities that the MOAA/S score was less than or equal to that category.

where “x” represents a particular score and Baseline(x) describes the probability of achieving that score in the absence of the drug, I_{MAX} , the maximal achievable probability, and η represents any intersubject variability in the overall effect. The logits were converted into cumulative probabilities [e.g., P(x)] by means of the logistic transformation:

The probabilities for each category were obtained by subtraction from the cumulative probabilities, with the probability of achieving a score of ≤5 being unity.

Simulations
Preliminary Monte-Carlo simulations of the pharmacokinetics and pharmacodynamics of remimazolam after various IV infusion regimens were performed using ModelMaker (version 4; Pugh Computers, Aberystwyth, Wales) and the pharmacokinetic and pharmacodynamic parameters from the final models of remimazolam. The following protocol was observed for the simulations:

The initial loading dose was administered as a 1-minute IV infusion.
Maintenance doses, of either 15- or 60-second duration, were introduced if the MOAA/S score reached 4 during the 18-minute duration of the “procedure.”
If the initial loading dose had not reduced the MOAA/S score to <3, the first maintenance dose was given at that score.
If the MOAA/S score still did not decrease to 3, the subject was considered to be a dropout because of insufficient sedation for the initiation of the procedure.
A maximum of 4 maintenance doses was allowed, each no sooner than 2 minutes after the previous dose.
Subjects who had received 4 maintenance doses and whose MOAA/S scores had reached 5 during the 18 minutes of the procedure (excluding dropouts) were considered to be failures.
The variable parameters (i.e., body weight and the 5 pharmacokinetic and 4 pharmacodynamic parameters with between-subject variability [BSV] [Tables 1 –3 ]) were randomly selected from the same log-normal distributions as those found in the modeling. All 25 possible combinations of 5-, 6-, 7-, 8-, and 9-mg initial loading doses with 2.5-, 3.0-, 3.5-, 4.0-, and 4.5-mg maintenance dose regimens were examined and the preferred dosing regimen was selected on the basis of maintenance of a MOAA/S score of 2 to 4 in as many subjects as possible. One thousand subjects were simulated each time; the pharmacokinetic and pharmacodynamic parameters of the “subjects” were the same for each of the 25 dosing regimens.

Table 2: Summary of Final BIS Population Pharmacodynamic Models of Remimazolam and Midazolam

Table 3: Summary of Final MOAA/S Population Pharmacodynamic Models of Remimazolam and Midazolam

Monte-Carlo simulations of the pharmacokinetics and pharmacodynamics of remimazolam were then undertaken with the preferred initial loading dose/maintenance dose regimen using NONMEM. A dataset of 1000 “subjects,” which included an initial loading dose of 1-minute duration and time points of up to 40 minutes postdose, was first created. The plasma concentrations and pharmacodynamic scores for the 1000 putative subjects were simulated using the parameters from the final models, the variable pharmacokinetic and pharmacodynamic parameters being randomly sampled from the log-normal distributions found from the modeling. An ad hoc program in Excel's Visual Basic for Applications was used to identify when each subject's MOAA/S score had reached the value specified by the above protocol and an additional dose was inserted into the dataset. A second NONMEM simulation was then performed using the same individual pharmacokinetic and pharmacodynamic parameters as in the initial simulation. The process was repeated 3 times so that each subject could have received a total of up to 4 maintenance doses within the 18-minute hypothetical procedure. The final output was analyzed in Excel to estimate the distribution of MOAA/S and BIS scores at each time point.

Context-sensitive half-times were calculated for both remimazolam and midazolam by simulation using ModelMaker and the population mean pharmacokinetic parameters from the preferred models. Simulations of continuous 50 mg/h (remimazolam) and 0.075 mg/kg/h (midazolam) infusions of lengths between 1 minute and 8 hours were made, and the time for the concentrations to decrease to half their maxima was calculated. Context-sensitive half-times, defined as the time for the score to increase from its simulated minimum value to the average of the minimum and the resting figure, were also estimated for the BIS data. A similar calculation was made for the MOAA/S scores to estimate the relationship between the duration of the infusion and the time for the score to increase from 4 to 5.

RESULTS
Population Pharmacokinetic Analysis of Remimazolam and Midazolam
The raw data of each subject consisted of a mixture of 16 arterial (1 minute to 4 hours, n = 852) and 5 venous (2 to 8 hours, n = 184; all 12-hour samples were below the limit of quantification) concentrations, with 2 samples at the 2-, 3-, and 4-hour time points. A single clearly erroneous venous concentration was omitted from the analysis. Some demographic data of the subjects enrolled are included in the accompanying report.^{2} Inspection of the 2- to 4-hour remimazolam plasma levels and, to a lesser extent, those of midazolam showed that the venous concentrations were significantly higher (paired t test; P < 0.0001 for both compounds), the gradients of the plots of venous against arterial levels equaling 1.44 for remimazolam and 1.21 for midazolam (Fig. 1 ). Mammillary compartmental models (Fig. 2 ), which assume uniform distribution of the drug in a central compartment containing both arterial and venous plasma, were therefore flawed, and physiologically based recirculation models with separate venous and arterial compartments and a cardiac blood/pulmonary compartment for the dose (Fig. 2 ) were also examined. The volume of the cardiac/pulmonary compartment was assumed to be proportional to body weight (1 L/70 kg). Initial models also assumed that the arterial volume was proportional to body weight (0.65 L/70 kg), but significant improvements were obtained when this assumption was removed. Models with the elimination clearance of remimazolam from the central (Fig. 2 ), pulmonary, peripheral, and hepatic (volume = 1.5 kg/70 kg, blood flow = 26% cardiac flow) compartments were examined (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ).

Figure 1: Comparison of venous and arterial plasma concentrations of remimazolam and midazolam 2, 3, and 4 hours after IV administration to healthy volunteers. Gradient = venous:arterial concentration ratio.

Figure 2: Diagrammatic representation of typical (1) recirculation and (2) mammillary pharmacokinetic models of remimazolam and midazolam with clearance from the central compartment.

The best model, in terms of lowest objective function and meaningful pharmacokinetic parameters, was Run Rem-R13, which involved clearance from the peripheral compartment and BSV on elimination clearance and the arterial, venous, peripheral, and deep volumes (Table 4 , Fig. 2 ) (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ). Addition of BSV on a sixth parameter offered little or no improvement in objective function, and allowing both arterial and pulmonary volume to be determined by the modeling led to identifiability problems and physiologically unlikely estimates (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ). Three- and four-compartmental mammillary models, which used the same dataset and assumed that the observed differences between the arterial and venous levels were random errors, gave much poorer fits to the data (Table 4 ; e.g., ΔOFV = +370.1 for the 4-compartment model, Rem-C51, relative to Rem-R13).

Table 4: Comparison of Performance of Some Recirculation and Compartmental Pharmacokinetic Models of Remimazolam and Midazolam

The pharmacokinetics of midazolam, however, were not readily defined by the recirculation models, the covariance step of the run corresponding to the preferred model for remimazolam (Run Mdz-R13) failing, and both 3- and 4-compartment mammillary models gave considerably better fits (Table 4 ). For example, the objective functions of the successful recirculation model with a fixed arterial volume (Mdz-R07) and the 4-compartment mammillary model (Mdz-C51) had objective functions of 2262.0 and 2111.1, respectively (ΔOFV = +150.9, Table 4 ).

Covariate analysis of the possible effects of gender, body weight, and heart rate was examined for each remimazolam pharmacokinetic parameter. Physiologically reasonable findings that sex and body weight predicted cardiac output (on average, females weighed 10 kg less than males) were observed (ΔOFV approximately −20, Table 4 ) (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ). After the introduction of sex as a predictor of cardiac output, no other covariate was statistically significant at the 0.01 level (ΔOFV = −6.6). Body weight was not a significant predictor of elimination clearance (Run Rem-R22, ΔOFV = −4.9 relative to Rem-R13, online Appendix 2 [see Supplemental Digital Content 1, http://links.lww.com/AA/A339 ]) within the 65- to 90-kg weight range studied, and this can be seen by inspection of a graph of body weight against clearance, which consists essentially of a random distribution with a correlation coefficient for a power function of 0.09 (Fig. 3 ). The performance of the final models of both remimazolam and midazolam was good with average relative standard errors of the population mean pharmacokinetic parameters of 5.7% and 24.2% (Runs Rem-R13 and Mdz-C51, respectively, Table 1 ). Average between-subject variability was similar for the 2 compounds (32% and 29%) and residual errors were small (proportional errors: 12.9%, 7.0%; additive errors: 0.14 μg/mL, 2.5 μg/mL for remimazolam and midazolam, respectively). Goodness-of-fit plots for the remimazolam recirculation model confirmed the quality of its performance: plots of predicted and individually predicted concentrations lie close to the line of identity, whereas the individually weighted residuals are symmetrically distributed about zero (Supplemental Fig. 1, see Supplemental Digital Content 3, http://links.lww.com/AA/A341 ; see Appendix for figure legend). However, there is some positive bias among the later (venous) concentrations (Supplemental Fig. 2, see Supplemental Digital Content 4, http://links.lww.com/AA/A342 ; see Appendix for figure legend). A typical pharmacokinetic curve (subject #507, Fig. 4 ) shows the rapid decline in remimazolam concentrations and the well-fitted higher venous, compared with arterial, levels observed at the later time points. The calculated systemic exposure of each subject to remimazolam (Fig. 5 ) was similar to the results from noncompartmental analysis.^{2}

Figure 3: Plot of body weight against elimination clearance of remimazolam after IV administration at various doses (recirculation model R13). This graph is similar, but not identical, to that of the accompanying report,2 which was based on noncompartmental calculations.

Figure 4: Observed and predicted plasma concentrations of remimazolam after a 1-minute, 9 mg IV infusion to subject #507 (recirculation model Rem-R13).

Figure 5: Relationship of elimination exposure, calculated by noncompartmental methods to that estimated by recirculation model Rem-R13, of subjects to remimazolam after IV administration to healthy volunteers at various doses.

Population Pharmacodynamic Analysis of Remimazolam and Midazolam
Pharmacodynamic modeling, based on the parameters from the final pharmacokinetic models of remimazolam (Rem-R13) and midazolam (Mdz-C51), of both the BIS and MOAA/S data, was undertaken. The BIS observations were continuous data and sigmoid inhibitory effect models with BSV on all 4 pharmacodynamic parameters were readily fitted to the remimazolam and midazolam data (Table 2 ). A typical BIS curve (again for subject #507) is shown in Figure 6 ; although the fitted curves look quite realistic, the peak effect is not well modeled. A plot of individual-weighted residuals against time is consistent and indicates a general small negative bias at the 2- to 4-minute time points (Supplemental Fig. 3, see Supplemental Digital Content 5, http://links.lww.com/AA/A343 ; see Appendix for figure legend). The slope of the graph of predicted against observed BIS scores (0.77) also indicates bias, but the corresponding graph of individually predicted scores is symmetrical with a slope of 1.01 (Supplemental Fig. 4, see Supplemental Digital Content 6, http://links.lww.com/AA/A344 ; see Appendix for figure legend).

Figure 6: Observed and predicted Bispectral Index (BIS) scores after a 1-minute, 9 mg IV infusion of remimazolam to subject #507 (model Rem-R13-BIS01).

Covariate analysis of the effects of sex and age on IC_{50} , k _{e0} , and I_{MAX} , undertaken for remimazolam with model Rem-R13-BIS01, indicated that gender was a significant predictor of IC_{50} , with men apparently being twice as sensitive to the effects of the drug (IC_{50} : males = 0.21 μg/mL, females = 0.42 μg/mL, ΔOFV = −16.3, Table 5 ). No other covariate effects were statistically significant at the 0.05 level.

Table 5: Covariate Analysis of Pharmacodynamic Models of Remimazolam

The categorical MOAA/S data were modeled by means of ordered logistic regression. This introduced additional parameters corresponding to the baseline values of each category with the result that the models could be overparameterized. The MOAA/S models were initially explored using the pharmacokinetic parameters from remimazolam and the best then applied to the midazolam data. Model MOA02, with BSV on IC_{50} , I_{MAX} , k _{e0} , and the Hill coefficient, γ, was preferred for remimazolam (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ) and was also satisfactory for midazolam (Table 3 ). Covariate analysis, similar to that undertaken with the BIS data, used the simpler MOAA/S model MOA04 (BSV on IC_{50} , k _{e0} , and γ); it also indicated that sex was a predictor of IC_{50} , but the effect was smaller than for the BIS analysis and not statistically significant at the 0.01 level (Table 5 ).

The MOAA/S curves for subject #507 (Fig. 7 ) fitted the observed data well, and, overall, the percentage of exact correspondences in scores was 84.8% with 13.4% differing by 1 unit for model Rem-R13-MOA02 (Fig. 8 ). The corresponding figures for midazolam were 78.0% and 19.8% (Supplemental Fig. 5, see Supplemental Digital Content 7, http://links.lww.com/AA/A345 ; see Appendix for figure legends).

Figure 7: Observed and predicted Modified Observer's Assessment of Alertness/Sedation (MOAA/S) scores after a 1-minute, 9 mg IV infusion of remimazolam to subject #507 (model Rem-R13-MOA02).

Figure 8: Plot showing the correspondence between observed and predicted Modified Observer's Assessment of Alertness/Sedation (MOAA/S) scores after a 1-minute IV infusion of various doses of remimazolam (model Rem-R13-MOA02). The areas of the circles are proportional to the number of examples of each point.

Midazolam was apparently more potent than remimazolam with estimated population mean IC_{50} values of 0.07 and 0.08 μg/mL (BIS and MOAA/S models, respectively) compared with 0.26 and 0.42 μg/mL (Tables 2 and 3 ). However, the pharmacodynamic response was significantly slower as shown by the smaller midazolam k _{e0} values (0.053 and 0.050 min^{−1} for BIS and MOAA/S, respectively) compared with 0.14 and 0.25 min^{−1} for remimazolam (Tables 2 and 3 ). Very high population mean values of the Hill coefficient were apparent for midazolam, particularly for the BIS model (8.6, Table 2 ), but those for remimazolam were unremarkable (1.6, 1.4 for the BIS and MOAA/S models, respectively; Tables 2 and 3 ). The maximal effect of remimazolam was apparently larger than that of midazolam (BIS population mean I_{MAX} = 39.3 compared with 19.4 and MOAA/S logit = 21.9 compared with 9.5; Tables 2 and 3 ).

Simulations
Simulations of the arterial and venous plasma concentrations of remimazolam after administration to subject #507 (Fig. 9 ) show how the recirculation model initially gives increasing arterial levels with a peak being reached approximately 10 seconds after the end of the infusion and that the venous concentrations exceed the arterial ones at approximately 5 minutes. The concentrations in the 2 compartments gradually become closer, but the venous levels are always higher (Fig. 9 ).

Figure 9: Simulated arterial and venous plasma concentration profiles of remimazolam after a 1-minute, 9 mg IV infusion of remimazolam to subject #507 (model Rem-R13).

Context-sensitive half-times were simulated for both remimazolam (50 mg/h) and midazolam (0.075 mg/kg/h) using the population mean pharmacokinetic parameters (Table 1 ), after infusions of various lengths between 1 minute and 8 hours (Fig. 10 ). Because a mammillary model had been fitted to the midazolam plasma levels, there was no difference between the simulated arterial and venous concentrations. The context-sensitive half-times of midazolam increased with the duration of the infusion reaching a value of 60 minutes after an 8-hour infusion, in agreement with a literature figure of approximately 70 minutes.^{6} The recirculation model for remimazolam differentiated between arterial and venous plasma levels and the more relevant, as far as sedation is concerned, arterial context-sensitive half-times were simulated. The results were much shorter than those of midazolam and reached a steady-state value of between 7 and 8 minutes after an approximately 2-hour infusion (Fig. 10 ).

Figure 10: Comparison of context-sensitive half-times for remimazolam (arterial plasma levels, recirculation model Rem-R13) and midazolam (arterial and venous concentrations equal, mammillary compartmental model Mdz-C51). Midazolam doses = 0.075 mg/kg/h; remimazolam doses = 50 mg/h.

Context-sensitive half-times for the BIS scores were taken as the time for the score to increase from its minimum value to the average of the minimum and the baseline. Again, the simulated value for midazolam was long, at least 6 hours, after an 8-hour infusion (Fig. 11 ). In comparison, the maximal half-time for remimazolam was <1 hour and this figure was constant with infusions longer than 3 hours (Fig. 11 ). Because the MOAA/S data were categorical, a true context-sensitive half-time could not be estimated. However, the times for the score to increase from 4 to 5 after various lengths of infusion were considered useful parameters; a maximal value of 15 to 16 minutes was obtained with remimazolam infusions longer than 3 hours (Fig. 11 ). The 0.075 mg/kg/h dosing regimen for midazolam required more than an hour for any reduction in MOAA/S score. After 6- to 8-hour infusions, the time for the MOAA/S score to increase from 4 to 5 was 45 to 50 minutes.

Figure 11: Relationship between infusion length of remimazolam and midazolam and context-sensitive time for pharmacodynamic scores after typical doses. Half-time for Bispectral Index (BIS) scores equals time for score to increase from the minimum to mean of the minimum and baseline value. Midazolam doses = 0.075 mg/kg/h; remimazolam doses = 50 mg/h. MOAA/S = Modified Observer's Assessment of Alertness/Sedation.

Monte-Carlo simulations (1000) of a range of possible initial loading dose/maintenance dose regimens were undertaken as a guide for future studies. The distribution of the pharmacokinetic parameters was assumed to be the same as that found for recirculation model Rem-R13 and the pharmacodynamic parameters were obtained from MOAA/S model Rem-R13-MOA02. Some descriptive statistics of the variable pharmacokinetic and pharmacodynamic parameters from the Monte-Carlo simulations are shown in Table 6 . The geometric means of these parameters are approximately equivalent to the population means estimated from the modeling (Tables 1 , 3 , and 6 ).

Table 6: Summary of Pharmacokinetic and MOAA/S Pharmacodynamic Parameters Obtained from the Monte-Carlo Simulations of Remimazolam

The optimal regimen, in terms of minimizing the numbers of dropouts, failures, and subjects with MOAA/S scores of zero, seemed to be a 6-mg initial loading dose followed by 3-mg maintenance doses. This delivered sedation with MOAA/S scores between 2 and 4 to at least 70% of the subjects from 2 minutes after the end of the loading dose until the end of the 18-minute procedure (Fig. 12 ). Estimates of percentage of subjects suitably sedated ranged from 70.5%, 3 minutes after the start of the infusion (95% confidence intervals: 67%–74%), to >90% (95% confidence intervals: 92%–94%) between 4 and 8 minutes later. Sedation was rapid (Fig. 12 ) with all subjects predicted to attain their minimum MOAA/S score after a single dose within 2 minutes of the end of the infusion. Taking account of those simulated subjects who would need a second dose of remimazolam at the 2-minute time point before being sufficiently sedated for initiation of the procedure, and omitting the dropouts, minimum MOAA/S scores are predicted to be attained by the end of the “top-up” infusion in 85% of subjects (95% confidence intervals: 81%–89%). The 30% of the putative subjects whose simulated sedation was not optimal included dropouts (8.6% and 7.3% with 60- and 15-second maintenance doses, respectively) whose MOAA/S score did not reach 3 after 2 doses. The pharmacodynamic parameters of these subjects included either a large IC_{50} or a small I_{MAX} , i.e., the I_{MAX} /IC_{50} ratio for these subjects was small, the geometric mean equaling 22.3 compared with 50.2 for the whole population (Table 6 ). There were also failures (11.5% and 15.2% with 60- and 15-second maintenance doses, respectively) whose MOAA/S score increased to 5 before the end of the procedure despite receiving 4 maintenance doses; in these cases, the I_{MAX} /IC_{50} ratio was also small (geometric mean = 33.6), but larger than those of the dropouts. A few of the subjects (6.5%, 7.0% with 60- and 15-second maintenance doses, respectively) attained zero MOAA/S scores, i.e., loss of consciousness, usually 1 minute after the end of the initial loading dose; these subjects had higher than average I_{MAX} /IC_{50} ratios (geometric mean = 114.9). Simulated recovery from sedation was generally rapid (Fig. 12 ) with 89% (95% confidence intervals: 87%–91%) of the successfully treated subjects (i.e., excluding dropouts and failures) having attained MOAA/S scores of 5 within 16 minutes of the end of the procedure. The pharmacodynamics of the subjects who recovered more slowly were characterized by a small value of the Hill coefficient (γ); the geometric mean for the 42 subjects with scores of <5 at the 40-minute time point (i.e., 22 minutes after the end of the procedure) equaled 0.56 compared with 1.42 for the whole population (Table 6 ).

Figure 12: Monte-Carlo simulations of Modified Observer's Assessment of Alertness/Sedation (MOAA/S) scores during and after a 6-mg initial loading dose of remimazolam followed by up to four 3-mg maintenance doses (model Rem-R13-MOA02). Maintenance doses were given no closer than 2 minutes after the previous dose. They were introduced if the MOAA/S score reached 4 during the 18-minute duration of the “procedure.” If the initial loading dose did not reduce the score below 3, the first maintenance dose was given at that score. “‰” = frequency per thousand subjects.

DISCUSSION
Two population pharmacokinetic models were applied to the plasma levels of remimazolam and midazolam using the ADVAN6 subprogram of NONMEM: a conventional mammillary model with 3 or 4 compartments, and a more physiologically based recirculation model (Fig. 2 ).^{7} The assumption implicit in the simple compartmental models that the drug is instantaneously distributed throughout the central compartment is only valid with compounds whose distribution and clearance kinetics are relatively slow compared with blood flow. This was clearly not the case for remimazolam with its venous/arterial ratios being significantly >1 (Fig. 1 ). The physiologically based recirculation pharmacokinetic models overcome this deficiency by including both venous and arterial compartments as well as a compartment, consisting of cardiac blood and the lung, for the dose (e.g., online Appendix 1 [see Supplemental Digital Content 1, http://links.lww.com/AA/A339 ]). These models correctly fitted the observation that venous concentrations of remimazolam were significantly higher than arterial ones taken at the same time points between 2 and 4 hours after dosing and also indicated that the arterial levels will be higher up to approximately 5 minutes postdose (Fig. 9 ). This prediction has yet to be examined and it is probable that better models would have been found if both arterial and venous samples had been available within the first 15 minutes of the study. However, because of study logistics with a large number of procedures being repeatedly performed within a short time frame, venous samples were not started until the 2-hour time point. The mismatch between the venous and arterial plasma levels of remimazolam, which led to the recirculation model being investigated, was not anticipated, and it was assumed that early arterial concentrations would be sufficient to define the pharmacokinetics of both remimazolam and midazolam.

Although the recirculation model with clearance from the peripheral compartment displayed the lowest objective function, the absence of any venous data at early time points means that any conclusions about the site of the hydrolysis of remimazolam to CNS 7054 must be considered tentative. The degree of partitioning of remimazolam between red blood cells and plasma has not been determined, and “cardiac output” in the recirculation models (approximately 3.7 L/min) lies between the typical values of blood and plasma of approximately 5 and 2.5 L/min, respectively. The volumes of both arterial blood and the lungs/heart were initially assumed to be proportional to body weight to minimize any identifiability problems. However, it was found that allowing the arterial volume to be independent of body weight gave a significantly improved fit (models Rem-R07, Rem-R12, Rem-R13, online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ) with physiologically reasonable pharmacokinetic parameters (Table 1 : arterial volume approximately 1 L, BSV = 58%). On the other hand, if both the arterial and cardiac/pulmonary volumes were unfixed, physiologically unreasonable parameters resulted (online Appendix 2, see Supplemental Digital Content 2, http://links.lww.com/AA/A340 ) and fixing the volume of the lung to 1 L/70 kg was considered to be the best approach.

Whereas the quality of the recirculation models depended on the 2-, 3-, and 4-hour arterial concentrations being consistently smaller than the corresponding venous concentrations, the mammillary models required that they be similar. The relatively small, although statistically significant, difference between these midazolam plasma levels meant that the recirculation models did not minimize easily and the mammillary models (e.g., online Appendix 1 [see Supplemental Digital Content 1, http://links.lww.com/AA/A339 ]) were superior. The pharmacokinetic parameters obtained for midazolam were within literature ranges (elimination clearance: 0.30 ± 0.06 L/h/kg [0.25–0.54 L/h/kg^{5} ]; terminal half-life: 3.7 ± 0.3 hours [1.8–6.0 hours^{5} ]; volume of distribution: 1.59 ± 0.24 L/kg [1.0–3.1 L/kg^{5} ]).

The recirculation models were, however, much more appropriate for remimazolam with its greater elimination clearance. Thus, the difference in objective function between the same recirculation and mammillary models for remimazolam (Rem-R07/Rem-C03) and midazolam (Mdz-R07/Mdz-C03) were −326.4 and +156.0, respectively (Table 4 ). With its elimination clearance being 3 times that of midazolam (population means: 66.7 and 22.6 L/h) and its steady-state volume of distribution significantly less (population means: 89.0 and 121 L), both the terminal half-life and mean residence times of remimazolam were only 20% to 25% of those of its comparator (Table 1 ), indicating that it is likely to be a superior drug in terms of rapid recovery from sedation.

Other than a clear relationship between gender/body weight and cardiac output, no covariate effects were observed for the pharmacokinetic parameters of remimazolam. In particular, there was no significant relationship between body weight and elimination clearance (Fig. 3 ), which indicates that there will be no advantage in dosing healthy subjects by weight, rather than as a fixed dose, in terms of consistency of exposure to remimazolam within the weight range studied (65–90 kg).

No attempts were made to model the pharmacokinetics and pharmacodynamics simultaneously because the quality of the kinetic data was so much more reliable, and a 2-stage approach was used. The pharmacokinetic parameters for each subject were first estimated, using a mammillary compartmental model for midazolam and a recirculation model for remimazolam, and then used to simulate arterial concentration profiles for the pharmacodynamic modeling. Sigmoid inhibitory effect pharmacodynamic models were successfully fitted to the BIS and MOAA/S data resulting from treatment with remimazolam and midazolam. The BIS data are continuous and the models never exactly fitted the extreme values (i.e., E_{0} and E_{0} − I_{MAX} ). The MOAA/S scores, however, were categorical and the extreme values were often achieved by the models, particularly scores of 5 (full alertness). The effect of this is that the MOAA/S graphs appear much steeper (Figs. 6 and 7 ) even though the population pharmacodynamic parameters, k _{e0} , IC_{50} , and Hill coefficient, were broadly comparable (Tables 2 and 3 ).

Midazolam was apparently more potent than remimazolam with an IC_{50} of 0.07 (BIS), 0.08 (MOAA/S) μg/mL compared with 0.26 (BIS), 0.42 (MOAA/S) μg/mL for remimazolam (Tables 2 and 3 ). However, the Hill equation (see sigmoid inhibitory effect equation in Methods) assumes that the whole concentration/effect range is covered by the data. This can be assumed to be true for remimazolam, because complete anesthesia (MOAA/S scores of zero) was induced in 14 of 30 subjects given the higher doses (0.1–0.3 mg/kg) but not for midazolam, which was administered as a single 0.075 mg/kg dose (1 of 18 subjects with MOAA/S scores <1). As a result, the maximal effect (I_{MAX} ) appeared to be approximately 2-fold greater for remimazolam; higher midazolam doses would have defined its pharmacodynamic parameters better, increasing I_{MAX} and probably changing IC_{50} significantly. There was a marked difference in the values of k _{e0} , the rate constant for equilibration of the drug between plasma and the effect site; the population mean values for midazolam were approximately 0.05 min^{−1} compared with 0.14 min^{−1} (BIS) and 0.25 min^{−1} (MOAA/S) for remimazolam. The larger figure for remimazolam combined with its shorter pharmacokinetic half-life suggests that recovery from sedation will be considerably more rapid.

A limited covariate analysis of the pharmacodynamic parameters suggested that sex might be a predictor of sensitivity to remimazolam for both the BIS and MOAA/S models (Table 5 ) with males having approximately 40% smaller values of IC_{50} than females. It cannot be concluded from these fairly limited data that sex is a genuine predictor of sensitivity to remimazolam, but the possibility should be borne in mind in future studies.

The simulated population mean context-sensitive half-times also enable some comparisons to be made between remimazolam and midazolam. The relatively short pharmacokinetic half-time for remimazolam of 7 to 8 minutes compared with approximately 60 minutes for midazolam (Fig. 10 ) indicates that recovery from sedation or anesthesia will be much more rapid. In addition, the context-sensitive half-time of remimazolam appears to be relatively insensitive to the duration of the infusion, reaching its maximum after a 2-hour infusion compared with at least 8 hours for midazolam. A similar picture is provided by analogous simulations of the BIS and MOAA/S data (Fig. 11 ).

Preliminary Monte-Carlo simulations examined all 25 possible combinations of 5 initial loading doses (5–9 mg) followed by 5 smaller maintenance doses (2.5–4.5 mg) using ModelMaker, simulation software that could be programmed to recognize when a maintenance dose was needed and introduce it during an individual run. Additional Monte-Carlo simulations were then performed with NONMEM, using the combination that appeared to deliver the most consistent sedation of 2 to 4 on the MOAA/S scale. Single 1-minute, 6-mg infusions, followed by a 3-mg dose 2 minutes later if necessary, are predicted to lead to rapid sedation within 3 minutes of the start of treatment in >80% (95% confidence intervals: 81%–89%) of subjects but with loss of consciousness (MOAA/S = 0) in a small minority (approximately 7%, Fig. 12 ). The simulations also indicated that lower doses (e.g., 3 mg given over 15 or 60 seconds) at intervals no shorter than 2 minutes are adequate to maintain sedation in most subjects. The rapid onset of the effects of remimazolam means that maintenance doses will be able to be given more accurately than those of slower-acting drugs such as midazolam. The simulations highlighted subjects, characterized by a low I_{MAX} /IC_{50} ratio, who would be relatively insensitive to the effects of remimazolam. These included dropouts (approximately 8%) whose minimum MOAA/S score was >2 and failures (11%–15%) whose sedation was not controlled for the total length of the procedure with only four 3-mg maintenance doses. Because the choice of procedural length and number of maintenance doses was somewhat arbitrary, the latter category could be managed by increasing the number of allowable maintenance doses. The problem of the dropouts could probably also be overcome by permitting the first maintenance dose to equal the 6-mg initial loading dose if the MOAA/S score were still >3 at the 2-minute time point. Recovery from sedation resulting from treatment with remimazolam appears to be rapid, with almost 90% (95% confidence intervals: 87%–91%) of the successfully treated subjects predicted to be fully alert with MOAA/S scores of 5 within 16 minutes of the end of the procedure (Fig. 12 ).

CONCLUSIONS
Population pharmacokinetic and pharmacodynamic models developed for remimazolam and midazolam fitted the observed data well. Simulations based on these models show that remimazolam delivers extremely rapid sedation, with its maximal effect being reached within 3 minutes of the start of treatment. This property will enable maintenance doses to be given more accurately than slower-acting drugs. No covariate effects considered to be clinically relevant were observed, suggesting that dosing by body weight may offer no advantage over fixed doses in terms of consistency of exposure to remimazolam within the weight range studied (65–90 kg).

DISCLOSURES
Name: Hugh R. Wiltshire, PhD.

Contribution: This author performed the modeling and helped to prepare the manuscript.

Attestation: Hugh R. Wiltshire has approved the final manuscript.

Conflicts of Interest: Hugh R. Wiltshire has consulted for PAION UK Ltd.

Name: Gavin J. Kilpatrick, PhD.

Contribution: This author helped design the study.

Attestation: Gavin J. Kilpatrick has approved the final manuscript.

Conflicts of Interest: Gavin J. Kilpatrick was an employee and director of PAION at the time of conduct of the study. Dr. Kilpatrick owns shares and share options in PAION.

Name: Gary S. Tilbrook, PhD.

Contribution: This author helped analyze the data.

Attestation: Gary S. Tilbrook has approved the final manuscript.

Conflicts of Interest: Gary S. Tilbrook was an employee of PAION at the time of conduct of the study.

Name: Keith M. Borkett, BSc.

Contribution: This author helped design the study and prepare the manuscript.

Attestation: Keith M. Borkett has seen the original study data, reviewed the analysis of the data, approved the final manuscript, and is the author responsible for archiving the study files.

Conflicts of Interest: Keith M. Borkett was an employee of PAION at the time of conduct of the study.

This manuscript was handled by: Tony Gin, MD, FRCA, FANZCA.