Dynamics of target-mediated drug disposition: characteristic profiles and parameter identification

In this paper we present a mathematical analysis of the basic model for target mediated drug disposition (TMDD). Assuming high affinity of ligand to target, we give a qualitative characterisation of ligand versus time graphs for different dosing regimes and derive accurate analytic approximations of different phases in the temporal behaviour of the system. These approximations are used to estimate model parameters, give analytical approximations of such quantities as area under the ligand curve and clearance. We formulate conditions under which a suitably chosen Michaelis–Menten model provides a good approximation of the full TMDD-model over a specified time interval.


Introduction
The interaction of ligand and target in the process of drug-disposition offers interesting examples of complex dynamics when target is synthesised and degrades and when both ligand and ligand-target complex are eliminated. In recent years such dynamics has received considerable attention because it is important in the context of data analysis, but also, more generally, in the context of system biology because this model serves as a module in more complex systems [1].
Based on conceptual ideas developed by Levy [2], the basic model for target mediated drug disposition (TMDD) was formulated by Mager and Jusko [3]. Earlier studies of ligand-target interactions go back to Michaelis and Menten [4]. We also mention ideas about receptor turnover developed by Sugiyama and Hanano [5]. Mager and Krzyzanski [6] showed how rapid binding of ligand to target leads to a simpler model, Gibiansky et al. [7] studied the related quasi-steady-state approximation to the model and Marathe et al. [8] conducted a numerical validation of the rapid binding approximation. Gibiansky et al. [9] also pointed out a relation with the classical indirect response model. For further background we refer to the books by Meibohm [10] and Crommelin et al. [11], and to the reviews by Lobo et al. [12] and Mager [13].
In practice, the Michaelis-Menten model is often used when ligand curves exhibit TMDD characteristics (see e.g. Bauer et al. [14]). Recently, Yan et al. [15] analysed the relationship between TMDD-and Michaelis-Menten type dynamics. We also mention the work by Krippendorff et al. [16] which studies an extended TMDD system which includes receptor trafficking in the cell.
The characteristic features of TMDD dynamics were first studied in [3] under the condition of a constant target pool, i.e., the total amount of target: free and bound, was assumed to be constant in time. Under the same assumption, a mathematical analysis of this model was offered by Peletier and Gabrielsson [17]. This assumption was made, in part for educational reasons, because it makes a transparent geometric description possible, in which qualitative and quantitative properties of the dynamics can be identified and illustrated. In a recent paper Ma [18] compared different approximate models under the same assumption of a constant target pool.
In the present paper we extend this analysis to the full TMDD model and do not make the assumption that the target pool is constant. This means, in particular, that we shall now also be enquiring as to how the target pool changes over time and how it is affected by the dynamics of its zeroth order synthesis and first order degeneration.
The analysis in this paper consists of a combination of numerical simulations based on a specific case study, and a detailed mathematical analysis of the set of three differential equations that constitute the full TMDD model. Our main objective will be to gain insight in such issues as: 1. Properties of concentration profiles When the initial ligand concentration is larger than the endogenous receptor concentration, the dynamics of target-mediated drug disposition results in a characteristic ligand versus time profile. In Fig. 1 we show such a profile schematically. One can distinguish four different phases in the dynamics of the system in which different processes are dominant: (A) a brief initial phase, (B) an apparent linear phase, (C) a transition phase and (D) a linear terminal phase. We obtain precise estimates for the duration of each of these phases and for each of them we obtain accurate analytical estimates for the concentration versus time graphs of the ligand, the receptor and the ligand-receptor complex.
On the basis of ligand concentration versus time curves we will develop instruments for extracting information about the target and the ligand-target complex versus time curves.  [7,8,15]. In this paper we point out how the full model and the reduced Michaelis-Menten model differ significantly in the initial second order phase and in the linear terminal phase, in that the terminal rate (k z ) of ligand in the full model is much smaller than that in the Michaelis-Menten model.

Parameter identifiability
In this paper we focus on the classical TMDD model, as presented in for instance [3] and shown schematically in Fig. 2. However, since we will focus on the typical features of the interaction between the ligand L, the receptor R and their complex RL, we drop the peripheral compartment (V t , Cl d ) of the ligand. In mathematical terms, the model then results in the following system of ordinary differential equations:  In the first phase (A) drug and target rapidly equilibrate, in the second phase (B) the target is saturated and drug is mainly eliminated directly by a first order process, in the third phase (C) the target is no longer saturated and drug is eliminated directly, as well as in the form of a drug-target complex, and in the final, fourth phase (D) the drug concentration is so low that elimination is a linear first order process with direct as well as indirect elimination (as a drug-target complex) Fig. 2 Schematic description of target-mediated drug (or ligand) disposition. The ligand L binds reversibly (k on /k off ) to the target R to form the ligand-target complex RL, which is irreversibly removed via a first order rate process (k e(RL) ), and in addition is eliminated via a first order process (k e(L) = Cl (L) /V c ) The quantities L, R and RL are assumed to be concentrations, k f = In L /V c denotes the infusion rate of the ligand (here In L denotes the infusion of ligand and V c the volume of the central compartment), and k on and k off denote the second-order on-and first-order off rate of the ligand. Ligand is eliminated according to a first order process involving the rate constant k e(L) = Cl (L) /V c , where Cl (L) denotes the clearance of ligand from the central compartment. Ligand-target complex leaves the system according to a first order process with a rate constant k e(RL) . Finally, receptor synthesis and degeneration are, respectively, a zeroth order (k in ) process and a first order (k out ) process.
In the absence of a zeroth-order infusion of ligand, i.e., when k f = 0, the steady state of the system (1) is given by Thus, this is the situation when there is no free or bound ligand. The receptor concentration then satisfies a simple turnover equation involving zeroth order synthesis and first order degeneration: The TMDD system can be viewed as one in which two constituents, the ligand and the target, or receptor, are interacting with one another whilst each of them is supplied and removed, either in their free form, or in combination in the form of a complex. We shall find that the total quantities of ligand and receptor will play a central role. Therefore, we put We deduce from the system (1) that their behaviour with time is given by the following pair of conservation laws: Note that in this system the on-and off rates of ligand to receptor no longer appear. We investigate the dynamics of the system (1) that is generated by two types of administration of the ligand: (i) Through a bolus dose. Then k f = 0. We denote the initial ligand concentration by L(0) = L 0 = D/V c , where D is the dose and V c the volume of the central compartment.
(ii) Through a constant rate infusion In L . Then k f = In L / V c [ 0 and L(0) = 0.
When we assume that prior to administration the system is at baseline, the initial values of ligand, receptor and ligand-receptor complex will be In this paper we focus on the situation when k e(RL) , k off and k e(RL) are small compared to k on R 0 . Whereas in our earlier investigation [17], the total amount of receptor was constant, because it was assumed that k e(RL) = k out , here this assumption is no longer made and, generally, R tot will vary with time. However, we show that there exists an upper bound for the total amount of receptor in the system, free and bound, which is independent of the amount of ligand supplied and holds for both types of administration. Specifically, we prove that, starting from baseline, For the proof of this bound we refer to Appendix 3.
Anchoring our investigation in a case study in which ligand is administered through a series of bolus doses, we dissect the resulting time courses of the three compounds, L, R and RL and identify characteristic phases, Phases A-D shown in Fig. 1. We associate these phases with specific processes and show, using singular perturbation theory [19,20,21], that individual phases may be analysed through appropriately chosen simplified models, yielding accurate closed-form approximations. They offer tools which may be used to compute critical quantities such as residence time, and to verify whether different approximations to the full TMDD model, such as the rapid binding approximation [6] and the quasi-steady-state approximation [7,18] are valid in the different phases. These issues are discussed at the conclusion of this paper.
Much of the mathematical analysis underpinning the results presented throughout the text is presented in a series of Appendices at the end of the paper.

Case study
The ligand, target and complex concentration-time courses used throughout this analysis, were simulated to mimic real experimental observations obtained on a monoclonal antibody dosed to marmoset monkeys. Simulated data were generated with a constant coefficient of variation (2 %) of the error added to the data. Synthetic data were used for pedagogic and proprietary reasons in order to answer the question: ''To what extent is the parameter precision affected by including/not including target (R) and complex (RL) data. The model used for generating the data is shown in Fig. 2 and the actual parameters are stored in Table 1. WinNonlin 5.2, with a Runge-Kutta-Fehlberg differential equation solver, was used for both simulating and regressing data. A constant CV (proportional error) model was used as weighting function. All dose levels (Concentration-time courses) were simultaneously regressed for the ligand (L), ligand and target (L, R) and ligand, target and complex (L, R and RL) data analyses, respectively. Kinetic data of high quality-as regards spacing in time and concentration and very low error level-were used intentionally to demonstrate (a) improved precision when using two or more sources of chemical entities; (b) as a check that one gets back approximately the same parameter estimates as used for generating the original dataset(s).

Data analysis
The case study is based on three sets of simulated concentration versus time data (I, II and III), each set of data obtained after following four rapid intravenous injections of the ligand or antibody (L). These datasets are shown in Figs. 3 and 4. They are increasing in richness: the first set (I) contains ligand profiles only, the second set (II) contains ligand as well as target or receptor profiles and the third set (III) contains profiles of all three compounds: the ligand, the receptor and the ligand-receptor complex.
The purpose of this study is to demonstrate the possibility of fitting the eight-parameter model shown in Fig. 2, to three different sets of high quality data with increasing richness, and show how precision of the estimates of the model parameters increases when successively information about target (II) and target and complex (III) is added. We use this data set for two purposes: (i) for data analysis and (ii) for highlighting critical features of the temporal behaviour of the three compounds.
Simulated data from three sources (ligand, target and complex) were intentionally used. We have experienced that data of less quality gave biased and imprecise estimates as well as biased and imprecise predictions of ligand, target and complex.
The central volume V c was assumed to be equal to 0.05 L/kg and fixed. The other parameters are then re-estimated   Fig. 2 to the data. Estimates are given in Table 2 0.001 Left semi-logarithmic graphs of simulated plasma concentrations of L (red discs) and R (blue squares) versus time (Data set (II)) and on the right the same, but also semi-logarithmic graphs of RL (green triangles) (Data set (III)), taken after administration of four rapid intravenous injections D of 1.5, 5, 15 and 45 mg/kg, respectively. V c for these doses was fixed at 0.05 L/kg. The dots are simulated data and the solid curves are obtained by fitting the model sketched in Fig. 2 to the data. Estimates are given in Table 2 for datasets I-III (Table 2). Thus, we know a priori what values they should have. Dataset I is made up from simulated concentration-time profiles covering five orders of magnitude in concentration range and from 0 to 500 h. Dataset II contains the same simulated ligand (L) profiles as in dataset I as well as target (R) concentration-time profiles obtained at each dose level. Dataset III includes dataset II but is enriched by four simulated time-courses of the ligand-target complex (LR) as well.
The four doses are D = 1.5, 5, 15 and 45 mg/kg. The volume V c of the central compartment being 0.05 L/kg, this yields the following initial ligand concentrations L 0 = D/ V c = 30, 100, 300 and 900 mg/L.
The parameter values given in Table 1 yield the following values for the dissociation constant K d and the constant K m related to K d : and Here, K d is a measure of affinity between drug (ligand) and target, whereas K m is more of a conglomerate of affinity (K d ) and irreversible elimination of the ligand-target complex (k e(RL) ) and used for comparisons to the Michaelis-Menten parameter K M of regression model Eq. (45). Unless the removal of the ligand-target complex is fully understood, one should be careful about the interpretation of an apparent K m -value. K m can be very different from the affinity K d . Here, solid biomarker (physiological or disease markers) data on effective plasma concentrations may be a practical guidance.
Summarizing we may conclude that: • Dataset I-which involves L-allows the prediction of robust ligand concentration-time profiles within the suggested concentration and time frame. We see that if only ligand data are available, the majority of parameters except for k on , k off and k e(RL) are estimated with high precision. The latter three parameters are still highly dependent on information about the time courses of either target and/or complex. Since k kon , k off and k e(RL) have low precision (high CV%) we would discourage the use of these parameters for the prediction of tentative target and complex concentrations. • Dataset II-which involves L and R-still gives good precision in all parameters except k off and k e(RL) , which will also be highly correlated. Since we also have experimental data of the target we encourage the use of this model for interpolation of target concentrationtime courses, but not for concentration-time courses of the complex. • Dataset III-which includes L and R, as well as RL- gives high precision in all parameters. Since we also have measured the complex concentration-time course with high precision we obtained k off and k e(RL) values with high precision. We doubt that the practical experimental situation can get very much better than this latter case where we have simultaneous concentration-time courses of L, R and RL with little experimental error due to biology and bio-analytical methods. Dataset III is an ideal case; the true experimental situation seldom gets better.
We also doubt the practical value of regressing too elaborate models to data. Models that capture the overall trend nicely but result in parameters with low precision and biased estimates may be of little value.
The volume of the central compartment V c ought to fall somewhere in the neighbourhood of the plasma water volume (0.05 L/kg) for large molecules in general and antibodies in particular. In our own experience of antibody projects this has been the case when data contained an acceptable granularity within the first couple of hours after the injection of the test compound. Therefore we assumed V c to be a constant term (0.05 L/kg) in this analysis and not part of the list of parameters to be estimated. We think this increases the robustness of the estimation procedure and is biologically viable.

Critical features of the graphs
The graphs in Figs. 3 and 4 exhibit certain characteristic features and so reveal typical properties of the dynamics of the TMDD system.
(a) Initially all the ligand graphs in Fig. 3 exhibit a rapid drop which increases in relative sense as the ligand dose decreases. Over this initial period, which we refer to as Phase A, (cf. Fig. 1), R(t) exhibits a steep drop that becomes deeper as the drug dose increases. (b) After the brief initial adjustment period, the graphs for large doses reveal linear first order kinetics over a  (2) R 0 mg/L 12 (4) 12 (1) 12 (1) J Pharmacokinet Pharmacodyn (2012) 39:429-451 433 period of time (Phase B) that shrinks as the drug dose decreases. At the lowest dose the linear period has vanished and the graph exhibits nonlinear kinetics. (c) For the larger doses, there is an upward shift of the linear phase that appears to be linearly related to the ligand dose; the slope of this linear phase appears to be dose-independent. (d) The point of inflection in the log(L) versus time curvethe middle of Phase C-which we observe in the graphs for L 0 = 100 and 300, moves to the right as the initial dose increases, but stays at the same level. This is clearly seen in Fig. 3 in which the baseline value R 0 and the value of K d and K m are also shown. (e) For the lower doses we see that the log(L) versus time curve eventually becomes linear again, with a slope that is markedly smaller than it was in the nonlinear Phase C that preceded it. This part of the graph corresponds to Phase D in Fig. 1.
Summarising, in the ligand graphs of Fig. 3 we see for the higher drug doses the different phases A-D that were pointed out in Fig. 1. In the following analysis we explain these features and quantify them in that, for instance, we present estimates for the upward shift referred to in (c) and the right-ward shift of the inflection point alluded to in (d).

Dynamics after a bolus administration
Here we assume that ligand is supplied through an intravenous bolus administration and that there is no infusion, i.e., k f = 0. Thus, we focus on the system In order to obtain a first impression of typical concentration versus time courses for the three compounds, we carry out a few simulations of the system (9). We then present a mathematical analysis in which we delineate and discuss the four phases A, B, C and D in the characteristic ligand versus time graph shown in Fig. 1.

Simulations
We use the same initial doses as in the case study, i.e., the initial ligand concentrations are L 0 = 30, 100, 300 and 900 mg/L. The parameter values are given in Table 3, which is the same as Table 1, except that (i) the parameters V T and Cl d are absent because the tissue compartment has been taken out and (ii) the elimination rate has been reduced to k e(L) = 0.0015 h -1 so that the different phases in the ligand versus time graphs can be distinguished more clearly.
In Fig. 5 we show three ligand versus time graphs. In the figure on the left, L is given on a logarithmic scale and in the two figures on the right L is given on a linear scale, with the one on the right a blow-up of the initial behaviour when L 0 = 900 mg/L. (a) A rapid initial adjustment (see the blow-up on the right). (b) A first linear phase with a slope which is independent of the dose, and which shifts upwards as the drug dose increases. (c) A transition phase which shifts to the right as the drug dose increases, but maintains its level. (d) A final linear terminal phase with a slope k z that is again independent of the drug dose. For the parameter values of Table 3 we find that k z & k e(RL) = 0.003.
Since in Figs. 3 and 4 the time was restricted to 500 h, the terminal phase (D) has only just begun for the lowest dose, and not even started for the higher doses.
In Fig. 6 we present the receptor dynamics: concentration versus time profiles for, respectively, R, RL and R tot when L 0 [ R 0 . Whereas in [17] it was assumed that k out = k e(RL) , and hence the receptor pool was constant (R tot : R 0 ), for the parameters in Table 3 we have k out [ k e(RL) , so that R tot is no longer constant.
We make the following observations: (i) The total amount of target R tot (t)-free and bound to ligand-increases from its initial value R 0 = 12 to a maximum value R * , and then drops off again towards its baseline value R 0 . A similar observation was made in [15,22].
We shall prove that Thus, if k e(RL) \ k out , as is the case with the parameter values of Table 3, then R * [ R 0 and the total target pool increases before it returns to the baseline value R 0 . Alternatively, if k e(RL) [ k out , then R * \ R 0 and we show that the total target pool first decreases before it returns to the baseline value R 0 .
(ii) As the drug dose increases, R(t) & 0 for an increasing time interval and the graphs of RL(t) and R tot (t) trace-for the same increasing time interval-a common curve C in the (t, R tot )-plane (cf. Fig. 9). This curve C is monotonically increasing and tends to the limit R * as t ! 1: If k e(RL) [ k out , we show that an analogous phenomenon occurs along a curve C; which still tends to R * , but is now decreasing.
In Fig. 7 we present graphs of R -R 0 , RL and R tot -R 0 on a logarithmic scale. We note two conspicuous features: (i) The three graphs exhibit a kink (a sharp angle) which shifts to the right (increasing time) as the drug dose increases.
(ii) R(t) -R 0 tends to zero as t ! 1 in a bi-exponential manner, whilst RL(t) and R tot (t) -R 0 converge to zero in a mono-exponential way.

Low dose graphs
We conclude these simulations with a comparison of highdose and low-dose graphs. We do this by adding simulations for initial ligand concentrations which are smaller than R 0 . Specifically, we add the values L 0 = 0.3, 1, 3 and 10 mg/L to the graphs shown in Figs. 5 and 6. Figure 8 shows simulations of the ligand, target and complex concentration-time courses after eight different intravenous bolus doses. The initial drop will be difficult to capture unless that is taken care of experimentally within the very first minutes or so when the second-order process occurs. We see that in the low-dose graphs (L 0 \ R 0 ) the signature shape of Fig. 1 is no longer present. Instead, as L 0 decreases, the ligand curves become increasingly bi-exponential, and condensed into an apparent monoexponential decline. Still the terminal slope after the highest and the lowest dose are the same. As the ligand doses decrease, the target profile becomes less affected in terms of intensity (depth) and duration below the baseline concentration. In fact, one can show that if R(0) = R 0 and RL(0) = 0, then The proof of this upper and lower bound is given in Appendix 3. It is immediately clear that when L 0 ? 0 the upper as well as the lower bound converges to R 0 . Therefore, we may conclude that for every t [ 0, When we replace R(t) by R 0 in the system (9) the resulting system is linear, involving only L and RL. This explains the bilinear character of the log(L) versus time graphs (see also [17]).

Mathematical analysis
We successively describe the dynamics in the four phases: A-D. Throughout we assume that the ligand has large affinity for the receptor, that the elimination rates are comparable, and that the bolus dose is not too small. Specifically we assume: where M is a constant which is not too large, i.e., eM ( 1: These three assumptions were inspired by the parameter values given in Table 3 and the initial values of L and R. We also mention a similar model used for the study of Interferon-b 1a in humans [23] fitted with comparable parameter values. Phase A Ligand, receptor and receptor-ligand complex quickly reach Plateau values (L; R; RL) (see the right graph of Fig. 5). Since L 0 [ R 0 by Assumption B, the supply of free receptor is quickly exhausted, so that these plateau values are approximately given by Note that this is confirmed by the initial portion of the graph of L(t) shown in Fig. 5: we see that L drops by  Table 3. The dashed lines indicate the baseline R 0 and K d , and the dotted line the level R * approximately R 0 over a time span of about 0.04 h. In Appendix 2 we give details of the dynamics in this initial phase and the approach to the plateau values. We show that it takes place over a time interval (0, T 1 ), where 1 When L 0 = 900 mg/L this yields a half-life t 1/2 of about 0.01 h, in agreement with what is shown in Fig. 5. For a detailed study of Phase A we refer to [17] and to Aston et al. [24].
Phase B Over the subsequent time span when L(t) ) K d , say over the interval T 1 \ t \ T 2 , the three compounds are in quasi-equilibrium. This means that R, RL and R tot = R ? RL are-approximately-related to L through the expressions: Since in this phase L(t) ) K d , they may actually be approximated by the simpler expressions With these equalities, the system (4) can be reduced to the simpler form Note that for R tot we have obtained a simple indirect response equation (see also [9]).
Phase C When L(t) = O(K d ), say over the interval T 2 \ t \ T 3 , the approximation (17) is no longer valid. Applying a scaling argument appropriate for this regime, we show in Appendix 5 that to good approximation so that in this phase the rapid binding assumption is approximately satisfied [18].
Phase D When L(t) ( K d , i.e., beyond T 3 , the ligand concentration is so small that the dynamics is linear again. The critical times T 1 , T 2 and T 3 provide a natural division of the dynamics in four phases: A, B, C and D, as was done in Fig. 1. In Phase A (0 \ t \ T 1 ), ligand, receptor and complex reach quasi-equilibrium, in Phase B (T 1 \ t \ T 2 ), the bulk of the ligand is eliminated from the system while most of the receptor is bound to ligand and in quasi-equilibrium. Phase C (T 2 \ t \ T 3 ) is a nonlinear transitional phase in which L exhibits a steep drop, and finally, in Phase D (T 3 \t\1) the three compounds converge linearly towards their baseline values.
Receptor graphs: Phase B In Phase B, which extends over the interval T 1 \ t \ T 2 , the system (18) holds. Since the second equation only involves R tot as a dependent variable, it can be solved explicitly to yield Because T 1 is small by the estimate (15), it follows from the system (4) and the initial conditions that R tot (T 1 ) & R 0 . Thus, to good approximation, we may put T 1 = 0 and R tot (0) = R 0 , and so simplify the above expression for R tot (t) to where we recall from Eq. (10) that R * = k in /k e(RL) . Plainly, if T 2 were infinite, then We denote the graph of the function R tot (t) by C: In Fig. 9 we see how for the different drug doses, the simulations of R tot (t) follow the graph C up till some time, when they suddenly depart from C: Remark We recall from the approximation (17) that R(t) & 0 for T 1 \ t \ T 2 and hence that RL(t) & R tot in this phase of the dynamics, as we see confirmed in Fig. 6.
We conclude with a bound of the target pool. It is evident from the receptor graphs in Figs. 6, 7, 8 and 9 that for the data of the case study, the total receptor concentration R tot is not constant. However, it remains bounded for all time and does not keep on growing. In fact it is possible to prove the boundedness of R tot under very general conditions on the data, and actually obtain a sharp value for the upper bound. Specifically, we have The proof is given in Appendix 3.
Ligand graphs: Phase B Over the interval (T 1 ,T 2 ) in which it is assumed that L ) K d , we can use (17) to simplify the first equation from the system (18) as follows: Subtracting the second equation of (18) from (23) we obtain Equation (24), together with the initial value LðT 1 Þ ¼ L; can be solved explicitly. In Appendix 4 we show that the solution is given by where T 1 , being small, has been put equal to zero. We recall from Eq. (14) that L % L 0 À R 0 : In Fig. 10 we compare the numerically computed ligand versus time graphs L(t) with the analytic approximation L approx (t) given by the expression (25). It is evident that the two curves are very close until L has become so small that it is comparable to K d , i.e., until the end of Phase B, where the nonlinearity pitches in.
Since the function L approx (t) is decreasing, L approx (0) [ 0 and L approx (t) ?k in /k e(L) \ 0 as t ! 1; it follows that there exists a unique time T * [ 0 at which L approx (t) vanishes. We readily conclude from the definition of L approx (t) that T * is given by Remark It is clear from Eq. (24) that for larger values of L 0 , we may estimate k e(L) from the slope of log(L) and in Fig. 10 we see that T * yields a good estimate for T 2 , the end of Phase B. Information about k e(L) and T * combined yields an estimate for k in from the expression which can be derived from (26).
When K d is very small, as is the case for the parameter values of Table 1, the approximation L approx (t) given by (25) is valid for most of the range of L and can therefore be used to obtain an estimate for the area under the ligand curve AULC. An elementary computation yields Details are given in Appendix 4. This expression for AUCL yields for the clearance It is interesting to compare these approximate expressions for AUCL and CL with the corresponding expressions for mono-exponential ligand elimination. They are, respectively, D/(k e(L) V c ) and k e(L) V c . Thus we see that they are related by a factor which only depends on two dimensionless critical numbers: (i) the ratio l of R 0 and L 0 and (ii) the ratio j of the direct elimination rates of receptor, k out , and ligand, k e(L) . Ligand and receptor graphs: Phases C and D In Phase C the ligand concentration L has become comparable to K d and, as is shown in Appendix 5, we have to good approximation This suggests using a different scaling of L. For definiteness we assume that Phase C comprises the time interval in which L drops from 10 9 K d to 0.1 9 K d , and denote the times that L(t) reaches 10 9 K d and 0.1 9 K d by, respectively, T 2 and T 3 . Thus, L(T 2 ) = 10 9 K d and L(T 3 ) = 0.1 9 K d . In Phase C we use the approximation (30) in the ligand conservation law in (4): We now introduce K d , which-like L, R and RL, has the dimension of a concentration-as a reference variable for L, and introduce the dimensionless variables Using these variables in Eq. (31) we obtain Since e ( 1; and v = O(1) (cf. Fig. 6), we may neglect the term eu in the left-and the right-hand side of this equation and so obtain In the simulations shown in Figs. 5 and 6 we see that in Phase C, L drops rapidly from O(10 9 K d ) to O(0.1 9 K d ), i.e., by a factor 100, whilst R tot stays relatively close to R * and changes by no more than a factor 1/7 & 0.15. This suggests making the following assumption: Thanks to this assumption we may view v as a constant, which we may divide out and thus eliminate from the equation. We end up with a simple nonlinear equation for u, which is valid in Phase C: Since LðT 2 Þ ¼ 10 Á K d ; it follows that u(T 2 ) = 10. Equation (34) can be solved explicitly and we find for its solution: uðtÞ ¼ Ae Àk eðRLÞ ðtÀT 2 Þ 1 À Ae Àk eðRLÞ ðtÀT 2 Þ for t ! T 2 ; Returning to the original variables we obtain for the large time behaviour logfLðtÞg $ logðK d Þ þ logðAÞ À k eðRLÞ ðt À T 2 Þ as t ! 1 The asymptotic expression (36) yields estimates for (i) the terminal slope k z TMDD (k e(RL) ); (ii) the intercept of the asymptote of log{L(t)} of the ligand graph in the terminal Phase D with the vertical line {t = T 2 }.
The approximate identities in (30) imply that in Phases C and D, when L = O(K d ), we have to good approximation so that This is consistent with the value of k z found for L(t).
We also find that to good approximation so that This confirms what we see in Fig. 7: that for t [ T 2 the receptor concentration R(t) tends to R 0 in a bi-exponential manner, in contrast to the way RL(t) tends to zero, which is mono-exponential. For completeness we also compute the terminal slope by means of a standard analysis of the full TMDD system. This is done in Appendix 7. It is found that for the parameter values in Table 3, the terminal slope k z of all the compounds, is givento good approximation-by k z TMDD = k e(RL) . This confirms the limit in (36) and the exponent in (38).

Comparison with Michaelis-Menten kinetics
In many studies involving TMDD, models are employed that combine linear and saturable Michaelis-Menten type elimination (e.g. see [14]) of the form in which k, V max and K M are empirical parameters. The underlying assumption is that the MM-term can replace the combined first and second order processes of buildup and elimination via the complex in the TMDD model within a certain ligand concentration range. In light of (24), fitting the data for large values of the ligand concentration would yield k = k e(L) and V max = k in . Putting K M = K d , Eq. (41) then becomes Alternatively, we can take (31) as point of departure.
Following [15] we assume that R tot K d ( (K d ? L) 2 . Then the left hand side of (31) reduces to dL/dt. Assuming that R tot & R * and remembering that k e(RL) R * = k in , we may replace the factor k (e(RL) R tot in the right hand side of (31) by k in and so arrive at (42). Fitting to data of low ligand concentrations (L ( K d ), Eq. (42) reduces to the linear equation which yields a terminal slope k z MM given by This terminal slope is quite different from the value k z TMDD obtained in (36) and (38). For the parameter values of Table 3, we find that k z MM ) k z TMDD . Thus, the TMDD-model and the Michaelis-Menten (MM)-model exhibit very different terminal slopes, unless one also includes a non-specific peripheral volume distribution term in the MM-model. Michaelis-Menten model with peripheral compartment. Adding a peripheral compartment to the MM-model makes it possible to capture the slow terminal elimination that is typically seen in TMDD data. Figure 11 shows the regression of a 2-compartment model with parallel linear (Cl (L) ) and Michaelis-Menten (Cl = V max /(L p ? K M )) elimination: Fitting this model to the data shown in Fig. 3 results in the parameter estimates which, together with their precision (CV%), are given in Table 4. The reduced model mimics the concentration-time data for the two highest doses reasonably well, whereas the two lower doses display systematic deviations between observed and predicted data.
Since the reduced model has two parallel elimination pathways (linear and nonlinear) it has the intrinsic capacity of exhibiting linear first-order kinetics at low and at high concentrations. In the concentration-range in between it behaves nonlinearly. For higher concentrations the MM- Table 4 Parameter estimated from fitting the Michaelis-Menten model (Eq. (45)) to the data shown in Fig. 11. As in the Case Study,  Table 4. The dashed line in the middle of the plot indicates the estimated value of K M . Notice how far away it is from the original value of K m -marked by the thin drawn line-which was estimated by the TMDD model route is saturated and the linear elimination pathway dominates so that the system behaves linearly. However, the typical concentration-time pattern for ligand seen in a true TMDD system (cf. Figs. 5, 8, 10), cannot be fully described by the parallel linear-and MMelimination model. The reduced model displays typical biexponential decline (which is expected from a two-compartment model) at lower concentrations. That is generally not the case with the full TMDD model.
A clear distinction of the two models occurs initially, immediately after dosing (Phase A), when the second-order reaction between ligand and circulating target forms the complex. This process cannot be captured by the reduced model, which may cause biased estimates (too large) of the central volume.
In Table 5 we summarise these results. It shows that the MM-models (41) and (45) may be fitted successfully to the first part of the ligand versus time graph, although they miss the initial drop in Phase A. They catch the first part of Phase C, but the first model fails to catch the second part, where the graph joins up with the terminal Phase D. The second MM-model is an improvement, but still shows significant deviations for lower ligand-concentrations.
In this table, a plus (?) means that the corresponding phase can be adequately explained, whilst a minus (-) means that it cannot.

Constant rate drug infusion
We assume that the drug is administered through a constantrate infusion over a finite period of time, and we are interested in elucidating which parameters are critical in determining the time to steady state, the extent of steady state ligand, target and ligand-target complex concentrations and the dynamics after washout. Assuming that the infusion rate reaches its constant value k f in a negligible amount of time and that washout at time t washout is also instantaneous, we consider the following variant of the system (1) : in which H(t) denotes the Heaviside function: We assume that initially there is no ligand in the system, i.e., L 0 = 0. When the infusion lasts long enough, i.e., when t washout is large enough, the concentrations will converge towards their steady state values L ss , R ss and RL ss . Then, at washout, they will return to their pre-infusion values: L = 0, R = R 0 and RL = 0.
We first derive expressions for the steady state values. Then we carry out a series of simulations subject to the same assumptions as those made in (13), except that we replace Assumption C by As we shall see, L ss is an increasing function of k f , so that we require here that k f does not drop below a threshold value for which L ss (k f ) = R 0 . We discuss features of the dynamics exhibited in these simulations, especially the time to steady state after onset of infusion, and after washout.

Steady state concentrations of L, R and RL
For the steady state concentrations L ss , R ss and RL ss of the system (46) we find the following expressions. For the ligand-receptor concentration we obtain in which q = (k e(L) /k e(RL) )k out K m . In light of the conservation laws for ligand and target, we then obtain for the ligand and for the target, The expressions (47)-(49) are derived in Appendix 6. The formula (48) for the ligand shows that L ss will be smaller than expected from the ratio of ligand infusion rateto-clearance (In L /Cl (L) = k f /k e(L) ), due to the removal of ligand as part of the complex RL ss . The same reasoning may be used to explain why the circulating target concentration R ss given by (49) is smaller than the baseline concentration R 0 = k in /k out . Due to the removal of target by means of the complex, the target concentration R ss will drop further as the infusion rate increases and RL ss increases accordingly (cf. Eq. (49)).  Figure 12 shows graphs of L ss , R ss and RL ss as functions of the infusion rate k f for the parameter values of Table 3. Note that at large infusion rates, L ss increases approximately linearly, except for a downward shift. Indeed, this is confirmed analytically: by expanding the expression (48) for large values of k f , we obtain The reason is that for large infusions, elimination of ligand occurs primarily via the extra-target elimination route (k e(L) ). For small infusion rates the removal of ligand is also seen to be first order, but clearly at a much lower rate. In Fig. 12 the receptor concentration R ss is seen to decrease and the complex concentration RL ss is seen to increase as k f increases. However, it is interesting to observe that RL ss reaches an upper bound in spite of increasing levels of k f . This is also confirmed analytically: letting k f tend to infinity in (47) and (49) we obtain the limits R ss ðk f Þ ! 0 and RL ss ðk f Þ ! R Ã as k f ! 1 ð51Þ Thus, the upper bound is found to be R * = k in /k e(RL) . It is interesting to note that when plotted on a logarithmic scale, the elimination of target mirrors the elimination of ligand for larger values of the infusion rate (cf. Fig. 12 on the right). This can be understood from the relation obtained from (46). When we divide Eq. (52) by k on and take the logarithm, we obtain logðL ss Þ þ logðR ss Þ ¼ logðK m Þ þ logðRL ss Þ For larger values of k f we have RL ss & R * (cf. (51)), and hence logðL ss Þ þ logðR ss Þ % logðK m Þ þ logðR Ã Þ ð 53Þ which establishes the symmetry which is evident in the graphs for L ss and R ss shown on the right in Fig. 12.
Using the parameter values in Table 3, we obtain 1/ k e(L) = 667, k in /k e(L) = 73, R * = k in /k e(RL) = 36.7 and log(K m ) ? log(R * ) = 0.45. We see that these values are confirmed by the numerically obtained graphs shown in Fig. 12.
We note that the expressions (47)-(49) show that the steady state concentrations do not depend on the on-and off rates k on and k off individually, but only as part of the constant K m .

Simulations
We show simulations of concentration versus time graphs of the system (46) when drug is supplied through a constantrate infusion over a period of 5000 h. Four infusion rates are considered: k f = 0.12, 0.18, 0.30 and 0.54 (mg/L)/h. In Fig. 13 we show ligand graphs, on a linear and on a semilogarithmic scale, and in Figs. 14 and 15 we show graphs of the concentration of R, RL and R tot , first on a linear scale and then also on a logarithmic scale. In each of these figures we include the build-up phase as well as the washout phase.
In Fig. 13 we see that for larger values of k f the time-tosteady state of the ligand L is more or less independent of the infusion rate. The amounts of ligand are now so large that the receptor is quickly saturated and ceases to play an important role in the dynamics. What remains is the linear clearance of ligand, so that the ligand dynamics is described to good approximation by the equation which can be solved explicitly. Evidently which is consistent with (50) and  Fig. 12 The steady state concentrations L ss , R ss and RL ss graphed versus the infusion rate k f , on a linear scale (left) and on a log-log scale (right) for parameter values taken from Table 3 This results in a time to steady-state of 3-4 9 t 1/2 is &1850 h for the parameter values given in Table 3. The washout dynamics is very similar to the dynamics after a bolus dose, as described before: Phase A (0,T 1 ) now covers the infusion period so that T 1 coincides here with the time of washout. Phases B-D are plainly evident in the post-washout dynamics.
The dynamics of the receptor R and the receptor-ligand complex RL are shown in Figs. 14 and 15. As in Phase A in the bolus administration, the pre-dose receptor pool (R 0 ) quickly binds to the ligand. We see that the speed of receptor depletion increases with increasing infusion rate k f , consistent with the half-life estimate (15) after a bolus administration.
In due course, additional receptor is formed, albeit more slowly-and binds immediately to the ligand-resulting in an increase of RL as was observed after a bolus dose (see also Fig. 6). The dynamics is very similar-compare Eq.    Table 3. Note that the convergence of target R to R 0 is bi-exponential and that the decline of complex RL to zero, and the convergence of total target R tot to R 0 are mono-exponential J Pharmacokinet Pharmacodyn (2012) 39:429-451 443 This results in a time-to-steady-state of 3-4 9 t 1/2 is &924 h (cf. Table 3), which we see confirmed in Figs. 14 and 15.
Eventually RL(t) levels off at the steady state value RL ss , which is close to R * for the larger infusion rates (cf. Eq. (51)). After washout, when k f is large, R(t) & 0 and RL(t) & R * for a while before they abruptly return to their baseline values. It is evident from Fig. 15 that, as in Fig. 7, initially the slope of log(R 0 -R) is steeper than that of log(RL). This is in agreement with the analysis presented in Appendix 5, where it is shown that over that period the half-lives of R(t) -R 0 and RL(t) are, respectively, O(1/k out ) and O(1/k e(RL ).

Discussion and conclusion
We have shown how the concentration profile of ligand, receptor and ligand-receptor complex in the TMDD model can be divided into four different phases and how for each of these phases closed-form approximations can be derived. Inspired by a specific case study, the following assumptions were made about the parameter values (see also (13) and (13*)): for some moderate constant M [ 0, and L 0 [ R 0 . When ligand is administered through a bolus dose, L 0 [ R 0 , and the conditions in (58) are satisfied, four phases can be distinguished in the ligand elimination graph: a brief initial Phase A, a slow linear Phase B, a rapid nonlinear Phase C and then again a slow linear terminal Phase D (cf. Fig. 1). Thanks to accurate analytical approximations for these four phases as shown in the Eqs. (14) for Phase A, In Table 6 we list the parameter values which play a central role in the different phases of TMDD graphs. In Phase A : the drop (R 0 ) and the duration (O(1/(k on L 0 ))), in Phase B : the slope (k e(L) ) and the receptor input (k in ), in Phase C : the depth (K d ) and in Phase D : the terminal slope (k e(RL) ). In brackets we have included the parameters which can be estimated when results from earlier phases are used. Thus, since Phase A yields an estimate for R 0 = k in /k out and Phase B an estimate for k in , an estimate for k out follows.
It should be noted though that estimating R 0 may be difficult, since often there are no data for the first phase because it is over very quickly.
In Fig. 16 we summarise these results and show in the schematic ligand versus time graph (recall Fig. 1) the parameters which may be estimated from the different phases.
The four phases identified in the ligand elimination graph after a bolus administration are reflected in the structure of the receptor versus time graphs (receptor, receptor-ligand complex, and total amount of receptor). During Phase A the receptor pool is quickly depleted, and it remains so during Phase B. Then, during the Phases C and D it climbs back to the terminal baseline level R 0 .
The analytic approximations obtained for L, R and R tot may be used to verify the validity of the assumptions which underpin different approximations to the full TMDD model when the Assumptions A, B and C (or C * ) regarding the parameters are satisfied.
1. The rapid binding model [6,18], in which it is assumed that where K d is defined in (7). In Phase B, which is characterised by R(t) & 0, we have dR/dt & 0 so that, according to the second equation of the system (9), we have approximately which disagrees with (59). In contrast, in Phases C and D the identity (59) is satisfied according to the results established in Appendix 5 (cf. (102)), and Appendix 7 where it was shown that k z = k e(RL) . In Fig. 17 we show how the quantity D ¼ L Á R À K d RL varies with time and how D rapidly jumps from k in /k on down to zero at the transition of Phases B and C.
2. The quasi-steady-state model [7,18] in which it is assumed that where K m is defined in (8). Evidently, this assumption is not valid during Phases C and D, but it is during that part of Phase B in which RL(t) & R * . In that interval dRL/dt & 0 (cf. Fig. 6) and hence, by the third equation of the system (1), condition (61) is approximately satisfied.
Anchored on the data of the Case Study, the analysis in this paper is based on the Assumptions A, B and C (or C * ). The question arises whether the characteristic features of the ligand elimination curves, such as shown in Fig. 1, are still present when these assumptions are not met.
In general, the behaviour of nonlinear systems such as (1) is very sensitive to the values of the parameters and initial data involved. However, a number of features of the ligand versus time graphs is quite robust in that they may survive if e.g. Assumption A is not satisfied and K d and R 0 are comparable. Thus, the estimate (15) for T 1 suggests that the initial Phase A will remain short relative to typical times over which the other processes develop when K d /L 0 is small. We refer to [17] and [18] for a detailed analysis of this situation.
The approximate expressions for L and R tot in Phase B (Eqs. (25) and (20)) are still valid provided that R(t) & 0.
This will still be the case when Assumption A is replaced by K d ( L 0 [17,18].
In contrast, the analysis of the dynamics in Phase C that is carried out in Appendix 5 depends critically on Assumption A. It will be interesting to study the dynamics beyond Phase B when Assumption A does not hold, as it will be interesting to see how the value of a, b and c affects the dynamics.
We have selected a set of data (ligand and circulating target and complex) with low experimental variability, concentration-time courses at four ligand doses given as bolus injections, and well-spaced data in time that captures the necessary phases and shapes of a typical TMDD system. Based on this approach and the mathematical/analytical analysis, we can draw conclusions about the identifiability of the model parameters and appropriate system. When data are less precise and information rich, or, when target and/or complex are less accessible, the a priori expectations of parameter accuracy and precision will be lower.
Acknowledgments The constructive criticisms of the manuscript of the reviewers improved the quality of the text, and are highly appreciated.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.  Dimensionless quantities x, y, z L/L 0 , R/R 0 , RL/R 0 ; dimensionless concentrations x; y; z L=L 0 ; R=R 0 ; RL=R 0 ; dimensionless concentrations u, v, w L/K d , R/R 0 , RL/R 0 ; dimensionless concentrations s k on R 0 t; dimensionless time e The symbol A is defined by the expression B f(x) * g(x) as x ? ? f(x)/g(x)?1 as x? ?
There exist constants K [ 0 and

Appendix 2: Short time analysis
We present a mathematical analysis for Phase A and show how ligand, receptor and ligand-receptor complex rapidly converge to the new Plateau values L; R and RL given in (14), where we assume that L 0 [ R 0 . In order to identify terms in the system (9) which dominate the dynamics in this initial phase, we define the dimensionless concentrations and the dimensionless parameter l = R 0 /L 0 . Substituting them into the system (1) we obtain Introducing the dimensionless time s = k on R 0 t, then results in the dimensionless system where we have introduced the following dimensionless constants (see also (58)): We multiply the third equation of (64) by l and add the resulting equation to the first equation to obtain Similarly, we add the second and the third equation to obtain d ds ðy þ zÞ ¼ e bð1 À yÞ À cz f g ð 67Þ Since e ( 1 and ea ( 1; eb ( 1 and ec ( 1 by Assumption B (cf. (13)), it follows that to good approximation, xðsÞ þ lzðsÞ ¼ 1 and yðsÞ þ zðsÞ where we have used the fact that x(0) = 1, y(0) = 1 and z(0) = 0. Note that in light of the definition of s and e; the assumption es ( 1 is equivalent to k off t ( 1. Remark These approximate equalities state that the total amount of ligand (L ? RL) and the total amount of receptor (R ? RL) in the system remain more or less constant during this initial phase. Using the expressions from (68) to eliminate x and y from the equation for z in (64), we obtain where we have dropped the term elð1 þ cÞz because elc\ec ( 1 so that its impact is negligible. The right hand side of Eq. (69) has two zeros: z = 1 and z = 1/l [ 1. Plainly f(z) [ 0 for 0 B z \ 1 and for z [ 1/l, whereas f(z) \ 0 for 1 \ z \ 1/l. Since z(0) = 0 it follows that z(s) increases and approaches the lower zero z ¼ 1 as time tends to infinity: When we linearise Eq. (69) at z = 1, write zðsÞ ¼ 1 þ fðsÞ and omit the terms involving f n for n [ 1 we obtain the equation Therefore, Thus, for es ( 1; i.e., for k off t 1/2 to be small, as required in (68), the parameters and the initial ligand concentration need to be such that k off /(k on (L 0 -R 0 ) ( 1 or L 0 -R 0 ) K d .
In light of (68), the limit for z(s) obtained in (70) implies the following limits for x(s) and y(s): Returning to the original variables and writing L ¼ L 0 x; R ¼ R 0 y and RL ¼ R 0 z; we conclude that over a time span of O(1/(k on (L 0 -R 0 )) we have Appendix 3: Bounds of the receptor concentration We establish two a-priori bounds: one for the total receptor concentration by R tot = R ? RL and one for the free receptor concentration R. Both bounds are global in time, i.e., they hold for all time.
Lemma 1 Suppose that R tot ð0Þ ¼ R 0 ¼ def k in =k out : Then we have: Proof Adding the last two equations of the system (1), we obtain We now consider two cases: (I) k e(RL) C k out and (II) k e(RL) \ k out Case I (k e(RL) C k out ). We write Eq. (76) as and conclude that since, by assumption, k e(RL) C k out , we have Define the function q 1 ðtÞ ¼ def R tot ðtÞ À R 0 : Then, transforming Eq. (78) and the initial value of R tot to this new variable, we obtain dq 1 dt À k out q 1 ; q 1 ð0Þ ¼ 0 ð79Þ This implies that d dt e k out t q 1 ðtÞ À Á 0 and hence q 1 ðtÞ q 1 ð0Þ ¼ 0 Therefore, as asserted. This completes the proof of the first inequality. h Case II (k e(RL) \ k out ). We now write Eq. (76) as dR tot dt ¼ k in À k eðRLÞ R tot þ ðk eðRLÞ À k out ÞRL \k in À k eðRLÞ R tot Put q 2 ðtÞ ¼ def R tot ðtÞ À R Ã ; where R * = k in /k e(RL) , and note that by assumption, q 2 (0) = R tot (0) -R * \ 0. Proceeding as in the previous case we conclude that q 2 (t) \ 0 for t C 0 and hence that This completes the proof of the second inequality, and hence of the lemma. h

Remarks
(1) The upper bound in Lemma 1 is independent of the amount of ligand that is supplied to the system and it is valid whether drug is given through a bolus administration or by means of a constant rate infusion. (2) Lemma 1 generalises an observation made before in [3,15,17], which states that R tot ðtÞ R 0 if and only if k eðRLÞ ¼ k out ð83Þ Next, we prove bounds for R.

Lemma
2 Suppose that L(0) = L 0 , R(0) = R 0 and RL(0) = 0. Then Proof From the conservation law (4) we deduce that L tot (t) is a decreasing function of t. Therefore, LðtÞ þ RLðtÞ L 0 for all t ! 0 ¼) RLðtÞ L 0 ð85Þ In order to prove the upper bound, we write where Plainly R(0) = R 0 \ R ? , and we claim that R(t) \ R ? for all t C 0. Suppose to the contrary that there exists first a time t 0 [ 0 such that R(t 0 ) = R ? , i.e., RðtÞ\R þ for 0\t\t 0 and Rðt 0 Þ ¼ R þ We see that R(t) approached R ? from below as t % t 0 ; so that dR/dt C 0 at t = t 0 . However, we deduce from (86) that dR/dt(t 0 ) \ 0. Therefore, we have a contradiction and we may conclude that there exists no time t 0 [ 0 for which R(t) attains the value R ? . This proves the claim and thereby the upper bound. h To prove the lower bound we write where Proceeding as with the upper bound, we establish the lower bound.
Corollary 1 Since R ? (L 0 ) ? R 0 and R -(L 0 ) ? R 0 as L 0 ? 0, it follows that for any time t C 0, Appendix 4: Derivation of the approximation L approx (t) and expressions for AULC and CL We wish to solve the initial value problem dL dt ¼ Àk eðLÞ L À k in ; Lð0Þ ¼ L ð92Þ where L % L 0 À R 0 is the ligand concentration right after Phase A. We multiply this equation by e k eðLÞ t to obtain d dt ðe k eðLÞ t LÞ ¼ e k eðLÞ t k in Integration over (0,t) yields the desired expression for the solution L approx (t) of Eq. (92): L approx ðtÞ ¼ ðL þ mÞe Àk eðLÞ t À m; m ¼ k in k eðLÞ ð93Þ Plainly, L(t) ? -m \ 0 as t ! 1: Since L approx (t) is decreasing and L(0) [ 0, it follows that there exists a unique time T * [ 0 for which L approx (t) vanishes. An elementary computation yields the following expression for T * : For the three largest initial ligand concentrations used in case study this formula yields, as L 0 decreases, respectively, T In order to derive the approximation for AULC(L 0 ), we write Z T Ã 0 L approx ðtÞdt ¼ 1 k eðLÞ Z k eðLÞ T Ã 0 fðL þ mÞe Às À mgds ¼ 1 k eðLÞ fðL þ mÞð1 À e Àk eðLÞ T Ã Þ À mk eðLÞ T Ã g When we now substitute the expression for T * from (94) into (95), we obtain the desired formula: from which (28) follows. Applying the definition CL(D) = D/AULC(D) we obtain an expression for CL(D). dn dt ¼ Àðk on R 0 þ k eðLÞ Þn þ k off f dg dt ¼ Àk on R 0 n À k out g þ k off f df dt ¼ k on R 0 n À ðk off þ k eðRLÞ Þf 8 > > > > < > > > > : ð111Þ when higher order terms are omitted. It is convenient to write this system in vector-and matrix notation: and A is the coefficient matrix of the linear system (111): We assume that the matrix -A has three distinct eigenvalues k 1 , k 2 and k 3 and that Y 1 , Y 2 and Y 3 are the corresponding eigenvectors. The General Solution of Eq.
(112) then takes the form where C 1 , C 2 and C 3 are arbitrary constants. The eigenvalues k i (i = 1, 2, 3) are the roots of the equation detðA þ kIÞ ¼ 0 Thus, we find that k 1 ¼ k out and k 2 and k 3 are the roots of the quadratic equation where a ¼ k on R 0 þ k eðLÞ þ k off þ k eðRLÞ ¼ k on R 0 1 þ ea þ eð1 þ cÞ f g b ¼ k on R 0 k eðRLÞ þ k eðLÞ k off þ k eðLÞ k eðRLÞ  (13)), we have to good approximation, This implies that b=a 2 ¼ ec ( 1; so that we conclude that k 1 ¼ k out ¼ 0:0089; k 2 ¼ a ¼ k on R 0 ¼ 0:11 and i.e., the terminal slope k z is given by k z = k e(RL) = 0.003. This agrees with the numerical value found in Figs. 5, 6, 7 and 8. For an earlier derivation we refer to [25].
Remark Naturally, when k e(RL) = 0 then ligand will still leak out of the system, albeit more slowly, thanks to direct elimination of ligand. For this case we obtain the following terminal slope: