Optimising Graphical Techniques Applied to Irreversible Tracers
Yasser Alzamil
1,2
,Yulia Hicks
1
, Xin Yang
1
and Christopher Marshall
3
1
Institute of Health, Technology and the Digital World, Cardiff University, School of Engineering,
Queen’s Buildings,The Parade, CF24 3AA, Cardiff, U.K.
2
Diagnostic Radiology Department, University of Hail, Applied Medical Sciences College,
P.O. Box 2440, Hail, Saudi Arabia
3
Wales Research & Diagnostic Positron Emission Tomography Imaging Centre PETIC, Cardiff University,
School of Medicine, UHWMainBuilding, HeathPark, CF14 4XN, Cardiff, U.K.
Keywords: PET, Patlak Graphical Analysis,
18
F-FDOPA Quantification, Image Analysis, Dpetstep Simulation.
Abstract: Graphical analysis techniques are often applied to positron emission tomography (PET) images to estimate
physiological parameters. Patlak analysis is primarily used to obtain the rate constant (K
i
) that indicates the
transfer of a tracer from plasma to the irreversible compartment and ultimately describes how the tracer binds
to the targeted tissue. One of the most common issues associated with Patlak analysis is the introduction of
statistical noise that affects the slope of the graphical plot and causes bias. In this study, several statistical
methods are proposed and applied to PET time activity curves (TACs) for both reversible and irreversible
regions that are involved in the equation. A dynamic PET imaging simulator for the Patlak model was used
to evaluate the statistical methods employed to reduce the bias introduced in the acquired data.
1 INTRODUCTION
This article describes the first experiment in the
optimisation process for 18F- FDOPA quantification
in Parkinson’s disease (PD) images. [18F]FDOPA
dynamic images of the brain are used and Patlak
graphical analysis will be applied first to the PET
data, then several statistical techniques will be
implemented to Patlak equations in terms of reducing
bias and noise as well as improving the accuracy of
the result. The calculation of goodness of fit of the
slope and the standard errors of the regression line
are the reference for the evaluation of result
improvement. Graphical techniques are considered
simple for the PET data analysis. Patlak equations for
irreversible tracers were used to establish another
graphical technique for reversible binding tracers.
These techniques are common and preferred due to
the independence from compartmental modelling,
which is harder compared to the graphical methods.
One issue usually generated with this type of analysis
is the bias introduced, which is caused by the
sensitivity to statistical noise (Logan 2003) that is
usually introduced by using the ordinary least
squares (OLS) method. This could be reduced by
selecting and applying several other methods that
include the feasible generalized linear least squares
(FGLS), total least square (TLS), robust fitting
regression (RFR).
In this article, these methods first will be applied on
simulated phantom PET data generated from
dPETSTEP simulation tool (Haggstrom et al. 2016)
to see how different statistical methods will affect the
image analysis result and how they behave with
different levels of noise. The equations, materials,
algorithms and tools used in the experiment will be
described and the result will be discussed and
evaluated in this article.
2 THEORY
The main thing that distinguishes between reversible
and irreversible tracers is the experiment length; in
other words, a tracer can be reversible over a long
period, but it could be irreversible during experiment
or scanning time (Logan 2000). For irreversible
tracers, the Patlak reference model is used to
calculate the [
18
F]FDOPA uptake in ROI, and it is the
same standard model using the blood data as an input
function except that the blood activity is replaced by
the activity of a reference TAC. The reference must
Alzamil, Y., Hicks, Y., Yang, X. and Marshall, C.
Optimising Graphical Techniques Applied to Irreversible Tracers .
DOI: 10.5220/0006513700170026
In Proceedings of the 11th International Joint Conference on Biomedical Engineering Systems and Technologies (BIOSTEC 2018) - Volume 2: BIOIMAGING, pages 17-26
ISBN: 978-989-758-278-3
Copyright © 2018 by SCITEPRESS Science and Technology Publications, Lda. All rights reserved
17
be devoid of targeted receptors (Patlak et al. 1983;
Patlak and Blasberg 1985). The Patlak method can
be modelled as two-tissue compartments (Figure 1).
Figure 1: Two-tissue compartmental model where Cp(t) is
the plasma tracer concentration at time (t). C1 is the free
tracer concentration in tissue and C2 is the trapped tracer
concentration. K1, k2 and k3 are unidirectional rate
constants of the tracer between plasma and tissues.
With FDOPA tracer, the regions used to generate
graphics include one reversible region (cerebellum)
and one irreversible region (striatum). TACs of
blood plasma usually calculated and used an input
function for Patlak analysis but a reference tissue
region replaces plasma measurement here. Patlak
plot equation with reference tissue as an input
function, which can be described by
C
T(t)
C
T
'(t)
=K
C
T
'
τ
t
0
C
T
'(t)
+V
(1)
K=slope=
k
2
k
3
(k
2
+k
3
)
(2)
Where CT(t) is the TAC values of the ROI tissue
(striatum), CT’(t) is the reference tissue TAC
(cerebellum) at time t, K is the unidirectional net
uptake rate constant and V is the blood volume
fraction . The equation works for t > t*, t* is the
equilibration time where the radioactivity ratio
between reference and ROI tissue becomes
reasonably constant (Ikoma et al. 2008). K represents
the tracer net uptake calculated from the regression
slope, and V is the intercept, which is equal to the
volume fraction of blood in ROI tissue at time 0. This
means measured activity in ROI (striatum) is divided
by the reference tissue activity that represents (y-
axis), and plotted against the integral of the reference
TAC from the injection time divided by the reference
activity, which represents (x-axis). For
18
F-FDOPA
tracer that targets brain receptors in striatum, the
model plot will result in a straight line after t* (Patlak
and Blasberg 1985) as illustrated in (Figure 2).
Figure 2: PET activity measured from ROI is divided by
the reference tissue activity that represents (y-axis), and
plotted against the integral of the reference TAC from the
injection time divided by the reference activity, which
represents (x-axis). The model plot resulted in a straight
line after t* = 30 min in this analysis.
2.1 Ordinary Least Square
OLS is one of the simplest methods of linear
regression and it is frequently used to analyse both
experimental and observational data. It aims to
closely fit a function with the data by minimizing
the sum of squared differences between the observed
responses in the given dataset and those predicted by
a linear function. On a plotted graph, this is seen as
the sum of the squared vertical distances between
each data point in the set and the corresponding point
on the regression line (RL). The smaller the
differences, the better fitted the model of the data
(Hayashi 2000).
In the case of a model with p explanatory variables,
the OLS regression model is expressed as:
 

 
(3)
Where Y is the dependent variable, β
0
is the intercept
of the model, X
j
corresponds to the j
th
explanatory
variable of the model (j=1 to p), and ε is the random
error with expectation 0 and variance σ². With n
observations, the estimation of the predicted value of
the dependent variable Y for the i
th
observation is
given by:
 


(4)
In a model, this leads to the following estimators of
the parameters:
β = (X
T
DX)
-1
X
T
Dy σ² = 1/(W p*)

i
(y
i
- y
i
)
(5)
BIOIMAGING 2018 - 5th International Conference on Bioimaging
18
where β is the vector of the estimators of the βi
parameters, X is the matrix of the explanatory
variables preceded by a vector of 1s, y is the vector
of n observed values of the dependent variable, p* is
the number of explanatory variables to which we add
1 if the intercept is not fixed, w
i
is the weight of the
i
th
observation, W is the sum of the w
i
weights, and
D is a matrix with the w
i
weights on its diagonal. The
vector of the predicted values can be written as
follows (Freedman 2009):
y = X (X
T
DX)
-1
X
T
Dy
(6)
2.2 Generalized Linear Least Squares
Generalized linear least squares (GLLS) was applied
in different PET quantitative model equations
(Logan et al. 2001) to remove bias, and it was
developed originally by Feng et al. (1993). GLLS is
used for estimating the unknown parameters in
models based on a linear regression when calculated
data shows a certain degree of correlation in the
residuals. This method was applied to the data in two
parts: for times 0 to T1 and from T1 to the end time.
The parameters generated a curve used as input to the
linear regression analysis (Logan et al. 2001). The
two types of GLLS are called weighted least squares
(WLS) and feasible generalized least squares
(FGLS). WLS can be applied when all the off-
diagonal entries of the covariance error matrix (W)
are 0. In FGLS, the opposite occurs when the
covariance of errors are unknown (Strutz 2010).
With FGLS, the calculation progresses in two steps:
first, the residuals are obtained by using OLS to
establish the errors covariance matrix that shows a
consistency in estimation. Second is implementing
the idea of GLLS, which is to divide the data given
into two sides, and the part that has the low variance
is given more weight than the other sides to generate
a more accurate fitted line. With finite samples, an
estimator’s accuracy with FGLS can be improved by
an iteration process where residuals are used to
update the errors covariance estimator and
consequently FGLS estimation is updated (Long and
Trivedi 1992; Freedman 2008; Gujarati 2009). The
FGLS estimator may or may not be unbiased in small
samples but if (
) is a consistent estimator of (W),
then the FGLS estimator is asymptotically unbiased,
efficient, and consistent. Monte Carlo studies have
shown that the FGLS estimator generally yields
better estimates than the OLS estimator (Kennedy
2008).
The general linear regression model is defined by
the following set of assumptions:
Linearity in parameters is the functional
form
y = X +
(7)
Error term has mean zero
E() = 0
(8)
Errors are nonspherical
Cov() = E(T) = W
(9)
Where W is any nonsingular TxT
variance-
covariance matrix of disturbances.
Error term has a normal distribution
~ N
Error term is uncorrelated with each
independent variable.
Cov (,X) = 0
There are two types of nonspherical errors: First is
when an error term does not have constant variance,
and this is called heteroscedasticity. In this type of
error, the disturbances are drawn from probability
distributions that have different variances and the
error term has non-constant variance; the variance-
covariance matrix of disturbances is not given by a
constant multiplied by the identity matrix (i.e., W
2
I). Second, the errors are correlated, which is
called autocorrelation or serial correlation where
disturbances are correlated with one another. It
occurs when using time-series data. When the
disturbances are correlated, the variance-covariance
matrix of disturbances is not given by a constant
multiplied by the identity matrix (i.e., W
2
I). This
is because the elements off the principal diagonal of
W, which are the covariance of the disturbances, are
non-zero numbers (Granger 1994).
In a general linear regression model stated in matrix
format, the sample of T multivariate observations
(Y
t
, X
t1
, X
t2
, …, X
tk
) is generated by a process
explained below:
y = X + , ~ N(0, W)
or
y ~ N(X, W)
An FGLS estimator uses the sample of data to obtain
an estimate of W where the true W is replaced with
its estimate
. The FGLS estimator is given by the
rule:
FGLS
= (X
T
-1
X)
-1
X
T
-1
y
Optimising Graphical Techniques Applied to Irreversible Tracers
19
The variance-covariance matrix of estimates for the
GLS estimator is
Cov(
) = (X
T
-1
X)
-1
The FGLS estimator is also a WLS estimator. The
WLS estimated is derived as follows. Find a TxT
transformation matrix P such that μ* = , where μ*
has variance-covariance matrix Cov(μ*) = E(μ*
μ*T) = σ
2
I. This transforms the original error term μ
that is nonspherical into a new error term that is
spherical. Use the matrix P to derive a transformed
model (Granger 1994; Kennedy 2008):
Py = PXβ + Pμ
or
y
*
= X
*
β + μ
*
Where y
*
= Py, X
*
= PX, and μ
*
= Pμ. The FGLS
estimator is the OLS estimator applied to the
transformed model, which is considered a
computational device only to obtain efficient
estimates of the parameters and standard errors of the
original model of interest.
2.3 Robust Fitting Regression
Robustness denotes the solidity of conclusions and
how their differences from assumptions are assigned
to a certain model. This means that small changes in
the data distribution do not cause large changes in the
variance of the estimates (Western 1995). RFR
provides an alternative to OLS when assumptions are
invalid within the model. RFR provides much
improved regression coefficient estimates when data
noise or outliers are existent. The influence of
outliers is down-weighted by making the outlying
residuals larger and simpler to detected, plus
performing an iterative procedure to identify outliers
and reduce the impact on the coefficient estimates.
Robust regression implements its own residual
analysis and reduces or completely removes
numerous data points, so a decision should be made
whether these observations are essential in the
analysis (Hintze 2001; Kutner et al. 2004). The most
common general method of robust regression is a
class of techniques called M-estimators that discount
the impact of outlying observations (Fox 2002),
introduced by Huber (1964). Consider the linear
model
y
i
= α + β
1
xi
1
+ β
2
xi
2
+ . . . + β
k
x
ik
+ ε
i
= x’
i
β + ε
i
The fitted model for the ith of n observations is
y
i
= a + b
1
x
i1
+ b
2
x
i2
+ . . . + b
k
x
ik
+ e
i
= x’
i
b + e
i
The general M-estimator minimizes the objective
function
ρ
e
i
=
n
i=1
ρy
i
-x'
i
b
n
i=1
The function ρ gives the contribution of each residual
to the objective function and should have the
following properties:
ρ(e) 0
ρ(0) = 0
ρ(e) = ρ(−e)
ρ(e
i
) ρ(e
i’
) for |e
i
| > |e
i’
|
For least-squares estimation, ρ(e
i
)=e
i
2
, let ψ=ρ be
the derivative of ρ, differentiating the objective
function while considering the coefficients b, and
assume the partial derivatives to be 0, which
enervates a system of k+1 estimating equations for
the coefficients:
ψy
i
-
bx'=0
n
i=1

Assume the weight function w(e) = ψ(e)/e , and let
w
i
= w(e
i
), then the estimating equations can be
written as
w
i
y
i
-
bx'=0
n
i=1
An iterative solution called iteratively reweighted
least-squares (IRLS) is required due to the
dependency between weights, residuals and the
estimated coefficients. The iteration is performed
following these steps:
1. Select initial estimates b(0); for example, the least-
squares estimates.
2. At each iteration t, calculate residuals e
i
(t-1)
and related weights w
i
(t-1)
=w[e
i
(t-1)
]
from the previous iteration.
3. Apply for new weighted-least-squares estimates
b
t
=[X'W
(t-1)
X]
-1
X
'
(t-1)
y
BIOIMAGING 2018 - 5th International Conference on Bioimaging
20
Figure 3: Ordinary least squares (Left) and total least squares (Right) fits of a set of
m
= 20 data points in the plane. (---)
data points, 

, x- approximations 

, solid line () fitting model
, dashed lines (---) approximation errors.
X is the model matrix, with
as its i
th
row, and



is the current weight matrix.
These steps are repeated till the estimated coefficients
converge and covariance matrix of b is
Vb=
E(ψ
2
)
[E(ψ
'
)]
2
(X'X)
-1
(23)
Using
[ψ(e
i
)]
2
to estimate E(ψ
2
) and

ψ '
e
i
/n
2
to estimate
[Eψ']
2
produces the estimated
covariance matrix V
(b).
2.4 Total Least Squares
The TLS method, known also as error-in-variables
method or orthogonal regression method, is a general
approach that can be used in n-dimensional space
(Petras and Podlubny 2010). Many areas of
application use the TLS method such as signal
processing, image processing and economics. The
orthogonal distance, distance between data point and
fitted line, is the main category of TLS and
mathematically can be expressed by the following
relation (Petráš and Bednárová 2010):
R= |d
i
|
n
i=1
(24)
d is the orthogonal distance and the target is to find a
minimum of R; the TLS approach minimizes the sum
of the squared d from the data points to the fitting line.
With the TLS method, well-known mathematical
tools are usually used.
For linear regression model of the expression:
y = bx + a
(25)
The coefficients a and b can be derived from the
following relations:
a=
y
i
n
i=1
-b
x
i
n
i=1
n
=y-bx
B=
1
2
y
i
2
n
i=1
-ny
2
-(
x
i
2n
i=1
-nx
2
)
nxy-
x
i
y
i
n
i=1
b=-
B
2
+1
The OLS and TLS methods assess the fitting accuracy
in different ways: the OLS method minimizes the sum
of the squared vertical distances from the data points
to the fitting line, whereas the TLS method minimizes
the sum of the squared d from the data points to the
fitting line. Figure 3 shows OLS and TLS fitting lines
as well as the data approximation. In the least squares
case, the data approximation is obtained by correcting
the second coordinate only. In TLS, the data
approximation is obtained by correcting both
coordinates (Markovsky and Van Huffel 2007). This
method takes into account the noise in the
independent as well as dependent variables (Varga
and Szabo 2002).
3 MATERIALS AND METHODS
3.1 Scanner and Reconstruction
Parameters
The scanner available in the PET Imaging Centre at
the University Hospital of Wales is GE Discovery
690 PET/CT (General Electric Healthcare), and it is
made for clinical and research use. The reconstruction
method used with FDOPA imaging is a maximum
likelihood ordered subset estimation maximisation
(ML OSEM) with Vue Point FX algorithm with time
of flight (TOF) correction.
Optimising Graphical Techniques Applied to Irreversible Tracers
21
3.2 dPETSTEP and Simulation
Parameters
dPETSTEP is a dynamic PET simulator that can
generate PET brain images and allows full simulation
of kinetic modelling (Haggstrom et al. 2016).
Additionally, the dynamic PET data can then be
model-fitted to produce physiological parameter
estimates. dPETSTEP uses MATLAB as a platform
and works as an extension of another application
named PETSTEP, and they both share some
commands. Both applications are open source and
available on the GitHub© website. dPETSTEP is a
fast and simple tool and can be used as an alternative
to Monte Carlo (MC), as it is 8000 times faster than
MC. In terms of kinetic analysis, dPETSTEP is very
helpful for the evaluation of different processing
choices with dynamic and parametric images that
usually require a long time for image analysis where
each pixel is analysed. To run the simulation, at the
beginning, the settings of the simulation must be
adjusted to be consistent with the PET scanner
features that are targeted by the analysis or
evaluation. Table 1 shows all required simulation
settings to generate simulated dynamic PET data for
18
F-FDOPA images that are performed for a PD
patient. The simulated dynamic PET scan will be
acquired as 26 time frames over 94.5 min (1 × 30 sec,
4 × 1 min, 3 × 2 min, 3 × 3 min, and 15 × 5 min) with
a dose of 111 MBq. The cerebellum region from each
time frame will be used as the reference region tissue
(input function) for the quantitative analysis. The
output after running the simulator is a 4D matrix with
a reconstructed dynamic ordered subset expectation
maximization (OSEM) image with point spread
function (PSF) correction, which is similar in features
to the dynamic PET data produced by the real
scanner. The noise level of the 4D matrix is kept at a
minimum in its sinogram. The entire volume of
images is separated into single slices as 2D images
(Figure 4), and these are saved to a MATLAB version
4 MAT-file, which allows them to be loaded in any
MATLAB version with 2D double, character, sparse
arrays and without compression. This version allows
for each variable to have a maximum size of 100
million elements per array and 231 bytes per variable,
which allows all the data information to be kept
without losing any details. Scatter and decay
correction was applied during the reconstruction
process obtained from setting the parameters file
within the simulation codes. All slices introduced by
the 4D matrix are exported to PMOD 3.4 software to
perform further analysis in kinetic modelling with
specific ROI.
Figure 4: (Left) The entire volume of images is separated
into single slices as 2D images, which are saved to a
MATLAB version 4 MAT-file, which allows them to be
loaded to PMOD software (Right).
3.3 PMOD Software
PMOD software has been designed for researchers in
the molecular imaging field and provides suitable
tools for all quantitative data processing steps. This
will make researchers pay more attention to the
content and clinical data rather than programming
new tools from scratch. PMOD can process all type
of images such, as CT, MR and SPECT, in many
different imaging formats, from simple processing
tasks to sophisticated protocols and analysis all
packed into one graphical user interface application.
It is validated and refers to more than 1000 peer-
reviewed publication documents on kinetic modelling
and biomedical research (LLC 2017). In this
experiment, two applications are used for the
analysis, PBAS and PKIN, which allow import of the
simulated dynamic PET images and applying kinetic
analysis, respectively. PBAS is the main tool in
PMOD and can receive images in different proper
formats, one of which is v4 MAT-files, which should
be imported to PBAS as double SE data form. The
images imported must be in Bq/cc units, and the
orientation should be corrected if it looks different.
All corrections are already made by the simulator, and
PMOD should do nothing in terms of corrections. At
the beginning, all slices are merged into frames again,
and the times for each frame are consistent with the
protocol used with the simulator. After that, all
frames are saved as DICOM format, the Digital
Imaging and Communications in Medicine standard.
This format is data rich, and the header information
includes attenuation, scatter and decay correction,
normalisation, frame timing and reconstruction
parameters as well as the standard, required details of
matrix dimensions and pixel size. The computer
which was used for reconstruction of the images was
an Intel® Core™ I7-6700K processor 8MB Cache
4.00 GHz 4 Cores 8, with Windows 7 Enterprise 64-
bit OS and 32 GB RAM.
BIOIMAGING 2018 - 5th International Conference on Bioimaging
22
PBAS was used to draw and analyse the ROI (left
striatum) and reference tissue (cerebellum), the
reason for choosing the left striatum is that it is
assumed that in simulated images both the left and
right side will have the same amount of radioactivity.
The net influx rate constant (Ki) is usually calculated
over a range of 30-90 min for dynamic images, so all
frames between 30 and 90 minutes are averaged into
slices, and then the ROIs are drawn on all slices to
show them. Two volumes of interest (VOIs) at
minimum must be created; one VOI for the
cerebellum represents the reversible region, and
another one for the left striatum represents the
irreversible region. For regional analysis, the ROI is
drawn as one complete contour; ROI objects were
drawn freehand (Figure 5).
Figure 5: For regional analysis, the ROI is drawn as one
complete contour; ROI objects were drawn freehand;
simulated striatum (Left) and part of simulated cerebellum
(Right).
Each VOI includes contours that are drawn over the
simulated anatomical region required for the analysis,
and all VOIs are saved in a single file to be used again
as a template. In fact, even the ROIs contours are
drawn on the averaged dynamic frames; however, the
TACs at the end will be obtained from the original
dynamic PET volume. So, the images are averaged to
help guiding the drawing on all slices that show ROIs.
Then, the time activity curves (TACs) are generated
from the contours of regions required for the analysis.
TACs for all regions must be checked, and the tracer
in the reversible curve should show a clear washout.
PMOD has an algorithm to determine the best time
for t*. In this experiment, the Patlak reference tissue
analysis is based on t* = 30 min, and this is important
for slope calculations (Ki value). For regional
analysis, the result will be (Ki) value (slope). For each
specific ROI made, the intercepts that represent (V),
standard errors (SE) and Chi-square are calculated
within the analysis. The maximum standard error
accepted for (Ki) is 10%.
The simulation of the FDOPA dynamic PET data is
repeated 10 times and simulates different (Ki) values
within the normal healthy subjects, which are
between 0.012 and 0.014 according to literature,
considering that the imaging protocol of FDOPA for
scanning patients is based on giving both carbidopa
and entacapone. After finishing the kinetic modelling
analysis, Gaussian noise with a zero mean was
applied in 10 levels (0.05 0.5) for both the reference
tissue and ROI. All the TACs where then exported to
be analysed by our statistical calculation methods,
FGLS, TLS and RFR, which were explained before.
Codes were written in MATLAB software to analyse
all TACs obtained from the simulated images
processed by PMOD, and the (Ki) results were saved
to Excel sheets. The following section illustrates the
results for the simulated data, including statistics,
which were calculated to evaluate how the methods
affected the estimate of the final (Ki) values.
4 RESULTS
TACs for reference tissue (cerebellum) and ROI (left
striatum) that were obtained from the first simulation
data without noise are plotted as an example in Figure
6. The curves show that activity in the cerebellum
starts to increase significantly, reaching 12000
Kbq/cc, and then decreases (washout) to around 4000
Kbq/cc and starts what is known as the equilibration
stage, which is the time point most preferable for
Patlak analysis. The activity washout confirms the
reversibility of the cerebellum region. In contrast, the
left striatum TAC increased until it reached a stable
level for the remainder of the scan, which confirms
the irreversibility. For OLS, FGLS, TLS and RFR
calculation, the results include noisy and non-noisy
data. After fitting the data in the statistical models,
two numerical methods used to evaluate the goodness
of fit for all linear regression analyses were included:
the sum of squares due to error (SSE) and the standard
error (SE), known as the root mean squared error
(RMSE) as well. The numerical measures are more
closely concentrated on a specific aspect of the data
and often try to compress that information into a
single number. SSE measures the total deviation of
the response values from the fit, where a value closer
to 0 indicates a smaller random error component
within the model selected, and the fit will be more
useful for prediction. It is also known as the summed
square of residuals.
Optimising Graphical Techniques Applied to Irreversible Tracers
23
SSE= (y
i
-y
i
n
i=1
)
2
(29)
Where y
i
is the i
th
value of the variable to be
predicted,
is the predicted value of y
i
. SE is an
estimate of the standard deviation of the random
components in the data and is similar to the SSE. An
SE value closer to 0 indicates a fit that is more useful
for estimation. SE is described as
SE=
MSE,MSE=
SSE
v
,v=n-m
(30)
Where MSE is the mean square error or the residual
mean square and v indicates the number of
independent pieces of information involving the n
data points that are required to calculate the sum of
squares, which is calculated as the number of
response values n minus the number of fitted
coefficients m estimated from the response values.
Figure 6: TACs for the reference tissue (cerebellum) and
ROI (left striatum) that were obtained from the first
simulation data without noise are plotted as an example.
A one-way repeated measured analysis of variance
(ANOVA) was conducted to evaluate the null
hypothesis that there is no change in the simulations’
Ki value when calculated with different statistical
regression models in all four groups (N = 10). The
results of the ANOVA indicated a significant effect,
Wilks’ Lambda = .001, F(3, 7) = 41.17, p < 0.01, η2
= .75. Thus, there is significant evidence to reject the
null hypothesis. Follow-up comparisons indicated
that each pairwise difference was significant, p <
0.01, except between OLS and TLS, where p = 0.09.
There was a significant difference between the
statistical models used, suggesting that using a
different linear regression model reveals a different
Ki value, which is the main parameter of diagnosis. A
repeated measures ANOVA is performed when the
samples are considered related (dependent) and more
than two groups. Table 2 illustrates te SE and SSE for
each method applied to simulation data. With the
OLS method, the min SE in all simulations including
noisy data was 0.05680 and the max was 1.01070; the
average SE was 0.23551. For SSE, the min was
0.00029, max 0.09286 and the average 0.00741. The
fold change range was between 0.00041 and 2.09565.
FGLS calculations with the data showed the SE min
was 0.0012377, max 0.0105768 and the average
0.0041978. The SSE min was 0.0000001, max
0.0000102 and the average 0.0000020, with a general
fold change range of 0.0032366-2.3815366. For TLS
analysis, the min SE was 0.05337, max 18.79341 and
the average SE in all data calculations was 2.86733.
The SSE min was 0.00026, max 32.10839 and
average 1.58100. The fold change range was
0.00506-2.09567. The min SE in the RFR analysis
was 0.07026, max 0.92162 and the average 0.37945.
The SSE min was 0.00045, max 0.07722 and the
average 0.01797 with a fold change range of 0.00755-
2.83829.
Table 1: A one-way repeated measured analysis of variance
(ANOVA) result.
Test
Wilks’
Lambda
F(3, 7)
p -value
ANOVA
(N=10)
.001
41.17
p < 0.01
pairwise
difference
OLS and
TLS, pairwise
η
2
p < 0.01
p = 0.09
.75
Table 2: SE and SSE values for each method applied to
simulation data.
Stat. Model
SE
SSE
fold
change
OLS
Min
0.0568
0.0003
0.0003
-
0.7927
Max
1.0107
0.0929
Avg.
0.2355
0.0074
FGLS
Min
0.0012
0.1 

0.0049
-
2.5119
Max
0.0106
10.2 

Avg.
0.0042
2.0 

TLS
Min
0.0697
0.0004
0.0107
-
2.0719
Max
1.3071
0.1553
Avg.
0.4505
0.0237
RFR
Min
0.0703
0.0005
0.0196
-
1.6780
Max
0.9216
0.0772
Avg.
0.3795
0.0180
5 DISCUSSION
For most of the Ki values, it increases in proportion
to the level of noise. SE and SSE also increase with
higher noise level. In simulations 1-5, the data
behaved quite similarly with the noise effect
compared to simulations 6-10, which had more
scattering, and this could have been caused by
BIOIMAGING 2018 - 5th International Conference on Bioimaging
24
increasing the Ki amount for the data without noise.
A one-way repeated measurements ANOVA test
indicated that there is significant evidence to reject
the null hypothesis and confirm that those statistical
models reveal different results for the final Ki value.
Each statistical regression model dealt with the
simulation data in a different way, and based on the
goodness of fit evaluations, the best fitted regression
model can be chosen. Follow-up comparisons
indicated that each pairwise difference was
significant except between OLS and TLS. The p-
value between OLS and TLS shows that null-
hypothesis about equality can not be rejected at
confidence level 0.05 that is used with this test. The
significant difference between the statistical models
used, suggesting that using a different linear
regression model reveals a different Ki value, which
is the main parameter of Parkinson’s disease
diagnosis. The min SE and SSE were found with
FGLS, and this suggests that FGLS is the best of these
models to fit the noisy data. FGLS also had the second
lowest fold change rate, which measures the change
difference from the Ki value without noise. OLS and
TLS showed the lowest fold change rate, and that
indicates more resistance to the noise effect than the
other methods. RFR had the highest fold change rate
compared to the data without noise, which indicates
sensitivity to noise.
In the first simulation, the OLS Ki data had similar
values at all levels of noise, unlike the other data for
the rest of the simulations. For the first two
simulations, all statistical regression models behaved
similarly at all noise levels; thus, the Ki data points
were close in value. This could happen in link to the
original Ki values without noise where considered the
lowest among the other simulations. Other
simulations with low Ki values below 0.012, which
represents the diseased levels, could reveal more
details about how statistical models will behave. In
addition, Ki values lower than 0.012 could clarify
whether the noise could lead to moving the Ki values
to be considered normal. The big difference within
the SE from one noise level to another indicates the
high sensitivity to noise in linear regression analysis,
which confirms the previous results in the published
literature review. The experiment has contributed to
the knowledge in the field by suggesting using FGLS
as a linear regression method for Patlak graphical
analysis. From the above results, using FGLS could
provide better data fitting with low SE and SEE,
indicating more accurate results than other
approaches. Including lower Ki values without noise
could give more details and better confirms the best
fitting model from regression methods. Repeating the
experiment with different equivalent time (t*) points
could reveal more details and alter the accuracy of
those methods. It would be useful to apply those
models on clinical data obtained from patient
dynamic images and compare the results to the
outcome of this experiment.
6 CONCLUSION
Graphical techniques in PET data analysis from
reference tissue are considered simple means of
analysis to obtain physiological parameters. One
issue usually generated with this type of analysis is
the bias introduced, which is caused by the sensitivity
to statistical noise. The study shows there is a
significant difference between the statistical models
used, suggesting that using a different linear
regression model reveals a different Ki value, which
is the main parameter of Parkinson’s disease
diagnosis. Analysing the PET data with different
statistical regression approaches and evaluating each
approach graphically and numerically could improve
the final result for more accurate diagnosis. FGLS
regression model could provide better data fitting in
Patlak analysis with low SE and SEE, indicating more
accurate results than other approaches. In the next
work, these methods will be applied to clinical PET
data obtained from a clinical trial to see how the final
result will be affected.
REFERENCES
Feng, D., Wang, Z. and Huang, S.-C. 1993. A study on
statistically reliable and computationally efficient
algorithms for generating local cerebral blood flow
parametric images with positron emission tomography.
Medical Imaging, IEEE Transactions on 12(2), pp. 182-
188.
Feng, D. et al. 1996. An unbiased parametric imaging
algorithm for nonuniformly sampled biomedical system
parameter estimation. Medical Imaging, IEEE
Transactions on 15(4), pp. 512-518.
Fox, J. 2002. Robust regression. An R and S-Plus
companion to applied regression.
Freedman, D. A. 2008. On regression adjustments to
experimental data. Advances in Applied Mathematics
40(2), pp. 180-193.
Freedman, D. 2009. Statistical models: theory and practice.
Cambridge: Cambridge University Press.
Granger, C. W. J. 1994. Learning and practicing
econometrics: Griffiths, We, Hill, Rc, Judge, Gc.
Journal of Economic Literature 32(1), pp. 115-122.
Gujarati, D. N. 2009. Basic econometrics. Location: Tata
McGraw-Hill Education.
Optimising Graphical Techniques Applied to Irreversible Tracers
25
Haggstrom, I., Beattie, B. J. and Schmidtlein, C. R. 2016.
Dynamic PET simulator via tomographic emission
projection for kinetic modeling and parametric image
studies. Med Phys 43(6), p. 3104.
Haggstrom, I. et al. 2016. A fast dynamic PET simulator for
improved kinetic modeling estimation. Journal of
Nuclear Medicine 57(supplement 2), pp. 477-477.
Hayashi, F. 2000. Econometrics. Princeton, New Jersey:
Princeton University Press.
Hintze, J. 2001. NCSS and Pass. Number Cruncher
Statistical Systems, Kaysville.
Huber, P. J. 1964. Robust estimation of a location
parameter. The Annals of Mathematical Statistics
35(1), pp. 73-101.
Ikoma, Y. et al. 2008. PET kinetic analysis: error
consideration of quantitative analysis in dynamic
studies. Annals of Nuclear Medicine 22(1), pp. 1-11.
Kennedy, P. 2008. A guide to econometrics. 6th ed.
Malden, Massachusetts: Blackwell Pub.
Kutner, M. H., Nachtsheim, C. and Neter, J. 2004. Applied
linear regression models. Location: McGraw-
Hill/Irwin.
LLC, P. T. 2017. PMOD Biomedical Image Quantification.
[Format type]. Available at:
http://www.pmod.com/web/?page_id=476 [Accessed:
date].
Logan, J. 2000. Graphical analysis of PET data applied to
reversible and irreversible tracers. Nuclear Medicine
and Biology 27(7), pp. 661-670.
Logan, J. 2003. A review of graphical methods for tracer
studies and strategies to reduce bias. Nuclear Medicine
and Biology 30(8), pp. 833-844.
Logan, J. et al. 2001. A strategy for removing the bias in the
graphical analysis method. J Cereb Blood Flow Metab
21(3), pp. 307-320.
Long, J. S. and Trivedi, P. K. 1992. Some specification tests
for the linear regression model. Sociological Methods
& Research 21(2), pp. 161-204.
Markovsky, I. and Van Huffel, S. 2007. Overview of total
least-squares methods. Signal Processing 87(10), pp.
2283-2302.
Patlak, C. S. and Blasberg, R. G. 1985. Graphical
evaluation of blood-to-brain transfer constants from
multiple-time uptake data: generalizations. J Cereb
Blood Flow Metab 5(4), pp. 584-590.
Patlak, C. S., Blasberg, R. G. and Fenstermacher, J. D.
1983. Graphical evaluation of blood-to-brain transfer
constants from multiple-time uptake data. J Cereb
Blood Flow Metab 3(1), pp. 1-7.
Petráš, I. and Bednárová, D. 2010. Total least squares
approach to modeling: a MATLAB toolbox. Acta
Montanistica Slovaca 15(2), p. 158.
Petras, I. and Podlubny, I. 2010. Least Squares or Least
Circles? A comparison of classical regression and
orthogonal regression. Chance 23(2), pp. 38-42.
Strutz, T. 2010. Data fitting and uncertainty: a practical
introduction to weighted least squares and beyond.
Location: Vieweg and Teubner.
Varga, J. and Szabo, Z. 2002. Modified regression model
for the Logan plot. Journal of Cerebral Blood Flow &
Metabolism 22(2), pp. 240-244.
Western, B. 1995. Concepts and suggestions for robust
regression analysis. American Journal of Political
Science 39(3), pp. 786-817.
BIOIMAGING 2018 - 5th International Conference on Bioimaging
26