Modeling and Simulation to Support Phase Ib/IIa Dose Selection for WBP216, A Long Half-Life Fully Human Monoclonal Antibody Against Interleukin-6

WBP216 is an innovative IL-6 antibody, presenting high affinity to IL-6 and a long half-life (40–60 days). To optimize the dosage regimen for future clinical trials, pharmacokinetics (PK) and pharmacodynamics (PD) of WBP216 would be firstly characterized in Chinese rheumatoid arthritis (RA) patients. PK, CRP and DAS28 data of WBP216 were collected from 26 RA patients in a single ascending dose study. Non-linear mixed effects modeling was used for a population PK/PD analysis. A two-compartment model with a sequential zero-first order absorption and a first order elimination best described PK behavior of WBP216. Apparent systemic clearance was 0.015 L/h, central volume was 8.04 L. CRP as the fast-decreasing endpoint and DAS28 as the slow-reacting endpoint were both fitted well through an indirect response model. The baseline of ALT and free IL-6 were found associated with PK/PD parameters during covariates exploration. Simulation results confirmed that a loading dose regimen either of administration at weeks 0, 2, and 6 or doubling the maintenance dose level, followed by maintenance dosing of 75–150 mg every 8 weeks, was expected to provide a best risk/benefit ratio in future clinical studies. We hope this first PK/PD study of WBP216 in Chinese RA patients will help in the clinical development of WBP216 in future and provide a reference to the dosage optimization of similar antibodies with long half-life. Clinical Trial Registration: CTR20170306


INTRODUCTION
WBP216 is an IgG1 antibody that binds and neutralizes IL-6 and is designed as a biologic antirheumatic drug. Rheumatoid arthritis (RA) is featured by progressive articular disability, systemic inflammation, and high morbidity, which stem from a complex interaction between various inflammatory cells and cytokines (McInnes and Schett, 2011;Smolen et al., 2016). Biological therapies are recommended for treating RA refractory to synthetic chemical drugs (Newsome, 2002). IL-6 is one of major pro-inflammatory cytokines involved in RA pathogenesis. Inhibition of the IL-6 signaling pathway helps to reduce inflammation and pain in patients with RA (Ohsugi and Kishimoto, 2008;Raimondo et al., 2017). Tocilizumab, as an anti interleukin 6 receptor (IL-6R) monoclonal antibody, is launched in 2009 for the indication of RA, whose successful clinical use has proved the key role of IL-6 in RA pathogenesis (Choy et al., 2002).
WBP216 is a fully humanized monoclonal antibody that has showed its significant efficacy in treatment of RA in phase Ⅰ trial. WBP216 has a strong affinity for interleukin 6 (IL-6), with an affinity constant in the picomolar range. WBP216 prevents the interaction of IL-6 and its receptor, thus reducing proinflammatory activity. Unlike regular IgG1 antibodies, WBP216 includes YTE mutations in the fragment crystallization (Fc). These mutations increase the ability of the Fc fragment to bind the neonatal Fc receptor (FcRn) (Oganesyan et al., 2009), which protect WBP216 from intracellular degradation, and thereby extending its half-life to >40 days and bringing some different pharmacokinetic (PK) characters. Thus, WBP216 has the potential to relieve RA effectively and for a sustained period.
Due to the blockade of IL-6 signaling, WBP216 can directly inhibit the hepatic production of C-reactive protein (CRP) directly (Bastida et al., 2019) and the erythrocyte sedimentation rate (ESR). Thus, the levels of inflammatory markers (CRP and ESR) rapidly decrease after the initiation of IL-6 blocker treatment, even prior to any improvement in swollen or tender joint counts. Therefore CRP alone is not sufficient to assess efficacy. According to the American College of Rheumatology (ACR) response criteria, the disease activity score based on 28 joint (DAS28) in its two versions employing ESR or CRP are most frequently used in clinical trials and also in clinical practice. DAS28 includes assessment of swollen/tender counts for 28 joints and patient/evaluator physician global assessment, and usually decreases slowly after RA treatments compared to CRP and ESR (Smolen et al., 2003). Hence, rapidly and slowly decreasing pharmacodynamic (PD) biomarkers should be assessed together to comprehensively evaluate RA medications.
Population PK/PD modeling and simulation has proven to be a useful approach in facilitating drug development (Budha et al., 2015). We report here the population PK/PD analysis of WBP216 PK and serum CRP and DAS28-ESR data in RA patients from phase Ⅰa, and model-based simulation results. The goal of this modeling and simulation is to optimize dose levels, dosing intervals and loading dose regimens for such a long half-life antibody in Phase Ⅰb/IIa, a multiple ascending dose (MAD) study, which will help to improve the drug development efficiency. Moreover, a population PK/PD analysis in this phase can lead to a better understanding of PK/PD characteristics of WBP216. PK/ PD data obtained directly from RA patients instead of healthy volunteers allows to recommend more accurate dosing regimens that can provide significant CRP and DAS28-ESR reduction while minimizing the frequency of subcutaneous injections. Therefore, we aimed to establish a population PK/PD model of WBP216 using phase Ia PK/PD data in RA patients and then use the model to make phase Ib/IIa dose regimen decision in this study. To our knowledge, this is the first report of the PK/PD of WBP216 in Chinese RA population. Our study will provide valuable information for the development and application of WBP216 in clinic.

Study Population and Data
Briefly, phase Ⅰa (CTR20170306) was a multi-center, randomized, double-blinded, single ascending dose study of WBP216 in RA patients. The study protocol was designed in full compliance with Good Clinical Practice and the Declaration of Helsinki and approved by the independent ethics committee of Peking Union Medical College Hospital and Beijing hospital (Beijing, China). Subjects eligible to take part in this study were RA patients diagnosed over 6 months, treated with basic RA medications (e.g., methotrexate and leflunomide) stably for at least 28 days, aged from 18 to 70 years and with a body mass index (BMI) of 19-30 kg/m 2 . All subjects signed the Informed Consent Form before their participation.
A total of 36 subjects were enrolled into 5 dose cohorts (10, 30, 75, 150 or 300 mg) respectively and randomized to receive WBP216 or placebo subcutaneously, with 3 subjects receiving active WBP216 and 1 subject receiving placebo in the lowest dose cohort (i.e., 10 mg) and with 6 subjects receiving WBP216 and 2 subjects receiving placebo in other each dose group. Serum samples were collected at pre-dose and 0. 083, 0.5, 1, 3, 7, 14, 21, 28, 42, 56, 84, 112, 168 days post-dose to obtain total WBP216 concentration and CRP data. Total WBP216 was assayed by validated methods whose linear calibration ranges were 39.1-10,000 ng/ml for WBP216. The CRP and ESR measurements were conducted by the clinical laboratory center of the hospitals. The LLOQ for CRP and ESR are 0.01 mg/L and 1 mm/h, respectively. DAS28-ESR and ACR scores were estimated by doctors at baseline and on days 7, 28, 56, 84, 112, 140 and 168, with 112/140 days measurements only for the 75-300 mg groups. If any anti-drug antibody (ADA) measurement was positive in any subject during the trial period, this subject was described as ADA positive.

Characterization of WBP216 PK Properties
The YTE mutations in the Fc fragment of WBP216 may result in non-standard extended PK properties, which is one of our study focus. In order to explore the PK characteristics, PK parameters were firstly calculated based on individual plasma concentrationtime-profiles using non-compartment analysis (NCA) by validated Phoenix WinNonlin version 8.1 software (Pharsight Corporation, Mountain View, CA, USA). PK parameters through NCA were then analyzed and explored using WinNonlin or PRISM (GraphPad 8.0.1, San Diego, CA, USA). A power model (Eq. 1) proposed by Gough et al. was applied to assess the dose proportionality (Kevin Gough, 1995).
in which Y denotes PK parameters such as area under curve (AUC) or maximum concentration (C max ). This approach assumes that the underlying relationship between log(Y) and log(dose) is linear. β 1 indicates total dose proportionality. In our study, a less stringent criterion was used (Hummel et al., 2009), given the small sample size in phase Ⅰ across multiple dose groups. The estimate of β together with CI falling completely within the range of [0.5-2] was quantified as dose proportionality.

Population PK and PD Model Development
The relationship between drug exposure and response was evaluated using nonlinear mixed effects models (Phoenix NLME, version 8.1, Certara). First-order conditional estimation, extended least squares method (FOCE-ELS) was used to estimate pop PK/PD model parameters. The final structural model was determined by the objective function value (OFV) and Akaike information criterion (AIC). A sequential modeling strategy was used for fitting the models to the phase Ⅰa data. The population PK model was first developed and then the individual post hoc parameters from the final PK model were used to predict the individual WBP216 concentrations to drive the drug effect on CRP or DAS28 time profiles using appropriate PD models. We tried one, two or three-compartment PK model with a first order, saturate or sequential zero-first order absorption compartment. The elimination phase was also analyzed by fitting to first order, saturate elimination or target-mediateddrug-disposition (TMDD) models. Based on the mechanism of action, an indirect-response model was chosen as a starting point for PD model development for CRP or DAS28. A linear model, an E max model or a sigmoidal E max model (Hill equation) were applied to characterize the relationship between WBP216 serum concentrations and those disease activity measures.
Inter-individual variability (IIV) was assumed to follow a lognormal distribution and was described using exponential model (Eq. 2).
where P ij represents the individual value of the parameter for the ith parameter in the jth individual, θ i depicts the population typical value for the ith parameter, and ηij represents random effect in jth individual sampled from a normal distribution with a mean of zero and variance of ω 2 . The residual unexplained variability of WBP216 concentration and PD observation was described by a proportional or additive error model, respectively (Eqs 3, 4).
where Y obs and Y pred are the observed and predicted serum concentration in plasma or disease activity measures. And ε 1 is the proportional and ε 2 is the additive component of the residual error model, respectively. Both of ε 1 and ε 2 are assumed to be normally distributed in the range of (0, σ 2 ). After collinearity diagnostics and correlation analysis, possible covariates including weight, age, sex, baseline serum albumin (ALB), alanine amiotransferase (ALT), creatinine clearance (CLcr), total bilirubin (TBIL), free IL-6, ADA (negative or positive) etc. were tested on both PK and PD parameters. Continuous covariates were described using the power function and categorical covariate were modeled by exponential function, see Eqs 5 and 6.
where Effect i is the multiplicative factor for covariate i, Cov ij is the covariate value for individual j, Cov median is the median covariate value, and θcov i is the exponent or parameter for covariate i model. Potential covariates were incorporated into the base model one by one using stepwise forward inclusion. When OFV decreased by 6.63 (at p < 0.01), the covariate was selected for inclusion to develop a full model, followed by the backward elimination. The covariates were subtracted one at a time in a stepwise manner as well once OFV increased above 10.83 (at p < 0.001, df 1), until all remaining covariates in model were statistically significant.

Model Evaluation
During the process of models building, the goodness of fit (GOF) of different models was compared on the basis of OFV and AIC. GOF was graphically evaluated by inspecting plots of the individual or population predicted vs. observed values, and plots of the conditionally weighted residuals (CWRES) vs. population predictions or time.
Models were also validated internally using predictioncorrected visual predictive checks (pcVPC) as well (Bergstrand et al., 2011). On the basis of 1,000 times pcVPC simulation, the 90% prediction interval (PI) was compared with the 90% interval of the prediction-corrected observations. Bootstrap analysis was also performed for the final model along with a total of 500 data sets resampling randomly from the original data set (Ette, 1997). We reported the calculated 90% confidence interval (CI) of model parameters from successfully minimized runs.

Simulations for Phase Ib Dose Selection
Simulation was conducted in Phoenix NLME (version 8.1, Certara, USA) based on a Monte-Carlo simulation approach. Up to 100 Phase Ib trials were simulated using the uncertainty distribution in parameter estimates. 27 patients in each simulated trial were simulated using the IIV log-normal distribution in both PK and PD parameters. The distribution of covariates still leveraged those of phase Ⅰa data set. Serum CRP levels and DAS28 were simulated for a range of maintenance doses (30-300 mg) under three different dosing frequencies: every 4 weeks (Q4W), every 8 weeks (Q8W), and every 12 weeks (Q12W). The duration of drug effect was simulated up to week 72 with weekly virtual PD sampling. It was assumed that the PK/PD relationship based on the Phase Ⅰ study lasting 24weeks could be extrapolated to longer term studies.
Tocilizumab has proven its successful clinical efficacy, so it was used as the reference for comparison of clinical endpoints. Since mean CRP is decreased by around 90% and mean DAS28 could be reduced by over 56.5% using the recommended dosage of tocilizumab (ACTEMRA ® HIGHLIGHTS OF PRESCRIBING WBP216 would take a long time to achieve steady state exposure and efficacy because of its slow elimination rate constants. Because RA patients would require a rapid pain relief, a loading regimen of WBP216 would be necessary. Hence, we simulated two categories of loading regimen to achieve steady state exposure: Firstly, WBP216 was given in a more intensive frequency at an initial three-administrations, including 0-4-8th week, 0-2-4th week or 0-2-6th week; The second simulation used a loading dose, that doubled the confirmed maintenance dose level.

Clinical Data Summary
Up to 36 RA patients took active medicine and placebo in a ratio of 3:1, respectively, wherein 27 patients received WBP216. One of the subjects in the 75 mg group showed a huge fluctuation of CRP level after 21 days from administration, very different from other subjects. The CWRES of there CRP samples were also greater than six during CRP model exploration. This subject was considered as an outlier and excluded from this PK/PD model analysis. Descriptive statistics (baseline values) of potential covariates of 26 patients tested in the PK/PD analysis were summarized in Table 1. A total of 391 PK samples were obtained during the phase Ⅰa and PD data consisted of 384 CRP samples and 241 DAS28 samples. Since fewer than 10% of PK samples were below the lower limit of quantification (LLOQ), they were handled by M1 method (Keizer et al., 2015). PK concentrations whose corresponding ADA was positive were all above the LLOQ, so they were not discarded.

Inspection of WBP216 PK Properties
We performed an NCA analysis for different dose groups before population PK model development to understand the PK properties fully because of the unique YTE mutations in WBP216. Figure 1 showed that mean apparent clearance (CL/ F) tended to increase with increasing dose, contrary to the clearance change pattern of TMDD, which generally has a high clearance in lower dose groups (Mager, 2006). Large variability of CL/F among individuals was observed in both 75 and 300 mg dose groups. High individual variability of apparent volume (V/F) was also observed in higher dose levels (75-300 mg). Like CL/F, V/F presented an increasing trend over dose levels. The phenomenon of CL/F and V/F changing with dose levels was speculated to result from the changes of either the actual increased CL and V or decreased bioavailability (F). To clarify the real reason, the half-life of WBP216 was analyzed. Half-life (t 1/2 0.693*V/CL) is considered not to be impacted by F. In this analysis of half-life, the value of t 1/2 distributed evenly in five dose groups, and the mean t 1/2 remained almost the same in different groups, except for the 150 mg group, likely caused by lower CL/F and higher V/F value in this group. Those changing trends of CL/F, V/F and t 1/2 suggested that the decrease of F with doses was probably the cause of the increase of CL/F and V/F. Mean t 1/2 is 40-60 days, indicating the potential for a long dosing interval in therapeutic use. The power model fitted ln (AUC) over ln (dose) data well ( Figure 1D). Parameter β was calculated as 0.652 with 90% confidence interval (CI), [0.43, 0.87], part of which fell outside the lower limit of prespecified range [0.5, 2] (Hummel et al., 2009). This dose-dependent study indicated little lower than proportional increases in exposure (AUC).

Population PK/PD Model Development
According to AIC value and goodness of fit (GOF) plots, the selected final proposed PK/PD model structure is shown in Figure 2. The WBP216 serum concentrations were best described by a two-compartment PK model with sequential zero-first order absorption and first order elimination (see Eqs 7-9). An indirect-response model with an drug E max inhibition of the CRP or DAS28 zero order rate production constant (K in ) best described the disease measures (see Eqs 10, 11). Residual variability was characterized by a proportional error model for serum concentration, CRP and DAS28. All PK/PD model parameters are summarized in Table 2.
Wherein A depot , A1, A2 representing WBP216 amounts in absorption depot, systemic central and peripheral compartments respectively, were equal to zero when t 0. C1 and C2 were the concentration in central and peripheral compartment, respectively and equal to zero before dosing. Response was CRP or DAS28, and was equal to its baseline value before drug administration. The final PK model fitted the observed concentration data better by assuming that after WBP216 was subcutaneously administered, drug was released from the injection site at a FIGURE 2 | The proposed PK/PD model structure. Parameters are abbreviated as follows: SC subcutaneous injection; Td time of zero order release; Ka first order absorption rate; C1 concentration in central compartment; V1 central volume; CL clearance from central compartment; C2 concentration in peripheral compartment; V2 volume of peripheral compartment; Q clearance between central and peripheral compartments; E max maximum drug efficacy; EC 50 concentration needed for 50% of E max ; K in zero order rate constant for response production; K out first order rate for response loss.
Frontiers in Pharmacology | www.frontiersin.org February 2021 | Volume 12 | Article 617265 zero order rate to a depot compartment, and was then absorbed to a central compartment. Td is the period of zero order release and Ka denotes the first order absorption rate. The reason of selecting such a more complex absorption model will be discussed later. Based on the mechanism of action, WBP216 neutralizes IL-6 and inhibits CRP production, and further slow down disease progression, so it is plausible to set the drug inhibition on response production rate, K in (Sharma and Jusko, 1998). As shown in Table 2, the relative standard errors (RSE%) for almost all fixed-effect parameters were ≤30.9%. The uncertainties for random-effect parameters were <27%.
Overall, the precision of parameters estimates was acceptable. Parameters showed various inter-individual variability, ranging from 17 to 128.4%. E max for CRP was set to 1 since its value was always very close to 1 in all tested runs. E max for DAS28 was also fixed to 1, thus leading to a straightforward convergence of the model. All applicable shrinkage was below 18.1%, which was smaller than reported cut-off value 30% and could assure accurate IIV estimates and avoidance of misleading diagnostic plots (Savic and Karlsson, 2009). The correlation diagnosis chart of various covariates was presented in Supplementary Figure S1. The potential impacts of demographics and laboratory data baseline on PK/PD parameters of WBP216 were tested using a stepwise covariate modeling procedure. Those statistical significant covariate effects were identified and retained in the final model: ALT on CL, baseline free IL-6 on EC 50,CRP and K in,CRP (see Eqs 12-14 and Supplementary Table S2). Continuous covariates were described using the power function, centered by the median value. Apparent clearance decreased with increasing ALT with exponent −0.833, which explained around 18% of CL/F inter-individual variability (IIV). Higher baseline free IL-6 levels could result in increased K in,CRP (exponent 0.695), while lead to its  Table 2.

Model Diagnosis and Evaluation
Goodness-of-fit (GOF) plots for the final PK model in serum are shown in Figure 3. Plots of the population-and individualpredicted concentration vs. observed concentrations demonstrate no major bias. The conditionally weighted residuals (CWRES) were symmetrically distributed about zero axis and most points laid within the acceptable range (−2 to 2), suggesting that little to no bias accompanied with concentration or time. The GOF plots for CRP and DAS28 model were presented in Supplementary  Figures S2 and S3, respectively, which also performed well in visual diagnostic. The predictive performance was evaluated internally by pcVPC. Plots of pcVPC were presented in Figure 4. We can see that the 5th, 50th and 95th percentiles of prediction-corrected observations and predicted data were fairly consistent, especially a better match for DAS28. The 95th percentile of predicted PK data through VPC is a slight under-prediction and the CRP model over-predicted drug inhibition at 5th percentile slightly, which will be discussed in the part of discussion. Despite these small deviations, the 90% prediction interval of simulated data covered most of the observations. The bootstrap results are also shown in Table 2. The typical values of parameters and IIV estimates in the final model were pretty close to the median values from bootstrap validation, and fell within 90%CI of bootstrap parameters completely, which indicated high stability and precision of the final model.
Overall, the good performance of GOF, pcVPC plots and bootstrap estimations reconfirmed that the final PK/PD model was adequately developed and the predictive performance was sufficient to capture PK/PD observations.    (Broeder et al., 2002;Wells et al., 2009) (see Supplementary Table S1): moderate responders were patients with an improvement of >0.6 and a present DAS28 score of <5.1, or an improvement of >1.2 and a present DAS28 score of >5.1. We predefined that once patients start having a moderate response to WBP216, they would begin to feel pain relief intuitively. Patients' time needed to reach at least moderate response are summarized in Figure 6. The results showed that over half of patients got pain relief within around 10 days, while 12/26 patients did not reach moderate response until at least 28 days from drug administration. WBP216 caused a gradual increase in clinical efficacy during the phase Ⅰa period, which was consistent with the patients' reports of slow onset of drug action.
Therefore, in order to allow patients to benefit from WBP216 therapy more quickly, loading dose regimens were designed. On  the basis that the optimal dosing frequency was Q8W during maintenance therapy, we simulated four kinds of "loading dose schedules", Figure 7 showed examples of varying loading schedules followed by a maintenance dose level of 150 mg Q8W: (1) Loading dose at weeks 0, 4 and 8 followed by Q8W still took more than 8 weeks to reach steady state exposure.
(2) Loading dose at weeks 0, 2 and 4 followed by Q8W caused a steep WBP216 concentration increase, exceeding the peak concentration achieved by the maximum dose level (300 mg) in phase Ⅰa, which may raise potential safety concerns.
(3) Dosing at weeks 0, 2 and 6 followed by Q8W would reach plateau exposure at the second drug administration. (4) A loading dose of 300 mg, doubling the maintenance dose level, resulted in steady state instantly at the first dosing. Hence, the last two loading dose schedules are highly recommended for future clinical trials of WBP216. The simulated DAS28 profiles under the optimized loading dose regimen (3) and (4) were presented in Supplementary Figure S4, showing that the target efficacy was reached three weeks in advance.

DISCUSSION
We developed a population PK/PD model that characterized the relationship between WBP216 serum concentrations and the changes in CRP and DAS28. This is the first report of the PK/ PD of WBP216 in Chinese RA patients. Rheumatoid arthritis is a complex autoimmune-mediated inflammatory disease involving both genetic and environmental factors (Smolen et al., 2016). Many factors in RA patients such as the disease process, complications, and concomitant medications could cause significant PK and PD differences between patients and healthy volunteers. Different from some IL-6 (R) antibodies that have been developed PK/PD models based on first-in-human studies with data from healthy volunteers (Xu et al., 2011), our modeling and simulation study in RA patients directly could avoid the aforementioned confounding factors and recommend more accurate dose regimens for phase Ib/IIa trials.
A two-compartment PK model with sequential zero-first order absorption and first order elimination fitted the PK data well. During the process of characterizing PK properties, we found the bioavailability of subcutaneous administration decreased gently with increasing dose levels. This atypical absorption may because large antibodies from the SC site are transported through the tissue interstitium and into the lymphatic system slowly accompanied by tissue metabolism and hydration etc (Richter et al., 2012;Jung et al., 2018). Adopting a base PK model structure with sequential zero-first order absorption (OFV 4,935) performed significantly better than a traditional twocompartment model with first order absorption rate (OFV 5,088) or a saturated absorption model of adding the Michaelis-Menten equation to bioavailability (OFV 5,026). Although the NCA analysis suggested a slight trend of absorption saturation, it seemed not obvious enough to fit a typical M-M equation best in absorption phase. Linear clearance was sufficient to describe WBP216 elimination rather than a linear plus Michaelis--Menten elimination of tocilizumab (Abdallah et al., 2017), probably because WBP216 binds to free IL-6 while tocilizumab binds the IL-6R involving complex internalization. The typical CL/F (0.015 L/h) from model estimates was slightly smaller than CL/F of sirukumab (0.019 L/h, corrected by F) (Xu et al., 2011), a IL-6 antibody, but the estimated V1/F (8.1 L) of WBP216 seemed to be larger than 4.1 L in sirukumab and 5.6 L in tocilizumab. Therefore WBP216 has a longer half-life. PK parameters showed high inter-individual variability, however, only ALT was identified as a covariate accounting for 18% of the variability of CL/F. CL/F was negatively correlated with ALT with exponent −0.833. Chunze Li etc. also reported that clearance of trastuzumab emtansine correlated significantly with baseline albumin and AST (Li et al., 2017). And the clearance of another IL-6 antibody, siltuximab, was found to be impacted negatively by ALT as well in its poplulation PK analyis containing 460 participants (Nikanjam et al., 2019). It is generally accepted that a therapeutic antibody is unlikely to be impacted by functional hepatic impairment (Dirks and Meibohm, 2010). We did not have sufficient data to explain this phenomenon. Dose adjustment for hepatic dysfunction was undetermined, and this would require a large target population to validate this point further. Some studies reported that weight was an important factor affecting CL/F (Abdallah et al., 2017) or V1/F (Li et al., 2018). Our study did not find this covariate, perhaps because of limited sample size with narrow weight range (61.5 ± 7.9 kg) in this phase Ⅰa trial and the stringent covariate entry/elimination criteria set. In addition, only 3/27 patients exhibited ADA positive (Table 1), so ADA was correspondingly identified as a nonsignificant covariate influencing pharmacokinetic behaviors. All ADA samples of the outlier subject were detected negative, which could not impact ADA conclusion after exclusion.
We did not adopt the strategy of modeling the PK and PD data concurrently because the high variability of CRP or DAS28 would affect the estimate precision of PK parameters. An indirect response model with inhibition on response production was applied for both CRP and DAS28 endpoints based on WBP216's mechanism of action. Fast decreasing CRP is able to reflect the binding IL-6 ability of WBP216 directly since hepatic production of CRP is mainly mediated by IL-6 (Vermeire et al., 2004), while slow changing DAS28 tracks RA-related clinical efficacy closely. Unlike some early clinical studies that only focused on a fast-decreasing biomarker (i.e., CRP) (Xu et al., 2011;Mayer et al., 2015;Li et al., 2018), we evaluated both fast-and slow-decreasing endpoints simultaneously to understand the efficacy of WBP216 comprehensively. As Table 2 shown, whatever K in,CRP or K out,CRP had greater value than those of DAS28, hinting high turnover rate of CRP. K in,DAS28 0.003 h −1 and K out,DAS28 6*10 −4 h −1 were very closed to reported values of tocilizumab (K in,DAS28 0.0037 h −1 ; K out,DAS28 7.2 × 10 −4 h −1 ) (Bastida et al., 2019). While Levi et al. reported K in,DAS28 0.011 h −1 and K out,DAS28 15.8 × 10 −4 h −1 using data from 4 phase III studies of tocilizumab, which was almost 2 fold of our estimated DAS28 parameters (Levi et al., 2013). The differences in above reports may because the DAS28 baseline of subjects was about 5.3 in our and Bastida's Frontiers in Pharmacology | www.frontiersin.org February 2021 | Volume 12 | Article 617265 studies instead of 6.8 in Levi's research. The EC 50,CRP (194.37 ug/L) for CRP was fairly smaller than that of DAS28 (1,576.3 ug/L), which indicated a higher concentration was needed for half-maximally decreasing DAS28. Covariates analysis showed that baseline free IL-6 was expected to be positively correlated with K in,CRP , since CRP production is mainly stimulated by IL-6 in body. However, no covariates were discovered affecting K in and K out of DAS28. Baseline of free IL-6 was negatively associated with the EC 50 for CRP. The addition of baseline free IL-6 in model was able to explain 59.2% and 55.6% of IIV for K in,CRP and EC 50, CRP respectively. Those limited sample size and narrow demographic or laboratory data may not provide fairly accurate covariate impacts on PK/PD parameters but they provided a reference for future covariate analysis in larger population.
The final population PK/PD model was evaluated by GOF plots, pcVPC and bootstrap. Slight deviations were observed in the 95th percentile of serum concentration and 5th percentile of CRP between respective predicted and prediction-corrected observations ( Figure 4). It was noteworthy that one of patients in 300 mg group had extremely high PK exposure, over three fold of other subjects, which raised the 95th percentile of the observed concentration significantly. WBP216 was able to inhibit CRP to a very low level, almost close to zero, while the lower limit of quantitation of CRP can only reach 0.01 mg/ L. Our model is expected to have a higher uncertainty near zero because of high IIV and residual errors, so the model predicted a lower CRP 5th percentile value than that of observed data. However, this slight deviation did not affect overall predictive ability since the 90% prediction interval of simulated data covered the majority of the observations. WBP216 is the IL-6 monoclonal antibody with the longest half-life (40-60 days) by far. By comparison, the half-life of tocilizumab is reported as 11-13days and siltuximab as around 21day, leading to a dosing regimen of once every three or four weeks in clinical practice (ACTEMRA ® and SYLVANT ® , HIGHLIGHTS OF PRESCRIBING INFORMATION; FDA). According to our simulation results ( Figure 5), the long halflife of WBP216 would allow it to be optimally administered once every 8 weeks. Dosing every 4 weeks did not offer any advantage since excess drug exposure would occur. Although the simulation results showed 150 mg Q12W seemed also acceptable, the Q8W dosing frequency was able to maintain a more stable change in CRP and DAS28. A score of DAS28-ESR<2.6 defines RA remissions (Anderson et al., 2012;Smolen et al., 2016). In our simulations with the baseline DAS28 around 5.3, when dosing 75 mg Q8W and 150 mg Q8W, approximately 81% and 92% of virtual patients were able to achieve DAS28-ESR <2.6 after 24 weeks therapy, respectively. So both of dose levels (75 and 150 mg) will deserve to be tested in future clinical studies.
A few subjects in phase Ⅰa complained of getting limited relief until one or two months after drug administration ( Figure 6). This prompted us to design four loading schedules used for simulation. An initial loading regimen of dosing at weeks 0, 2 and 6, followed by a maintenance regimen of Q8W, achieved steady state at the second administration, which was consist with dosage regimen of infliximab (RENFLEXIS, HIGHLIGHTS OF PRESCRIBING INFORMATION; FDA). Another good loading dose option is to double the maintenance dose level. The simulated DAS28 profiles under the two optimized loading dose regimen did show that three weeks were saved to reach the target efficacy with the maintenance dose of 150 mg Q8W (Supplementary Figure S4).
Although ACR20/50/70 endpoints were also measured in this early clinical study, there was no obvious dose-dependent relationship in the probability of achieving ACR20/50/70 efficacy in this phase Ia study. The combination of both early biomarkers (CRP and DAS28) helped to build the PK/PD relationship. A population PK/PD approach proved again useful in integrating all available PK/PD data during early clinical phases. Our study provided a relatively complete paradigm to accelerate clinical development for similar drugs as well.
The limitations of the model are as follows: (1) Baseline ALT and free IL-6 considered to be statistically significant were observed as covariates. Nevertheless, the limited number of patients in our study and the strict entry criteria set for the clinical trial may result in ambiguous covariate effects. Therefore, the confirmation of covariate effects should be kept in mind in future studies, which will have more data of target population added in. (2) In order to design dose regimen for a MAD study, it was assumed that the PK/PD relationship based on 36 patients in the Phase Ia study lasting 24weeks can be extrapolated to a larger target population and a longer term study. Attention should be paid to this model hypothesis when drawing conclusions from the simulation results. Patients with different disease states, disease progression, drug resistance and combination may invalidate this hypothesis.
In summary, a population PK/PD model was first successfully established for WBP216. Fast-decreasing (CRP) and slow-decreasing (DAS28) biomarkers were modeled concurrently to assess efficacy of WBP216 fully. For WBP216 with an exceptionally long half-life (40-60 days), two kinds of loading dose regimens are recommended for the next clinical studies. We expect that the modeling and simulation will be valuable for dose selection during future clinical trials, and provide a reference for the PK/PD studies of similar antibodies.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by correspondence authors, without undue reservation.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by peking union medical college hospital; Beijing 2016-I2M-1-010). Those funders were not involved in the study design, collection, analysis, interpretation of data, the writing of this article or the decision to submit it for publication.