Diagnosing Instability and Oscillation in Titration

Why a Titration Rule Can Oscillate

A titration algorithm can move in the correct direction after each measurement and still perform poorly. If its updates are too large, or if each update uses an old measurement, the dose may repeatedly cross the target instead of settling near it. The key is to analyze the feedback rule around its target: small deviations either shrink, persist, or grow from one update to the next.

The earlier material on titration rules written as state-update equations and proportional dose-adjustment rules provides the setting for this analysis. Here, the focus is the local behavior of a proportional rule when the response measurement is delayed. The analysis is deliberately linear: it describes small deviations near a target, not every feature of a full pharmacokinetic or pharmacodynamic model.

Definition: Let $D^*$ be a target dose and $e_k=D_k-D^*$ the dose error at update $k$. Suppose the response near the target is approximately $R_k-R^*=g e_k$, where $g>0$ is the local response slope. A proportional update with gain $K>0$, using a measurement that is $d$ update steps old, has the linearized error equation $$ e_{k+1}=e_k-ae_{k-d}, \qquad a=Kg. $$ Here $d$ is a nonnegative integer, and $a$ is the dimensionless effective gain. The equation assumes a fixed local slope, a correctly signed response adjustment, a constant gain, and a pure delay of $d$ updates. The starting values $e_{-d},\ldots,e_0$ specify the initial history.

The product $Kg$, rather than $K$ alone, governs local behavior. A dose-adjustment gain has units of dose per response unit, while $g$ has units of response per dose, so their product is dimensionless. A steep exposure-response relationship therefore makes the same dose-adjustment gain more aggressive. The delay $d$ represents feedback latency in update steps; it is not, by itself, a model of gradual pharmacokinetic equilibration or a changing response slope.

Characteristic Roots and Stability

To identify the behavior of the recurrence, try a component of the form $e_k=C\lambda^k$, where $C$ is a constant and $\lambda$ is a characteristic root. Substitution gives

C\lambda^{k+1}=C\lambda^k-aC\lambda^{k-d} \quad\Longrightarrow\quad \lambda^{d+1}-\lambda^d+a=0.

For $a>0$, the target is locally asymptotically stable if every characteristic root has magnitude less than one. In that case every initial error history produces errors that approach zero. A root outside the unit circle permits error components that grow; roots on the unit circle prevent convergence for some initial histories. Complex roots produce oscillatory components, while a negative real root produces a sign-alternating component.

Theorem: For the recurrence $e_{k+1}=e_k-ae_{k-d}$, with integer $d\geq 0$ and $a>0$, convergence to zero for every initial history occurs exactly when $$ 0<a<2\sin\!\left(\frac{\pi}{4d+2}\right). $$ At equality the system is not asymptotically stable. For $d\geq1$, a conjugate pair of characteristic roots lies on the unit circle. For $d=0$, the single characteristic root is $\lambda=-1$.

The bound shrinks as the delay grows. With no delay, gains below 2 are stable; with one update of delay, the bound is 1; with two updates of delay, it is about 0.618. Thus, a gain that is safe with immediate feedback may be unstable when measurements arrive later.

Why the Stability Bound Has This Form

The boundary of stability occurs when a characteristic root reaches the unit circle. Set $\lambda=e^{i\theta}$, for $0<\theta<2\pi$, in the characteristic equation. Rearranging gives $a=\lambda^d(1-\lambda)$. Since $a$ is positive and real, both the magnitude and phase must match.

|a|=|1-e^{i\theta}|=2\sin(\theta/2), \qquad d\theta+\frac{\theta}{2}-\frac{\pi}{2}=2\pi m, \quad m\in\mathbb{Z}.

The phase condition gives $\theta=(1+4m)\pi/(2d+1)$. The angles in $(0,2\pi)$ correspond to $m=0,\ldots,d$. Their half-angles run from $\pi/(4d+2)$ to $\pi-\pi/(4d+2)$; the smallest possible value of $2\sin(\theta/2)$ is therefore $2\sin(\pi/(4d+2))$. This is the first gain at which a root can reach the unit circle as $a$ increases.

For completeness, at $a=0$ the characteristic polynomial is $\lambda^d(\lambda-1)$. The root at 1 moves inside the unit circle for small positive $a$: implicit differentiation of $\lambda^d(\lambda-1)+a=0$ at $(\lambda,a)=(1,0)$ gives $d\lambda/da=-1$. The other roots start at zero and remain inside for sufficiently small $a$. Roots vary continuously with the polynomial coefficients, so stability is retained until a root first reaches the unit circle. The boundary calculation identifies that first gain. When $d=0$, the polynomial is $\lambda-1+a$, so $\lambda=1-a$: the boundary $a=2$ is the single root $-1$, not a conjugate pair.

Worked Example: No Delay at the Stability Boundary

Suppose a local titration rule has no measurement delay, so $d=0$, and its effective gain is $a=2$. Take $e_0=1$ mg as a hypothetical initial dose error. The recurrence is $e_{k+1}=e_k-2e_k=-e_k$, and its characteristic root is $\lambda=1-2=-1$.

The resulting errors are $e_0=1$ mg, $e_1=-1$ mg, $e_2=1$ mg, and $e_3=-1$ mg. Their magnitude does not decrease, so the trajectory oscillates persistently rather than converging. This is the boundary case: $a=2$ is not inside the stable range $0<a<2$. For comparison, at $a=1.4$, the root is $1-1.4=-0.4$, so the errors alternate in sign while shrinking by a factor of 0.4 per update.

Damped Oscillation Within the Stable Range

Stability does not mean that every error moves monotonically toward zero. A stable system can oscillate, provided the oscillations decay. Complex roots occur in conjugate pairs because the recurrence has real coefficients. A pair with magnitude $r<1$ contributes a sinusoidal pattern whose envelope shrinks approximately as $r^k$. A negative real root likewise creates alternating signs, with a shrinking envelope if its magnitude is below one.

Worked Example: One-Step Delay With Damped Oscillation

Consider a hypothetical local rule with one update of measurement delay, $d=1$, and effective gain $a=0.8$. The stability bound is $2\sin(\pi/6)=1$, so this gain is below the threshold. The characteristic equation and roots are

\lambda^2-\lambda+0.8=0, \qquad \lambda=\frac{1\pm\sqrt{1-3.2}}{2} =0.5\pm0.7416i.

Each root has magnitude $\sqrt{0.5^2+0.7416^2}\approx0.8944$, which is less than one, and the roots are complex. Thus the local error is a damped oscillation. For the hypothetical initial history $e_{-1}=1$ mg and $e_0=1$ mg, the first values are $e_1=1-0.8(1)=0.2$ mg, $e_2=0.2-0.8(1)=-0.6$ mg, $e_3=-0.6-0.8(0.2)=-0.76$ mg, and $e_4=-0.76-0.8(-0.6)=-0.28$ mg. The values do not decrease at every step, but the root magnitudes establish that the oscillatory components decay over time.

When Oscillation Becomes Divergence

At a gain above the stability threshold, at least one characteristic root has magnitude greater than one. If that root is complex, its contribution typically both oscillates and grows. The exact trajectory depends on the initial history: a particular initial condition can suppress a mode, but stability requires convergence for every history. In practice, a titration rule should not be called stable merely because one short trajectory looks well behaved.

Worked Example: One-Step Delay Above the Threshold

Now use the same one-step delay, $d=1$, but increase the hypothetical effective gain to $a=1.2$. Since the threshold is 1, this gain is outside the stable range. The characteristic roots are

\lambda^2-\lambda+1.2=0, \qquad \lambda=\frac{1\pm\sqrt{1-4.8}}{2} =0.5\pm0.9747i.

Their magnitude is $\sqrt{0.5^2+0.9747^2}\approx1.095$, greater than one. The error components therefore have an oscillatory pattern with a growing envelope. For instance, starting from $e_{-1}=0$ mg and $e_0=1$ mg gives $e_1=1$ mg, $e_2=-0.2$ mg, and $e_3=-1.4$ mg. These early values illustrate reversals; the root magnitude, rather than a few early points, establishes the instability.

Diagnostics, Alternatives, and Limitations

For a real titration algorithm, compare the observed error trajectory with the predictions of the linearized rule. Repeated sign changes that shrink are consistent with damped oscillation; persistent oscillations suggest proximity to a stability boundary; increasing oscillations are a warning of an unstable gain-delay combination. These patterns are clues, not proof: noisy or sparse measurements, varying response sensitivity, and changes in the dosing algorithm can all alter the observed trajectory.

A useful first diagnostic is to estimate or bound the local slope $g$, calculate $a=Kg$, and compare it with the threshold for the assumed delay. Check that the delay is expressed in the same update steps as the recurrence, and that the sign of the dose adjustment is corrective. If the response slope changes substantially across the dose range, a single local gain cannot characterize the whole trajectory. If doses hit bounds, or if adjustments are rounded, the actual system is no longer the unconstrained linear recurrence analyzed here.

The main alternative is to simulate the full state-update rule using a pharmacometric model, including response dynamics, observation timing, dose bounds, and measurement variability. Such simulation can represent features excluded by this derivation, but its conclusions depend on the model and scenarios chosen. An immediate-feedback calculation ($d=0$) is simpler and can be a useful benchmark, but it can be misleading when measurements are delayed. The linear root test gives a transparent local stability criterion; it does not replace model-based assessment when dynamics are nonlinear or the response is not near the target.

Common pitfall: Seeing that each update has the corrective sign and concluding that the rule must converge. With delay, updates are based on past errors, so a correction can continue after the current response has already crossed the target. The correct check is the magnitude of every characteristic root of the delayed recurrence. Also, equality at the stability threshold is not convergence: for $d=0$ the boundary root is the single value $-1$, while for $d\geq1$ the boundary includes a conjugate pair on the unit circle.

The criterion is local and deterministic. It assumes a fixed linear response slope, a constant gain and delay, and no bounds or noise. It says whether small errors decay under that model; it does not guarantee that a patient-specific trajectory will be safe, that the target is clinically appropriate, or that a dose will remain within feasible limits.

Check Your Understanding

Answer each question, then compare your reasoning with the concise solutions below.

For $d=0$ and $a=0.6$, what is the characteristic root, and does the error converge monotonically or with alternating signs?
For $d=2$, calculate the effective-gain stability threshold to three significant figures.
A one-step-delay rule has $a=0.9$. Is it stable? What does the stability result alone say about whether the trajectory oscillates?
Why can a gain that is stable with immediate feedback become unstable when the measurement delay increases?

Solutions. 1. The root is $\lambda=1-a=0.4$, so errors decay without sign alternation. 2. The threshold is $2\sin(\pi/(4(2)+2))=2\sin(\pi/10)\approx0.618$; gains strictly below it are stable. 3. The one-step threshold is 1, so $a=0.9$ is stable. Stability guarantees decay, but not monotonicity; determining oscillation requires examining the roots, not only their magnitudes. 4. A delayed update uses an earlier error, so corrections can continue after the current response has changed. Mathematically, increasing $d$ lowers the maximum stable effective gain $2\sin(\pi/(4d+2))$.