Surge functions as PET model input

Indicator dilution curves, after intravenous bolus administration, have been observed to be of the form of the gamma variate function (Thompson et al., 1964; Starmer & Clark, 1970), and this has been explained by modelling vascular vessels as series of dilution chambers (Schlossmacher et al., 1967; Davenport, 1983) and with convective dispersion model (Leonard & Jorgensen, 1974; Harpen & Lecklitner, 1984). The general one-parameter n-compartmental model approximation for the delayed response in blood compartment (transit-time model) describes the blood curve with function

(DiStefano III, 2013, p 193). With two compartments (n=2) this becomes similar to the gamma variate ‑based surge function

(with maximum at t=1/λ and 0=A/λ2). Applying logarithm to the equation gives

which can be fitted to data using non-iterative linear regression methods. Surge function can be used in simulations (Herscovitch et al., 1983), but as such it is too simplistic for fitting input curves from PET studies (examples are shown in Figure 1). Input TAC in DSC- or DCE-MRI and contrast-enhanced CT is often fitted using gamma variate functions. Gamma variate function and especially LDRW model function have been shown to fit well many kinds of indicator dilution curves (Mischi et al., 2008; Brands et al., 2011). Sum of two gamma distribution functions can fit the primary bolus and the first recirculation peak (Davenport, 1983). Example of sum of three surge functions is shown in Figure 1).

Examples of surge functions and their sum
Figure 1. Examples of curves that can be produced with the surge function.
The parameters of the curves are A=10, λ=0.1 (red); A=2, λ=0.05 (blue); and A=0.5, λ=0.02 (green). Purple curve is the sum of these three surge functions.


A common practise to model the recirculation phase in the plasma TAC is to use a function that is the sum of surge function and its integral function. The integral of surge function is

and the sum of surge and its integral function is

, where c is a scaling factor for the integral. The integral functions of the surge functions that were shown in Figure 1 are shown in Figure 2 (left panel), and the sums of the surge functions and scaled integral functions shown in Figure 2 (right panel). Similar and more versatile curve shapes can be achieved without the integral function if more surge functions are summed, but that increases also the number of function parameters which complicates the fitting. The same approach can be applied to any function, for example exponentials.

Examples of integrals of surge functions Examples of sums of surge and its integral function
Figure 2. On the left, the integrals of three surge functions (red, blue, and green), and their sum (purple). Surge functions are the same as in Figure 1. On the right, the sums of surge functions and their scaled integral functions (red, blue, and green), and their sum (purple). The scaling factors were given values 0.001, 0.005, and 0.002.

Essentially, the sum of surge function and its integral functions consists of an exponential and surge function. The exponential function has been commonly added to gamma variate function to account for the recirculation, for example in [13N]NH4+ and [15O]H2O studies (Golish et al., 2001; Lüdemann et al., 2006), and in MRI (Parker et al., 2006). Long acquisition time with DCE-MRI causes significant renal clearance of contrast agent, which can be taken into account with yet another exponential term (Duan et al., 2017). These formulations resemble the exponential input function with a pair of repeated eigenvalues which is based on a compartmentalized model of radiopharmaceutical behaviour in the circulatory system ("Model 2" by Feng et al. 1993a and 1993b).

Gamma variate function with a recirculation term can be fitted to input curves and TTACs using fit_gvar.

A simplified version of these functions has been used to fit the late part of input curves in order to reduce the number of blood samples in FDG studies (Phillips et al., 1995) and to allow usage of an analytic method (Bonson et al., 2000):

, where m and n were population means from a larger dataset, leaving only two parameters, a and b, to be fitted from an individual input curve.


The functions given above can represent the plasma curves after bolus injection. To simulate the radiotracer infusion, of duration Tin, a rectangular (boxcar) function can be convolved with the surge function (Eq 1). Since surge function is used as response function, its AUC should be unity, which can be accomplished by giving it in this form:

TAC can be convolved with surge function using convsurg. Convolution operation can be replaced by subtracting integrals of the surge functions:

, which is easy to implement in spreadsheet program. The function that takes into account the initial delay (radiotracer appearance time, Tap) and the duration of radiotracer infusion (Tin) can be calculated in three parts (examples are shown in Figure 3):

Examples of surge functions for infusion and their sum
Figure 3. Examples of surge functions calculated with infusion duration Tin=100. Otherwise the function parameters are the same as in Figure 1. Purple line is the sum of the three functions.

Exponential input function with a pair of repeated eigenvalues

Based on a compartmentalized model of radiopharmaceutical behaviour in the circulatory system, Feng et al. (1993a and 1993b) proposed a formulation of four exponential functions for fitting PET radiotracer PTACs which include both ascending and descending phase in the measured data. Two of the eigenvalues (time constants of the exponentials) can be paired (Feng & Wang, 1991; Wang & Feng, 1992), leading to equation:

, where Tap is the appearance time of radioactivity in the blood, caused by the delay of bolus in the circulation from the IV injection site to the arterial sampling site. A fast bolus injection is assumed. Examples are shown in Figure 4). One of the exponential terms (with A3 and λ3) can obviously be left out (Wang & Feng, 1992; Feng et al., 1993b) if the kinetics are relatively simple. The equation contains surge function and exponentials and bears likeness to the sum of surge function and its integral function.

Examples of Feng M2 functions
Figure 4. Examples of exponential input function with a paired of repeated eigenvalues.
The parameters of the curves are A1=100, λ1=0.1, A2=50, λ2=0.05, and A3=50, λ3=0.01 (red); A1=10, λ1=0.01, A2=20, λ2=0.005, and A3=30, λ3=0.001 (blue); A1=3, λ1=0.005, A2=100, λ2=0.0005, and A3=200, λ3=0.0001 (green). For all curves Tap=50.

This function can be fitted to PET TACs with program fit_feng. Note that the function can be negative in the beginning phase, which would be non-physiological, and redundant with the appearance time Tap. The derivative of the function (when t>Tap) is

To ensure that the derivative function is >0 at t=0, that is, f'(0)>0, this condition has to be met:

Note that all lambdas in this equation have positive values, when concentration is decreasing and approaching zero.

Incorporated with injection schedule

The exponential input function with a pair of repeated eigenvalues given above cannot fit well blood data, if radiotracer is given as a constant infusion (Wong & Feng, 2005). The duration of the infusion, Tin, can be incorporated in the previous function by convolving a rectangular function with it as a response function. Like in the case of sum of exponential functions and surge function, this can be accomplished by simple subtraction of the analytically integrated response function from its delayed (for time Tin) version of itself (Wong & Feng, 2005; Wong et al., 2006). The integral of the function is

The function that takes into account the initial delay (radiotracer appearance time, Tap) and the duration of radiotracer infusion (Tin) can be calculated in three parts (examples are shown in Figure 5):

Examples of Examples of Feng M5 functions
Figure 4. Examples of exponential input function with a paired of repeated eigenvalues, incorporating injection schedule ("M5" in Wong et al., 2006).
The parameters of the curves are A1=200, λ1=0.5, A2=100, λ2=0.1, A3=100, λ3=0.02, and Tin=30 (red); A1=5, λ1=0.02, A2=10, λ2=0.005, A3=20, λ3=0.0001, and Tin=30 (blue); A1=10, λ1=0.1, A2=30, λ2=0.02, A3=10, λ3=0.001, and Tin=120 (green); A1=0.6, λ1=0.01, A2=5, λ2=0.002, A3=3, λ3=0.00001, and Tin=120 (purple).
For all curves Tap=0.

Wong et al (2006) compared the performance of this function to several functions that are aimed for bolus injection studies. Functions were fitted to blood curves that were obtained from human FDG studies, performed with 3 min infusion, collected with arterialized venous sampling. The blood curves were best fitted with this function, and fits were good also with sum of exponential functions where the ascending phase was fitted using polynomial.

See also:


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

Feng D, Wang Z. A three-stage parameter estimation algorithm for tracer concentration kinetic modelling with positron emission tomography. Proceedings, 1991 American Control Conference, vol 2 (1991): 1404-1405. doi: 10.23919/ACC.1991.4791609.

Feng D, Huang S-C, Wang X. Models for computer simulation studies of input functions for tracer kinetic modeling with positron emission tomography. Int J Biomed Comput. 1993a; 32: 95-110. doi: 10.1016/0020-7101(93)90049-C.

Feng D, Wang Z, Huang SC. Tracer plasma time-activity curves in circulatory system for positron emission tomography kinetic modeling studies. IFAC Proc Vol 1993b; 26(2 Pt 3): 175-178. doi: 10.1016/S1474-6670(17)48708-8.

Lüdemann L, Sreenivasa G, Michel R, Rosner C, Plotkin M, Felix R, Wust P, Amthauer H. Corrections of arterial input function for dynamic H215O PET to assess perfusion of pelvic tumours: arterial blood sampling versus image extraction. Phys Med Biol. 2006; 51: 2883-2900. doi: 10.1088/0031-9155/51/11/014.

Wagner JG. Linear pharmacokinetic equations allowing direct calculation of many needed pharmacokinetic parameters from the coefficients and exponents of polyexponential equations which have been fitted to the data. J Pharmacokin Biopharm. 1976; 4(5): 443-467. doi: 10.1007/BF01062831.

Wong K-P, Feng DD. Generalization of a physiological model of input function for PET data analysis. J Cereb Blood Flow Metab. 2005; 25: S625-S626. doi: 10.1038/sj.jcbfm.9591524.0625.

Wong K-P, Huang S-C, Fulham MJ. Evaluation of an input function model that incorporates the injection schedule in FDG-PET studies. 2006 IEEE Nuclear Science Symposium Conference Record; 2006: 2086-2090. doi: 10.1109/NSSMIC.2006.354325.


Definite integrals of surge functions:

Tags: , ,

Updated at: 2021-12-09
Created at: 2016-08-08
Written by: Vesa Oikonen