Modeling Complex Pharmacokinetics of Long-Acting Injectable Products Using Convolution-Based Models With Nonparametric Input Functions.

The interest in the development and the therapeutic use of long-acting injectable (LAI) products for chronic or long-term treatments has grown exponentially. The complexity and the multiphase drug release process represent serious issues for an effective modeling of the PK properties of LAI products. The objective of this article is to show how convolution-based models with piecewise-linear approximation of the nonlinear drug release function can provide an enhanced modeling tool for (1) characterizing the complex PK profiles of LAI formulations with completely different drug release properties, and (2) addressing key questions supporting the optimal development of LAI products by simulating the PK time course resulting from different dosing strategies. Convolution-based modeling and simulation were implemented in NONMEM, and 3 case studies were presented to assess the performances of this new modeling approach using PK data of LAI products developed using different technologies and administered using different routes: microsphere technology and aqueous nanosuspension intramuscularly administered and biodegradable polymer subcutaneously administered. The performance of the convolution-based modeling approach was compared with the performance of conventional parametric models using a reference data set on theophylline. The results of the comparison indicated that the nonparametric input function provided a more accurate description of the data either in terms of global measure of goodness of fit (ie, Akaike information criterion and Bayesian information criterion) or in terms of performance of the fitted model (ie, the percent prediction error on Cma x and AUC0-t ).

The interest in the development and in the therapeutic use of long‐acting injectable (LAI) products for chronic or long‐term treatments has grown exponentially during the last decades. LAI products present a number of clinical advantages with respect to conventional oral therapies such as reduced dose frequency, enhanced adherence to treatment, avoidance of first‐pass metabolism, longer apparent body half‐life, and better control of the clinical response with overall improvement of the quality of life., The increasing prevalence of chronic disorders, such as schizophrenia, diabetes, cardiovascular diseases, HIV, and cancer, and the increasing demand for minimally invasive treatments have further fueled the growth of the LAI drug delivery products.

Different technologies have been used in the development of LAI products including injectable drug crystal suspensions, polymer‐based microspheres and polymer‐based or lipid liquid crystal in situ forming,, , , , biodegradable microsphere systems (eg, made of Poly Lactic‐co‐Glycolic Acid copolymer), and so forth. The common objective of an LAI product is to achieve optimal safety, efficacy, and patient compliance by controlling the drug delivery into the human body from 2 weeks to several months. The design of effective release characteristics is usually driven by the optimal drug pharmacokinetic (PK) profile identified by the exposure‐response relationship and by the minimum effective and maximal tolerated concentrations. To facilitate the development process of LAI products, it is fundamental to develop a pharmacometric framework aimed to (1) establish an in vitro/in vivo correlation (IVIVC) for supporting the optimization of the formulations, (2) estimate the expected exposure in a chronic treatment, (3) determine if a lead‐in oral treatment is required to reach the effective concentrations as quickly as possible, (4) determine the best switching strategy from a current to a new LAI treatment, (5) estimate the time necessary to clear the drug when the steady state has been reached, and (6) evaluate the risk of a drug holiday.

The availability of a model describing the time course of drug concentration profile is instrumental for implementing such a pharmacometric framework. The critical issue in this development is the high variable and irregular PK profile of the LAI products. Usually, this profile consists of an initial release phase (a burst release phase) when a small fraction of the drug rapidly moves from the injection site to systemic circulation, a lag phase with minimal drug release, and a main release phase when the large majority of the active ingredients is released from the injection site into systemic circulation.

The complexity and the multiphase in vivo drug release process represent a serious issue for an effective modeling of the PK properties of LAI products. Several pharmacokinetic strategies have been developed and applied to analyze such atypical absorption profiles, including the use of physiologically based pharmacokinetic approaches, double or triple Weibull in vivo release models,, parallel zero‐order immediate release and, after a lag time, first‐order release,, , transit compartments for delayed drug release, combination of immediate first‐order release and transit compartments, and inverse Gaussian density absorption. Recently, a convolution‐based modeling approach was shown to represent a powerful and flexible tool for modeling complex pharmacokinetics of extended‐release and LAI products, and for maximizing the benefit‐risk ratio of a treatment by optimizing the drug release properties using IVIVC and integrated PK/pharmacodynamic models.

Using this approach, the time course of the drug concentration can be described by convolving an input function with a disposition and elimination function when input and disposition functions are described by parametric models. A generalization of this method is now proposed for increasing the flexibility of the model and for providing an extended ability to characterize an irregular drug release process. The new implementation of a convolution‐based model includes a nonparametric description of the input function and a parametric description of the drug disposition and elimination processes. A piecewise linear approximation of the nonlinear input function is used to describe the drug release rate (ie, the first derivative of the cumulative drug release). By definition, the cumulative drug release is a monotonic increasing concave nonlinear function bounded to zero (no drug release at time zero) and to an upper asymptotic value (the total amount of drug released). The first derivatives of this function at different times can be approximated by a sequence of linear functions using a piecewise linear approximation. Piecewise linear models are widely used in diverse fields, such as circuit theory, image processing, curve fitting, system identification, economic and financial analyses, and so forth. The factors that presently motivate the use of these types of models are the simplicity of their structure, the extreme flexibility for characterizing highly variable and complex profiles, and the possibility to be efficiently implemented in nonlinear fitting programs.,

The objective of this article was to show how convolution‐based models with a nonparametric input function can be implemented in standard software such as NONMEM and to evaluate the performance of this modeling strategy by showing that (1) the performance of this model was as good as (if not better than) the performance of a conventional PK model in the population analysis of a reference data set, (2) this model was suitable to provide accurate fit of the PK data of different LAI products with completely different drug release properties, (3) this model was suitable to simulate the PK time course resulting from different dosing strategies, and (4) this model was suitable to address key questions for the optimal development of LAI products.

Methods

Convolution‐Based Model

The time course of the drug concentration resulting from an arbitrary dose can be described as a function of the in vivo drug release and the disposition/elimination processes defined by the unit impulse response according to the convolution integral: where τ is a dummy variable used for integration, Cp is the plasma concentration as a function of time t, f is the drug input rate, and UIR is the unit impulse response function.

The function characterizing the drug delivery f can be estimated as the first‐derivative of the cumulative drug release function r:

The convolution integral model (equation 1) can be represented in a more manageable form using a system of differential equations. In case of a simple disposition process (say, 1 compartment with first‐order process), the UIR function is characterized by the volume of distribution (V) and by the first‐order elimination rate constant (kel). The equation describing Cp(t) is: where Ap(t) is the amount of drug and F is the relative bioavailability of the current formulation with respect to the reference formulation (the one that provided an estimate of the UIR function). In this scenario, Cp can be estimated by numerically integrating equation 3. This model can easily be generalized to account for complex disposition processes including nonlinearity in the PK distribution and elimination processes. For example, in case of a drug presenting disposition and elimination characterized by a biexponential shape, the model defined in equation 3 can be replaced by equations 5 and 6 as: where k12 and k21 are the first‐order transfer rate constants from the central to the peripheral compartment and Ap and Ap1 are the amount of drug in the central and peripheral compartments, respectively.

The implementation of the convolution‐based model requires that one specify the submodel characterizing the r(t) function. This can be defined using a parametric or a nonparametric approach.

Parametric Input Function

In a parametric scenario, the structural form of r(t) is assumed to be described by a parametric function, such as exponential or single and dual Weibull functions with unknown parameters. In case of a double Weibull function, the r(t) function can be written as: where t is time, FF is fraction of the dose released in the first process, TD and TD1 are times to release 63.2% of the dose in the first and in the second processes, and SS and SS1 are sigmoidicity factors for the first and the second processes, respectively. The dr/dt function can be analytically estimated using the first derivative of equation 7 or can be approximated using a finite difference approach: where Δ is a sufficiently small number.

Nonparametric Input Function

In the case of a nonparametric scenario, the shape of dr(t)/dt can be directly approximated by a piecewise linear function using a sequence of parameters (pi, i = 1, n − 1) estimated on the n − 1 PK sampling times (ti) during the absorption process as

The pi values are constrained to be ≥0 as, by definition, r(t) is a cumulative function with positive first derivatives. The piecewise approximated input function can be theoretically defined by n − 1 parameters. In reality, not all parameters are relevant for an accurate description of dr/dt: some parameters may take the same value in 2 adjacent time intervals (ie, pi = pi + 1), indicating that the drug is released with the same rate in this time frame or some parameters (pi) may be equal to zero, indicating that no drug is released in the interval ti to ti+1. The modeling strategy is conducted using a top‐down approach: the full piecewise model with the n − 1 pi parameters will be initially implemented, and this model will be subsequently simplified according to the results of the nonlinear fitting. Therefore, some pi values may be fixed to 0, or the total number of the pi parameters may be reduced. The comparison of alternative models will be performed by inspecting the overall goodness‐of‐fit criteria.

The PK samples were assumed to be collected at the postdose times: 0, 0.25, 0.5, 0.75, 1.25, 1.5, 2, 3, 4, 5, 6, 8, 10, 12, 16, 20, and 24 hours. The nonparametric approximation of dr(t)/dt can be defined as: where p1 to p14 are the parameter values providing a piecewise local approximation of dr/dt (ie, the der values), dadt(.) is the symbol for derivative, and A(.) is the concentration of drug. The r(t) function is estimated as the integrated value of the second differential equation (ie, A (2)).

The values of A(2) normalized by the A(2) value at the end of the PK time course provided an estimate of the percent of the dose absorbed at the different times. The ‘der’ values provided an estimate of the rate of drug release and, in mean time, converted the amount of drug into predicted concentration (ie, the A(1) values). This is the reason why the volume of distribution (V/F), used as a scaling parameter to convert amount in drug concentrations in the parametric model, does not explicitly appear in this model.

The model can be equally applied to a single data set (ie, the individual or the mean observed PK) or to a population of subjects using a nonlinear mixed‐effects modeling approach. In this scenario, the full model will include the submodel characterizing the interindividual variability.

Modeling Repeated‐Dose Administration With a Convolution‐Based Model

The implementation of a convolution‐based model represented by equation 10 is only valid for single‐dose administration because the value of the dose appears explicitly in the code defining the model. In this scenario, NONMEM is unable to account for multiple doses using the default modeling options. Therefore, the cumulative concentrations resulting from administration of repeated doses need to be explicitly estimated by defining a superposition rule. This can be implemented using the DOWHILE functionality in NONMEM. An example of the implementation of a superposition rule using the DOWHILE instruction is presented below. This code assumes that the dosing history (the time of the dose intake and the dose at each intake time) is described in the input datafile using the “time” and ”amt” variables, respectively.

$input id time dv mdv amt

…

; Definition of the variables used in the DOWHILE code

; for implementing the superposition rule:

; dosetime (.) = time of dose intake

; dose (.) = dose value

$abbr declare dosetime(100), dose(100)

$abbr declare dowhile i

$abbr declare dowhile ndose

$abbr declare inpt

; End of the definition of the variables used in the DOWHILE code

$PK

; Evaluate the dosetime(.) and dose(.) variables

; by setting the time of the dose intake and the dose value

; from the value of the variable “time” and “amt” in the input data file

; when the value of the variable “amt” is >0

; The “ndose” variable represents the total number of doses administered

callfl = ‐2

if (newind < 2) ndose = 0

if (amt > 0) then

ndose = ndose+1

dosetime(ndose) = time

dose(ndose) = amt

endif

; End of the definition of the ndose, dosetime(.) and dose(.) variables

…

$DES

; Dowhile loop for estimating the cumulative input function (the “inpt” variable)

; resulting from the administration of each dose at a

; time tt (the relative time of dose intake)

; estimated as the difference between the cumulative time from

; the first dose intake (t) and the time from the i‐th dose (dosetime)

inpt = 0

i = 1

dowhile (i< = ndose)

tt = t‐dosetime(i)

if(tt> = 0) then

if (tt < 0.25) der = p1

if (tt ≥ 0.25 & time < 0.5) der = p2

if (tt ≥ 0.5 & time < 0.75) der = p3

…..

if (tt ≥ 10 & time < 12) der = p13

if (tt ≥ 12) der = p14

inpt = inpt+dose(i)*pp

endif

i = i+1

enddo

; End of the calculation of the cumulative input function

kkk = inpt

dadt(1) = kkk − kel · A(1)

dadt(2) = kkk

Comparison of the Model Performance With Parametric and Nonparametric Input Functions

The performance of the convolution‐based model with a nonparametric input function was compared with the performance of a conventional parametric model. A population modeling approach was used to analyze the theophylline data set using the 2 modeling scenarios. The data, accessible in R as “theoph” data frame, consists of 11 measurements of theophylline concentration in venous blood plasma collected in 12 subjects between 0‐ and 24 hours postdose.

The parametric model was a 1‐compartment model with first‐order absorption and elimination rate constants and with an absorption lag time. The model was defined by a set of 2 differential equations: where ka and kel are the first‐order absorption and elimination rate constants, respectively. The theophylline concentration (Cp) was estimated as Cp = A(2)/V.

The convolution‐base model with a nonparametric input function was implemented assuming that the drug disposition was described by a 1‐compartment model with a first‐order elimination rate constant. The nonparametric input function der = dr(t)/dt was locally approximated by 4 parameters:

The interindividual variability was assumed log‐normally distributed. The final model was defined by a single equation:

The theophylline concentration (Cp) was estimated as Cp = A1.

A global measure of goodness of fit is provided by the objective function value (OFV) based on the final parameter estimates. The OFV is estimated as minus twice the log likelihood of the data (−2LL) in NONMEM. Generally, models with more parameters are expected to better describe a given data set than models with a restricted number of parameters. Therefore, when comparing the performance of alternative models, it is necessary to adjust the comparison for the number of parameters used. The Akaike information criterion (AIC) and the Bayesian information criterion (BIC) were proposed for comparing model performance when the number of parameters and the structure of the models are different. These criteria are defined as follows: where np is the total number of parameters in the model and N is the number of the observations. Among 2 models, the most informative will be the one with the lowest AIC and BIC values.

In addition to these criteria, the performance of the 2 modeling approaches was evaluated by computing the percent prediction error (%PE) between the observed and individual predicted Cmax and AUC0‐t (area under the curve from time zero to the time of the last measurement estimated using a noncompartmental approach). The %PE was estimated using the following equation: where n is the number of subjects.

In the 2 modeling approaches, the interindividual variability was initially assumed to be lognormally distributed and the residual error proportional to the individual predicted values. Models were fitted to data using the first‐order conditional estimation with interaction method in NONMEM.

Case Studies: LAI Products for the Treatment of Schizophrenia

Several LAI formulations of the second‐generation antipsychotic drug risperidone are now broadly available for the treatment of schizophrenia. Risperidone is an atypical antipsychotic that possesses a high affinity for serotonergic 5HT2A and dopaminergic D2 receptors, and its pharmacodynamic profile allows for clinical efficacy in both positive and negative symptoms with a lower risk of extrapyramidal symptoms and tardive dyskinesia. Risperidone is mainly metabolized via CYP2D6 to 9‐hydroxy‐risperidone, which has pharmacologic activity similar to that of the parent compound. Therefore, the combined exposure of risperidone and 9‐hydroxy‐risperidone (the active moiety) has to be jointly evaluated for determining the effective/safe exposure resulting from the administration of an LAI product.

The time course profile of the active moiety resulting from the administration of 3 products (ie, Risperdal Consta, Invega Sustenna, and Perseris) was analyzed to illustrate the performance of the proposed modeling approach (Figure 1).

Figure 1

Typical pharmacokinetic profiles of the active moiety concentration resulting from the administration of 3 different LAI products for the treatment of schizophrenia.

Simulations were conducted for each LAI formulation to identify the best‐performing dosing regimen and to assess the need for either a lead‐in oral dosing period or for estimating loading LAI dose(s) to reach therapeutic exposure as quickly as possible. The target exposure for the active moiety was set to a range of 20 to 52 ng/mL, as estimated in a retrospective analysis conducted in 217 patients.

Risperdal Consta

Microsphere technology was used for this LAI formulation of risperidone. Risperidone was encapsulated in polymer microspheres of polylactide coglycolide (PLG) containing metabolic precursors of lactic acid and glycolic acid. PLG is a biodegradable medical‐grade polymer commonly used in sutures and medical devices. After intramuscular injection of the microspheres suspended in an aqueous diluent, the copolymer is gradually hydrolyzed and the microspheres progressively degraded, ensuring the slow but steady release of risperidone over several weeks. The mean PK samples generated after a single gluteal intramuscular injection of 50 mg of LAI risperidone to schizophrenic patients were used in the analysis. The PK samples were collected immediately before injection and at 1, 2, 4, and 8 hours on day 1, on day 2 (24 hours postdose), on day 3 (48 hours postdose), on day 5 (96 hours postdose), and on days 8, 11, 15, 18, 22, 25, 29, 32, 36, 39, 43, 50, 57, 64, 72, 78, and 85. The release profile of the active moiety showed a small initial release within the first 24 hours (<1% of the administered dose), followed by a lag time of about 3 weeks, during which very little drug was released from the microspheres. Therapeutic plasma concentrations were reached 3 to 4 weeks after injection, were maintained for an additional 2 weeks (throughout 6 weeks after injection), and subsided by 7 weeks after injection. The prescribing information for Risperdal Consta reported in the FDA labeling indicated that the recommended intramuscular dose is 25 mg every 2 weeks. Patients not responding to 25 mg may benefit from a higher dose of 37.5 or 50 mg. The maximum dose should not exceed 50 mg every 2 weeks. Oral risperidone (or another antipsychotic medication) should be given with the first injection of Risperdal Consta and continued for 3 weeks (and then discontinued) to ensure that adequate therapeutic plasma concentrations are maintained prior to the main release phase of risperidone from the injection site.

Invega Sustenna

Invega Sustenna (paliperidone palmitate [PP]) is the palmitate ester of paliperidone (9‐hydroxy‐risperidone). Paliperidone is a selective monoaminergic antagonist that exhibits the characteristic dopamine type 2 and serotonin (5‐hydroxytryptamine) type 2A antagonism of the second‐generation antipsychotic drugs. The LAI of PP for 4 weekly intramuscular injections has been approved in the United States as Invega Sustenna. The PP delivery system is based on an aqueous nanosuspension that slowly dissolves at the intramuscular injection site and releases paliperidone into systemic circulation over an extended period. The approved dosing is 234 mg on day 1, then 156 mg 1 week later (day 8), with a recommended maintenance dose of 117 mg intramuscularly once monthly, although some patients may require lower or higher dosages (monthly dose range, 39‐234 mg). The mean PK samples generated in a single‐dose, open label, randomized, parallel‐group study designed to evaluate the dose proportionality of 4 fixed doses of PP (25, 50, 100, 150 mg) following an intramuscular injection in the gluteal or deltoid muscle of schizophrenic patients were used in the analysis.

Perseris

Perseris is a once‐a‐month LAI formulation of risperidone subcutaneously administered using the Atrigel Delivery System. This biodegradable polymer drug delivery system uses an in situ poly‐dl‐lactide‐co‐glycolide implant formed by subcutaneous injection of a viscous liquid formulation that forms an implant on contact with tissue fluids. Risperidone is both dissolved and suspended in the Atrigel Delivery System and then slowly released from the injection site. Clinically relevant active moiety concentrations were reached after the first injection of Perseris without use of a loading dose or any supplemental oral risperidone. The mean PK samples generated after single subcutaneous doses of 60, 90, or 120 mg to schizophrenic patients were used in the analysis. Following injection, risperidone plasma concentrations increased rapidly with peak concentrations 4 hours postdose in all dose groups. After the first peak on day 1, mean plasma risperidone concentrations decreased through approximately day 3 and then increased again to reach a second peak between approximately day 11 and day 18. Following the second peak, risperidone plasma concentrations decreased gradually over time. The prescribing information for Perseris reported in the FDA labeling indicated that Perseris has to be administered as an abdominal subcutaneous injection at a dose of 90 or 120 mg once monthly. Neither a loading dose nor any supplemental oral risperidone is required.

Software

The data of Risperdal Consta, Invega Sustenna, and Perseris used in the analyses were extracted from the referred publications using ScanIt software, version 2.0. All simulations and parameter estimations were conducted using NONMEM software, version 7.4 (ICON Development Solutions, Hanover, Maryland). The data management and graphical presentation of the results were conducted using R language, version 4.0.0. The noncompartmental analysis was conducted using the validated open‐source R library “NonCompart.”

Results

Comparison of Model Performance With Parametric and Nonparametric Input Functions: Theophylline Case Study

The best‐performing model with a parametric input function was a 1‐compartment model with a lognormally distributed interindividual variability for all the parameters. The same 1‐compartment structural model was used for the convolution‐based model with a nonparametric input function. However, this model performed better when the interindividual variability of the pi parameters was assumed normally distributed (OFV, −304.963 with lognormal vs −306.025 with normal distribution). The residual error was best described by a proportional model in both scenarios. The parameter values estimated using the final models with parametric and nonparametric input functions are presented in Table 1.

Table 1

Parameter Values Estimated Using the Final Models With the Parametric and Nonparametric Input Functions

Parametric Input Function		Nonparametric Input Function
OFV = −287.321		OFV = −306.025
Parametera		Parametera
Fixed effect
kel (1/h)	0.0846 ± 0.0045 (5.30%)	kel (1/h)	0.0838 ± 0.0043 (5.1%)
Vd (L)	471 ± 19.3 (4.10%)	p1	2.4 ± 0.65 (27.1%)
ka (1/h)	2.3 ± 0.517 (22.50%)	p2	2.37 ± 0.457 (19.3%)
lag (h)	0.114 ± 0.0394 (34.60%)	p3	1.27 ± 0.304 (23.9%)
		p4	0.15 ± 0.0448 (29.9%)
Random effect
kel (1/h)	0.0292 ± 0.0231 (79.1%)	kel (1/h)	0.0339 ± 0.0235 (69.3%)
Vd (L)	0.0201 ± 0.0055( 27.3%)	p1	5 ± 3.36 (67.2%)
ka (1/h)	0.573 ± 0.246 (42.9%)	p2	2.13 ± 1.17 (54.9%)
lag (h)	0.441 ± 0.412 (93.4%)	p3	0.965 ± 0.496 (51.4%)
		p4	0.0287 ± 0.01 (34.8%)
Residual error
Proportional	0.0831 ± 0.0119 (14.30%)		0.0514 ± 0.0045 (8.8%)

Model estimates are reported as mean ± standard error (percent relative standard error).

The goodness‐of‐fit diagnostic plots for the 2 modeling approaches are shown in Figure 2. Overall, there was no apparent bias in these diagnostic plots, suggesting that either the model with a parametric or the model with a nonparametric input function provided an acceptable description of the data. The individual observations with the 2 model‐predicted curves are presented in Figure 3 in 4 selected subjects to illustrate the different PK time course resulting from the 2 modeling approaches.

Figure 2

Goodness‐of‐fit plots for the convolution‐based model with a parametric (a) and nonparametric (b) input function. The blue dotted lines represent the reference identity line (in the observed vs individual and population model‐predicted concentration value plots) and the reference zero intercepts horizontal lines (in the IWRES and conditional weighted residual plots). The red lines represent the smoothing function lines. The gray lines join the individual observations (circles).

Figure 3

Individual observations (dots) and model‐predicted curves in 4 selected subjects (red solid lines, parametric input function; blue solid lines, convolution model with nonparametric input function).

The overall performance of the 2 modeling approaches was assessed using the AIC and the BIC criteria. The AIC and BIC values for the parametric model and for the convolution‐based model with nonparametric input function were −269.321 and −284.025, respectively, for the AIC and −243.376 and −252.314, respectively, for the BIC. These values indicate that the convolution‐based model with nonparametric input function performed better than the parametric model. The estimated mean ± standard error %PE for AUC was 0.90 ± 0.18 and 0.60 ± 0.11 for the compartmental and convolution‐based model with nonparametric input function, respectively. The estimated mean ± standard error %PE for Cmax was 5.05 ± 1.18 and 2.51 ± 0.44 for the compartmental and convolution‐based model with nonparametric input function, respectively (Figure 4). The comparison of the %PE estimated in the 2 modeling approaches was conducted using a paired t test. The result of this analysis indicated a statically significant reduction of the %PE (P = .023) for Cmax and a borderline reduction (P = .052) for AUC, estimated using the nonparametric input function.

Figure 4

Estimated %PE on Cmax and AUC0‐t for the compartmental parametric model (Param) and for the convolution‐based model with nonparametric input functions (Non‐Param).

The disposition and elimination processes of the active moiety resulting from the administration of Risperdal Consta at the dose of 50 mg were characterized by a 1‐compartment model using the kel parameter defining the elimination rate constant. The shape of the input function (dr(t)/dt) was approximated by a piecewise linear function described by 14 parameters estimated between 0 and 40 days postdose (Table 2). The residual error was best described by a proportional model. The model provided an accurate description of the observed data (Figure 5A), and the estimated fraction of the available dose was characterized by a sigmoidal shape requiring ∼28 days to release half the administered dose (Figure 5B). The active moiety exposure simulated after the administration of 50 mg every 2 weeks indicated that the expected concentrations during the initial 30 days from the start of the treatment remained well below the therapeutic concentration range of 20 to 52 ng/mL (Figure 5C). Therefore, as recommended in the labeling of this product, a new simulation was conducted to supplement the 50‐mg LAI dose administered every 2 weeks with a daily oral dose of 5 mg administered during the initial 21 days. The results of the new simulation (Figures 5C, 5D) indicated that this new dosage regimen was suitable to provide a clinical benefit just on the first day of treatment. The time required for clearing the active moiety concentration (ie, the time for reaching an exposure below the quantification limit) was estimated at 70 days after the last drug intake on day 60, with the assumption of a quantification limit of 0.1 ng/mL.

Table 2

Parameter Values Estimated Using the Convolution‐Based Model With a Nonparametric Input Function for Risperdal Consta (50 mg), Invega Sustenna (100 mg), and Perseris (90 mg)

Risperdal Consta (50 mg)			Invega Sustenna (100 mg)			Perseris (120 mg)
Time Interval (Days)	Param	Value	Time Interval (Days)	Param	Value	Time Interval (Days)	Param	Value
< 0.042	p1	0.005	<2	p1	0.011	< 0.0.021	p1	0.731
0.042‐0.085	p2	2.300	2‐3	p2	0.042	0.021‐0.042	p2	0.881
0.085‐0.2	p3	0a	3‐5	p3	0.009	0.042‐0.083	p3	0.497
0.2‐2	p4	0a	5‐7	p4	0.005	0.083‐0.125	p4	0.306
2‐10	p5	0.006	7‐9	p5	0a	0.125‐0.167	p5	0.409
10‐12	p6	0a	9‐11	p6	0.002	0.167‐0.25	p6	0.061
12‐18.1	p7	0.009	11‐13.5	p7	0a	0.25‐0.333	p7	0.298
18.1‐22	p8	0.014	13.5‐15	p8	0.007	0.333‐0.5	p8	0a
22‐25	p9	0.034	15‐19	p9	0a	0.5‐1	p9	0a
25‐29	p10	0.115	19‐21	p10	0.008	1‐1.5	p10	0a
29‐32	p11	0.161	21‐28	p11	0a	1.5‐2	p11	0a
32‐36	p12	0.052	>28	p12	0.0003	2‐3	p12	0.035
36‐40	p13	0.026		kel (1/h)	0.022	3‐4	p13	0.019
>40	p14	0a		Error(Pr)	0.050	4‐7	p14	0.031
	kel (1/h)	0.165				7‐10	p15	0.046
	Error(Pr)	0.098				>10	p16	0a
							kel (1/h)	0.053
							k12 (1/h)	0.02
							k21 (1/h)	0.133
							Error(Pr)	0.048
							Error(Ad)	0.355

Error(Pr), proportional residual error; Error(Ad), additive residual error.

The pi parameters describe the piecewise linear approximation of the input function. All parameters were estimated with a precision (standard error) < 0.0001.

Ranges are ≥ to <.

Fixed.

Figure 5

(A) Mean observed active moiety concentrations (dots) with the predicted values by the convolution model (solid lines). The shaded areas represent the 95% prediction interval. (B) Model‐predicted percent of the available dose absorbed at the different times. (C) Simulated active moiety concentrations resulting from the administration of 50 mg every 2 weeks. (D) Two simulation scenarios: (1) simulated active moiety concentrations resulting from the administration of 50 mg every 2 weeks with the addition of a daily dose (4 mg) of oral risperidone during the initial 3 weeks of treatment (blue line); (2) simulated active moiety with discontinuation of the treatment after 5 administrations (green line); the vertical dotted line represents the time of the last administered dose on day 60. The horizontal dotted lines identify the range of effective active moiety concentrations (20‐52 ng/mL).

The disposition and elimination processes of the active moiety resulting from the administration of Invega Sustenna at doses of 25, 50, 100, and 150 mg were characterized by a 1‐compartment model. The shape of the input function (dr(t)/dt) was approximated by a piecewise linear function described by 12 parameters estimated between 0 and 28 days postdose. For brevity, only the parameters estimated at the dose of 100 mg are presented in Table 2. The residual error was best described by a proportional model. The model provided an accurate description of the observed data (Figure 6A), and the estimated fraction of the available dose was characterized by a sigmoidal shape requiring ∼7 days to release half the administered dose (Figure 6B). The active moiety exposure simulated after the administration of a monthly dose of 100 mg remained below the therapeutic concentration range of 20 to 52 ng/mL (Figure 6C). Therefore, as recommended in the labeling of this product, a new simulation was conducted at the dose of 234 mg on day 1, then 156 mg 1 week later (day 8) and with a maintenance dose of 117 mg once monthly. The results of this new simulation (Figure 6D) indicated that the new dosage regimen was suitable for providing a clinical benefit just on the first day of treatment and to maintain the clinical benefit during the treatment period. The time required for clearing the active moiety concentration (ie, the time for reaching an exposure below the quantification limit) is estimated at more than 360 days from the intake of the last dose on day 60, with the assumption of a quantification limit of 0.1 ng/mL.

Figure 6

(A) Mean observed active moiety concentrations (dots) with the predicted values by the convolution model (solid lines) by dose. The shaded areas represent the 95% prediction interval. (B) Model‐predicted percent of the available doses absorbed at the different times. (C) Simulated active moiety concentrations resulting from the administration of 100 mg every month. (D) Two simulation scenarios: (1) simulated active moiety concentrations at the approved dosing of 234 mg on day 1, then 156 mg 1 week later (day 8) followed by the maintenance dose of 117 mg intramuscularly once monthly (red line); (2) simulation as defined in point 1 but with discontinuation of the treatment on day 90 (blue line); the vertical dotted line represents the time of the last administered dose on day 60. The horizontal dotted lines identify the range of effective active moiety concentrations (20‐52 ng/mL).

The disposition and elimination processes of the active moiety resulting from the administration of Perseris at doses of 60, 90, and 120 mg were characterized by a 2‐compartment model with the following parameters: kel, elimination rate;, k12 and k21, transfer rate constants from the central to the peripheral compartments, respectively. The shape of the input function (dr(t)/dt) was approximated by a piecewise linear function described by 16 parameters estimated between 0 and 10 days postdose. For brevity, only the parameters estimated at the dose of 120 mg are presented in Table 2. The residual error was best described by a combined additive and proportional model. The model provided an accurate description of the observed data (Figure 7A), and the estimated fraction of the available dose was characterized by a sigmoidal shape requiring ∼4 days to release half the administered dose (Figure 7B). The active moiety exposure simulated after the administration of 120 mg every month indicated that the expected concentrations after the first dose reached values within the therapeutic concentration range of 20 to 52 ng/mL (Figure 7C). Therefore, as recommended in the labeling of this product, no loading dose was recommended in the current clinical practice. The time required for clearing the active moiety concentration (ie, the time for reaching an exposure below the quantification limit) was estimated at ∼100 days from the last dose intake on day 60, with the assumption of a quantification limit of 0.1 ng/mL.

Figure 7

(A) Mean observed active moiety concentrations (dots) with the predicted values by the convolution model (solid lines). The shaded areas represent the 95% prediction interval. (B) Model‐predicted percent of the available dose absorbed at the different times. (C) Two simulation scenarios: (1) simulated active moiety concentrations resulting from the administration of 120 mg every month (green line); (2) simulation as defined in point 1 but with discontinuation of the treatment on day 90 (red line); the vertical dotted line represents the time of the last administered dose on day 60. The horizontal dotted lines identify the range of effective active moiety concentrations (20‐52 ng/mL).

Discussion

Many factors are known to affect the ability to achieve the target clinical response of a therapeutic agent such as the identification of the safe and effective exposure, the characterization of the pharmacokinetic properties of the active ingredient, and the design of a dosage form suitable for delivering the active ingredient at the rate appropriate for quickly reaching and maintaining the targeted clinical response. Among these factors, the development of formulations using LAI technologies has been recognized as a strategic tool for the treatments of chronic and long‐term diseases. The main objectives of an LAI product are to achieve optimal safety and efficacy response to enhance adherence to treatment and to better control the clinical response with overall improvement of the quality of life. A critical step in the development of these products is the assessment of the relationship between exposure and clinical benefit and the evaluation of the therapeutic window: minimum effective and maximum tolerated concentrations. The most common routes of administration of LAI formulations are intramuscular, subcutaneous, and, less commonly, intravenous, intraocular, implant, and intra‐articular routes. The limitation of a subcutaneous injection is the limited dosing volume (ie, no more than 1 to 2 mL), whereas a larger injection volume can be administered for intramuscular (up to 2 to 5 mL) and intravenous (up to 100 mL) injections. The PK characteristics of a drug released by an LAI formulation is strongly affected by the physicochemical (solubility and stability) properties, the dose, the local absorption characteristics at the injection site, the injection volume, and the physiological properties associated with the diffusion of the drug from the administration site to systemic circulation.

All these factors significantly affect the development of models appropriate for characterizing the PK properties of LAI products. Recently, a convolution‐based modeling approach has been proposed as a powerful and flexible tool for modeling complex pharmacokinetics of LAI products. A generalization of this method is now proposed for increasing the flexibility of the model to characterize irregular drug release process resulting from the administration of LAI products via different routes. This new implementation is based on a piecewise linear approximation of the nonlinear input function used to model the drug release rate. In this article, 3 case studies are presented to assess the performance of this new modeling approach using PK data of LAI products developed using different technologies and administered using different routes: Risperdal Consta developed using microsphere technology and intramuscularly administered, Invega Sustenna developed using aqueous nanosuspension and intramuscularly administered, and Perseris developed using biodegradable polymer and subcutaneously administered.

To assess the overall benefit of this new modeling approach, the performance of the convolution‐based model with nonparametric input function was compared with the performance of a conventional parametric model using a reference data set on theophylline and a simple PK model. The results of the comparison indicated that the nonparametric input function provided a more accurate description of the data in terms of global measure of goodness of fit (ie, AIC and BIC criteria) or in terms of performance of the fitted model (ie, the %PE values on the Cmax and AUC0‐t).

These findings suggest that the approach based on the convolution‐based model with nonparametric input functions could represent a novel alternative methodology for fitting not only complex multiphase PK profiles but also more conventional PK profiles usually described by first‐order processes.

The availability of a general, flexible, and powerful modeling approach for characterizing the PK of LAI products is instrumental for supporting an informed development process. In particular, the model can be used for determining the best dosing strategy and to evaluate if and in which conditions an oral lead‐in treatment period should be supplemented to reach as quickly as possible efficacious effective concentrations, as illustrated by the simulations conducted for Risperdal Consta and Invega Sustenna. Furthermore, a PK model is instrumental for evaluating the timeline for loss of protective effect associated with drug holiday or discontinuation of LAI treatments. Clinical response in the presence of nonadherence to a treatment is difficult to evaluate in clinical trials for ethical reasons, but the resulting plasma levels can be simulated using model‐based tools. This analysis provides insight on the predicted plasma levels in clinical settings with the objective of aiding caregivers and patients to make informed decisions on treatment nonadherence. This is particularly critical for LAI products because of their mode of delivery and pharmacokinetics especially designed for providing a sustained and prolonged drug release. Complete discontinuation or interruption of treatment can have severe ramifications, including relapse and need for hospitalization. Therefore, it is important to estimate for how long the drug exposure remains above a minimal effective value after discontinuation.

A critical step in the development and in the optimization of LAI dosage forms is the assessment of an IVIVC as a toll for predicting, accurately and precisely, the PK time course of an LAI product from the dissolution data and for optimizing the development of formulations with optimal drug release properties. However, the assessment of IVIVC for LAI products is challenging because of the multiphase release characteristics and the lack of adequate in vitro release testing methods. Recently, the convolution‐based model implemented using parametric or nonparametric input function has been shown to represent a powerful and easy‐to‐use tool for assessing a time‐scaled level A IVIVC: a point‐to‐point correlation between in vitro dissolution and in vivo absorption. On these bases, the proposed convolution‐based model with nonparametric input function represents an effective tool first for estimating the fraction of the dose released in vivo and then for establishing an IVIVC conditional to the availability of the UIR function estimated using IV or immediate‐release data.

Evidence from clinical studies demonstrates potential clinical and economic benefits from early initiation of LAIs with respect to lower relapse rates, fewer hospitalizations, reduced illness‐related complications and comorbidities, and decreased medical resource use compared with oral products.

LAI products offer unique features and benefits, but also potential risks have to be considered within the context of each patient when selecting a specific LAI medication. Unmet needs remain for improved LAI formulations with optimal efficacy to support early use, less frequent injection for better patient comfort and convenience, and improved safety and tolerability profile by selecting the most effective doses and dosing strategy. In this context, the convolution‐based modeling approach has been shown to represent an integrated modeling framework for optimizing the clinical benefit of treatments by estimating the dosage regimen and the in vitro and in vivo drug release rates that maximize the expected benefit/risk ratio of treatments. An essential component of this approach is the availability of an accurate model describing the PK properties of LAI products. The nonparametric approach for describing the input function provides an easy‐to‐implement, flexible, and general tool for addressing these modeling issues that otherwise would require a trial‐and‐error approach to identify a parametric definition of the input function model. Finally, the traditional method for estimating clinical benefit is based on a sequence of trial simulation scenarios, each one aimed at exploring different features of the LAI formulations. The major limitation of this approach is that only a limited number of scenarios can be explored in simulations, which may be very time‐consuming. Furthermore, this empirical approach does not guarantee that the best‐performing scenario identified among the limited number of simulations will deliver the maximal possible clinical benefit.

In conclusion, the proposed modeling and simulation approaches have been shown to represent an effective framework for describing complex and multiphase PK of LAI products, for identifying the optimal dosing strategy, and for facilitating the development of LAI formulations. The limitations of the methodology remain associated with the model‐refining strategy that have to be implemented to properly characterize the piecewise linear approximation of the input function. In any case, the complex PK of LAI products requires consistent efforts for identifying the appropriate model suitable for describing the data whatever will be the modeling strategy selected.

Conflicts of Interest

The authors declare nothing to disclose and no conflicts of interest and that they did not receive any financial support for this study. R.G. and F.B.G. are consultants of Pharmacometrica.

Funding

This work was supported by PharmacoMetrica.

Data Sharing

The model codes can be obtained by emailing the corresponding author.