Compartmental model ODEs

Compartmental models for PET data can be described in terms of a set of linear, first-order, constant-coefficient, ordinary differential equations (ODEs). Change of tracer concentration in one of the tissue compartments is a linear function of the concentrations in all other compartments:

$\frac{dC_{i}(t)}{dt} = f_{i}(C_{0}(t), C_{1}(t), C_{2}(t), ...) \nonumber$

From the linearity follows that kinetic measurements are the convolution of the tracer input function and the response function of the system.

Convolution:
$a(t) \otimes b(t) = \int_{0}^{t} a(\tau) b(t-\tau) d\tau$

For example, for two-tissue compartment model with tissue compartments C₁ and C₂, input function C₀(t), and rate constants K₁, k₂, k₃, and k₄, can be described mathematically with a system of two ODEs:

$\frac{dC_{1}(t)}{dt} = K_{1} C_{0}(t) - (k_{2} + k_{3}) C_{1}(t) + k_{4} C_{2}(t)$ $\frac{dC_{2}(t)}{dt} = k_{3} C_{1}(t) - k_{4} C_{2}(t)$

The simulated tissue TACs are used to

test and validate analysis methods and software
estimate compartmental model rate constants (model fitting).

Analytic solutions using Laplace transformation

Analytic solutions for linear first-order differential equations can be derived using Laplace transformation. For this 2TCM example the solution is (Phelps et al., 1979):

$\begin{align*} C_1(t) &= \frac{K_1}{\alpha_2 - \alpha_1} \left[(k_4-\alpha_1)e^{-\alpha_{1}t} + \ (\alpha_2 - k_4)e^{-\alpha_{2}t}\right] \otimes C_{0}(t) \\ C_2(t) &= \frac{K_1 k_3}{\alpha_2 - \alpha_1} \left[e^{-\alpha_{1}t} - \ e^{-\alpha_{2}t}\right] \otimes C_{0}(t) \\ &\text{, where} \\ &\alpha_1 = \frac{k_2 + k_3 + k_4 - \sqrt{(k_2 + k_3 + k_4)^2 - 4 k_2 k_4}}{2} \\ &\alpha_2 = \frac{k_2 + k_3 + k_4 + \sqrt{(k_2 + k_3 + k_4)^2 - 4 k_2 k_4}}{2} \end{align*}$

, where ⊗ denotes the operation of convolution. For one-tissue compartmental model the solution is a lot simpler:

$C_1(t) = K_1 C_{0}(t) \otimes e^{-k_{2}t}$

For simplified reference tissue model the solution is:

$C_T(t) = R_1 C_{R}(t) + \left( k_2 - \frac{R_1 k_2}{1 + BP_{ND}} \right) C_{R}(t) \otimes \ e^{-\frac{k_2}{1 + BP_{ND}} t}$

Direct solutions for discrete-time data

Direct numerical approaches based on simple (Eulerian) integration are conceptually simpler and computationally more effective than analytic solutions requiring convolution operations (Graham, 1985). Adams-Bashforth Predictor requires present and past values of x (x^m and x^m-1) to predict x^m+1:

$x^{m+1} = x^{m} + \frac{\Delta t}{2}(3 f^m - f^{m-1})$

Δt is the sample time difference. Fourth-order Adams Predictor requires four previous values to predict x^m+1:

$x^{m+1} = x^{m} + \frac{\Delta t}{24}(55 f^m - 59 f^{m-1} + 37 f^{m-2} - 9 f^{m-3}) \\ \text{, where} \\ f^{m-i} = f(x^{m-i}, t_{m-i}) \quad \text{for}\ i = 1, 2, 3 \text{.}$

As a fourth-order algorithm its errors are similar to the fourth-order Runge-Kutta algorithm (DiStefano, 2013). The implicit Second-order Adams-Moulton Corrector requires only one value of the past:

$x^{m+1} = x^{m} + \frac{\Delta t}{2}(f^{m+1} + f^{m}) \\ \text{, where} \\ f^{m+1} = f(x^{m+1}, t_{m+1})$

The Second-order Adams-Moulton method has been applied to PET data modelling (Kuwabara et al., 1993). Solutions for the most common compartmental models using this method are shown in TPC modelling reports.

As an example, the solution of the one-tissue compartmental model ODE,

$\frac{dC_{1}(t)}{dt} = K_{1} C_{0}(t) - k_{2} C_{1}(t)$

starts by integrating the equation, assuming that at time zero all concentrations are zero:

$C_{1}(T) = K_{1} \int_{0}^{T} C_{0}(t)dt - k_{2} \int_{0}^{T} C_{1}(t)dt$

Integrals are calculated using trapezoidal method. Concentrations in tissue compartments are calculated starting at time zero when all concentrations and integrals are assumed to be zero, and then proceeding one input function sample at a time, using the results from previous sample to calculate the next ones. Definite integral of n th compartment is implicitly estimated (Kuwabara et al., 1993) based on second-order Adams-Moulton method with Euler integration:

$\int\limits_{0}^{T} C_{n}(t)dt = \int\limits_{0}^{T-\Delta t} C_{n}(t)dt + \ \frac{\Delta t}{2} C_{n}(T - \Delta t) + \frac{\Delta t}{2} C_{n}(T)$

In this example we have just one compartment (n=1). After substitution and rearrangement:

$C_{1}(T) = \frac{K_1 \int\limits_{0}^{T} C_{0}(t)dt \ - k_2 \left[ \int\limits_{0}^{T - \Delta t} C_{1}(t)dt + \ \frac{\Delta t}{2} C_{1}(\small{T - \Delta t}) \right] }{1 + \frac{\Delta t}{2} k_2}$

For simplified reference tissue model the solution is:

$C_{T}(T) = \frac{R_1 C_{R}(T) + k_2 \int\limits_{0}^{T} C_{R}(t)dt \ - \frac{k_2}{1+BP_{ND}} \left[ \int\limits_{0}^{T - \Delta t} C_{T}(t)dt + \ \frac{\Delta t}{2} C_{T}(\small{T - \Delta t}) \right] }{1 + \ \frac{\Delta t}{2} \frac{k_2}{1+BP_{ND}}}$

Using this method the TACs from compartmental models can be easily simulated even in spreadsheet software. For an example, you can download an Excel file where PET tissue data can be simulated using two-tissue compartmental model from measured plasma data and user definable rate constants. Excel files for simulating tissue data with three-tissue compartmental models are also available (compartments in parallel and in series).

Linearization

Some compartment model ODEs can be linearized into multilinear equations from which the model parameters can be solved using standard linear optimization methods. These methods may produce somewhat biased results with noisy data, but the computations are very fast, enabling pixel-by-pixel calculation of parametric images.

Comparison of solution methods

The convolution method provides precise results only if the sampling time interval is very short (relative to the rate constants). The alternative method is computationally faster, and produces correct results (as far as it is possible from discrete data) without any interpolation; the convolution method only approaches the correct result with shorter interpolated sampling intervals (that is, smaller Δt is required in the convolution method). As an example, convolution and the alternative Adams-Moulton method are compared in a simulation of one-tissue compartmental model.

With data sampled at discrete times some errors in simulation will always be present. The longer the sampling time frame durations are, the more important it is to appropriately account for the frame durations in calculating AUC (integral).

Literature

Budinger TF, Huesman RH, Knittel B, Friedland RP, Derenzo SE (1985): Physiological modeling of dynamic measurements of metabolism using positron emission tomography. In: The Metabolism of the Human Brain Studied with Positron Emission Tomography. (Eds: Greitz T et al.) Raven Press, New York, 165-183.

DiStefano III J. Dynamic Systems Biology Modeling and Simulation. Academic Press, 2013. ISBN: 9780124104112.

Graham MM. Model simplification: complexity versus reduction. Circulation 1985; 72(5 Pt 2):IV63-IV68. PMID: 4053329.

Gunn RN, Gunn RS, Cunningham VJ. Positron emission tomography compartmental models. J Cereb Blood Flow Metab. 2001; 21: 635-652. doi: 00004647-200106000-00002.

Huang SC, Phelps ME (1986): Principles of tracer kinetic modeling in positron emission tomography and autoradiography. In: Positron Emission Tomography and Autoradiography: Principles and Applications for the Brain and Heart. (Eds: Phelps,M; Mazziotta,J; Schelbert,H) Raven Press, New York, 287-346.

Kuwabara H, Cumming P, Reith J, Léger G, Diksic M, Evans AC, Gjedde A. Human striatal L-DOPA decarboxylase activity estimated in vivo using 6-[¹⁸F]fluoro-DOPA and positron emission tomography: error analysis and application to normal subjects. J Cereb Blood Flow Metab. 1993; 13: 43-56. doi: 10.1038/jcbfm.1993.7.

Phelps ME, Huang S-C, Hoffman EJ, Selin C, Sokoloff L, Kuhl DE. Tomographic measurement of local cerebral glucose metabolic rate in humans with [F-18]2-fluoro-2-deoxy-D-glucose: validation of method. Ann Neurol. 1979; 6: 371-388. doi: 10.1002/ana.410060502.

Tags: Modeling, Compartmental model, Simulation, ODE, Convolution

Updated at: 2024-01-15
Created at: 2004-06-03
Written by: Vesa Oikonen