Why a Terminal Run Is Informative
A sequence of below-quantification-limit (BQL) records late in a concentration-time profile is not a sequence of measured low concentrations. Each record says that the assay result fell below its applicable lower limit of quantification (LLOQ). The event contains information about the model-predicted concentration and residual variability, but it does not identify the exact concentration at that time.
The order and timing of the records matter. A run beginning shortly after a quantified sample can challenge a slowly declining profile; a run beginning later may be compatible with that same profile. The appropriate likelihood uses the probability of each recorded BQL event under the observation model. Some M-methods delete or replace BQL records according to a rule, which can change the information presented to the model; M3 instead uses the censored likelihood.
Likelihood Contribution of a Run
Let \(Y_{ij}\) be the assay result for subject \(i\) at time \(t_{ij}\), and let \(L_{ij}\) be the applicable LLOQ. Let \(\eta_i\) denote the subject-specific random effects and let \(f(t_{ij};\theta,\eta_i)\) be the model-predicted concentration, where \(\theta\) denotes population parameters. Under a continuous conditional observation distribution with conditional CDF \(F_{ij}(\,\cdot\mid\eta_i)\), a BQL record contributes
Continuity at the LLOQ is important to this equality. In general, for a distribution that might have a point mass at \(L_{ij}\), the correct expression is \(P(Y_{ij}<L_{ij}\mid\eta_i)=F_{ij}(L_{ij}^{-}\mid\eta_i)\), where \(F(L^-)\) is the left limit of the CDF. The normal observation model used in the examples below is continuous, so the CDF at the LLOQ is appropriate.
Suppose a run contains observations \(j=a,\ldots,b\). If observations are conditionally independent given the subject’s random effects and model parameters, its conditional likelihood contribution is
This product is conditional on the shared \(\eta_i\). The subject-level marginal contribution is obtained by integrating over the random-effects distribution:
where \(g(\eta_i\mid\theta)\) is the random-effects density. In general, this is not the product of the separate marginal BQL probabilities. The observations share the same subject-specific effects, so integrating each probability separately would discard their dependence induced by \(\eta_i\).
For example, under an additive normal residual model, \(Y_{ij}=f(t_{ij};\theta,\eta_i)+\epsilon_{ij}\), with \(\epsilon_{ij}\sim N(0,\sigma^2)\), the conditional BQL probability is \(\Phi((L_{ij}-f(t_{ij};\theta,\eta_i))/\sigma)\). Here \(\Phi\) is the standard normal CDF. Conditional independence of residuals, the specified residual distribution, and the LLOQ attached to each record are assumptions of this calculation; they should reflect the analysis model and assay data.
Worked Example: Likelihood of Three Consecutive BQL Records
Consider hypothetical terminal observations with a common LLOQ of \(1.00\) mg/L and additive normal residual standard deviation \(\sigma=0.300\) mg/L. For one subject at fixed random effects, suppose the model-predicted concentrations are \(0.8131\), \(0.6024\), and \(0.3306\) mg/L at three consecutive times. Assume the conditional residuals are independent and continuous normal.
The BQL probabilities are calculated from \(P(Y<L\mid\eta)=\Phi((L-f)/\sigma)\). For the first record, the standardized value is \((1.00-0.8131)/0.300=0.6230\), giving \(P_1=\Phi(0.6230)\approx0.7333\). For the second, \((1.00-0.6024)/0.300=1.325333\ldots\), giving \(P_2\approx0.9075\). For the third, \((1.00-0.3306)/0.300=2.231333\ldots\), giving \(P_3\approx0.9872\).
By conditional independence, the probability of observing the full three-record run at these fixed effects is
Thus, under these specified predictions and error assumptions, the run has a conditional likelihood contribution of about \(0.657\). It does not imply that any of the three concentrations equalled a particular fraction of the LLOQ. With random effects estimated from a population model, the subject-level contribution must also account for their distribution, as shown by the integral above.
How Timing Constrains a Declining Profile
A terminal-phase model often represents expected concentrations as declining over time. For illustration, take \(f(t)=A\exp[-\lambda(t-t_0)]\) for \(t\ge t_0\), with \(A>0\) and \(\lambda>0\). A larger \(\lambda\) means a faster exponential decline. The model does not make each assay result equal to this curve: residual variation means that a result can be BQL even when the prediction exceeds the LLOQ, or quantified when it is below it.
A short justification explains why a run can distinguish profiles. Under the continuous additive normal model, a BQL event at time \(t\) has probability \(\Phi((L-f(t))/\sigma)\). For fixed \(L\) and \(\sigma>0\), this probability increases as \(f(t)\) decreases, because the argument \((L-f(t))/\sigma\) increases and \(\Phi\) is increasing. Consequently, a series of BQL events can favor a profile that has declined sufficiently by those times over one that predicts substantially higher concentrations. This is probabilistic evidence, not a direct measurement of the decline rate.
Worked Example: Comparing Two Decline Rates Against a Run
Suppose three hypothetical samples at days 2, 3, and 4 are all BQL, with \(L=1.00\) mg/L and additive normal \(\sigma=0.300\) mg/L. Compare two fixed terminal profiles, both with \(A=4.00\) mg/L at day 1: a slower profile with \(\lambda=0.5\) per day and a faster profile with \(\lambda=1.0\) per day. Assume conditionally independent residuals and use these profiles only to illustrate the run contribution.
For the slower profile, predicted concentrations at days 2, 3, and 4 are \(4e^{-0.5}=2.426\), \(4e^{-1}=1.472\), and \(4e^{-1.5}=0.8925\) mg/L. The corresponding standardized LLOQ distances are \((1-2.426)/0.300=-4.754\), \((1-1.472)/0.300=-1.572\), and \((1-0.8925)/0.300=0.3583\). The BQL probabilities are approximately \(0.0000010\), \(0.0580\), and \(0.640\), so the run contribution is about \(3.7\times10^{-8}\).
For the faster profile, the predictions are \(1.472\), \(0.5413\), and \(0.1991\) mg/L. The corresponding BQL probabilities are approximately \(\Phi(-1.572)=0.0580\), \(\Phi(1.529)=0.9368\), and \(\Phi(2.670)=0.9962\). Their product is about \(0.0541\). Given these fixed profiles and assumptions, the observed run is much more compatible with the faster decline. In a population analysis, this comparison would be part of the full likelihood, including quantified observations and integration over random effects; it is not by itself an estimate of \(\lambda\).
Likelihood and M-Method Handling
The censored-likelihood approach retains each BQL record as the event that it represents. In the terminology of the M-methods, M3 uses this censored likelihood. It uses the LLOQ and the observation model to calculate a probability, rather than assigning an invented concentration. In contrast, M1 deletes BQL records, while substitution methods such as M2 and M4 enter a chosen value—commonly half the LLOQ or the LLOQ—as though it were an observed concentration. These procedures therefore do not contribute the same information to a fit.
| Handling | What enters for each BQL record | Implication for a terminal run |
|---|---|---|
| Censored likelihood (M3) | Probability of the event \(Y<L\) under the observation model | Retains timing and event information without claiming an exact concentration |
| Deletion (M1) | No contribution from the BQL record | Removes information about when the profile fell below the assay limit |
| Half-LLOQ substitution (M2) | A value of \(L/2\) treated as a concentration observation | Creates repeated numeric values that may form an artificial plateau |
| LLOQ substitution (M4) | A value of \(L\) treated as a concentration observation | Places each BQL result at the assay boundary rather than below it |
Worked Example: How Substitution Can Create a False Plateau
A hypothetical profile has quantified concentrations of \(8.0\) mg/L on day 0 and \(4.0\) mg/L on day 1, followed by BQL records on days 2, 3, and 4. Let the LLOQ be \(1.0\) mg/L. Compare what the handling rules imply about the recorded values; this simple calculation is not a full population-model fit.
Under M2, each BQL record is replaced by \(L/2=0.5\) mg/L. The apparent log-concentration slope from day 1 to day 2 is \((\ln(0.5)-\ln(4.0))/(2-1)=-\ln(8)\approx-2.079\) per day. Between days 2 and 4, the assigned values are both \(0.5\) mg/L, so their apparent slope is \((\ln(0.5)-\ln(0.5))/(4-2)=0\) per day. Thus the substitutions create a flat late segment that the BQL records did not establish.
Under M4, the assigned value is \(1.0\) mg/L at all three times. The apparent slope from day 1 to day 2 is \((\ln(1)-\ln(4))/(2-1)=-\ln(4)\approx-1.386\) per day, followed by a zero slope from day 2 to day 4. Under M1, the three terminal records are removed, so they provide no information about the decline after day 1. The censored likelihood instead contributes the probability of being below \(1.0\) mg/L at each scheduled time; it neither imposes a flat segment nor discards the timing.
Assumptions, Diagnostics, and Limitations
A run likelihood is only as credible as its observation and design assumptions. Check that each record’s LLOQ is correct, that the BQL status is encoded as an event rather than an exact value, and that the selected residual distribution matches the model being fit. A product of conditional probabilities also assumes residual independence over the run given the model and random effects. If residuals are correlated, the joint probability of the sequence must reflect that correlation rather than being calculated as a simple product.
Inspect the fit in the context of the quantified data near the run, the sampling times, and the assay limits. A model that predicts almost no chance of a long run may indicate a poor decline profile, a residual-error specification problem, a mistaken LLOQ, or an unaccounted-for dependence; it does not automatically identify which component is wrong. Conversely, a high predicted run probability does not prove that the model is correct. The model should also be checked against other information in the data and plausible parameter values.
The likelihood uses the fact that an observation is below the threshold, but it cannot recover its unobserved concentration. A long run can inform the rate and timing of decline only through the assumed concentration and residual models. If most relevant observations are BQL, different parameter combinations may explain the same events. Substitution may appear simpler, but it changes the data as if concentrations had been observed at chosen values; deletion avoids those invented values but discards the events altogether. Neither alternative preserves the event information in the way the censored likelihood does.
Common Mistake
Confirm which consecutive records are BQL, their sampling times, and the LLOQ applicable to each record.
State the conditional distribution, residual assumptions, and whether the distribution is continuous at the LLOQ.
Calculate the conditional probability of each \(Y<L\) event and combine them using the correct joint model.
Assess how the run supports or challenges predicted timing of decline, while retaining uncertainty about the unobserved concentrations.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- For a continuous conditional distribution, what probability does one BQL record contribute? What adjustment is needed if there may be a point mass at the LLOQ?
- Why is a run’s marginal likelihood generally not the product of each record’s marginal BQL probability?
- With \(L=2.0\) mg/L, \(\sigma=0.5\) mg/L, and a conditional predicted concentration of \(1.5\) mg/L under additive normal error, calculate the BQL probability.
- What information about a declining profile is lost by M1 deletion, and what unsupported assumption is introduced by M2 substitution?
- Does a run of BQL records prove that the concentration fell monotonically at every time in the run? Explain.
Solutions. (1) The contribution is \(P(Y<L\mid\eta)=F(L\mid\eta)\) when the conditional distribution is continuous at \(L\); generally use the left limit \(F(L^-\mid\eta)\). (2) The records share subject-specific random effects, so their conditional probabilities must be multiplied first and then integrated jointly over the random-effects distribution. (3) The standardized value is \((2.0-1.5)/0.5=1\), so the probability is \(\Phi(1)\approx0.8413\). (4) M1 removes the evidence about when results fell below the assay limit; M2 treats \(L/2\) as if it were the measured concentration. (5) No. BQL status is a threshold event affected by residual variability; it does not establish exact concentrations or monotonicity between sampling times.