PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Phys Med Biol. Author manuscript; available in PMC 2010 April 28.
Published in final edited form as:
PMCID: PMC2860864
NIHMSID: NIHMS77389

Image Covariance and Lesion Detectability in Direct Fan-Beam X-Ray Computed Tomography

Abstract

We consider noise in computed tomography images that are reconstructed using the classical direct fan-beam filtered backprojection algorithm, from both full and short-scan data. A new, accurate method for computing image covariance is presented. The utility of the new covariance method is demonstrated by its application to the implementation of a channelized Hotelling observer for a lesion detection task. Results from the new covariance method and its application to the channelized Hotelling observer are compared with results from Monte Carlo Simulations. In addition, the impact of a bowtie filter and x-ray tube current modulation on reconstruction noise and lesion detectability are explored for full-scan reconstruction.

1. Introduction

Since noise is an important component of image quality, accurate estimation of image noise is crucial for the evaluation of image reconstruction algorithms. Two complementary ways that image noise may be statistically quantified are through image variance and image covariance. In x-ray computed tomography (CT), much research has focused on the calculation of image variance. For the case of parallel-beam filtered backprojection (FBP) reconstruction, some notable examples from the literature that consider image variance are (Tanaka and Iinuma, 1975; Tanaka and Iinuma, 1976; Brooks and Di Chiro, 1976; Barrett et al., 1976; Chesler et al., 1977; Gore and Tofts, 1978; Riederer et al., 1978; Faulker and Moores, 1984; Kak and Slaney, 1988; Davis, 1994). In the case of fan-beam FBP reconstruction, image variance was examined by (Huesman et al., 1977; Herman, 1980; Bennett et al., 1984; Bennett and Byer, 1984; Bennett and Byer, 1986; Pan, 1999; Pan and Yu, 2003; Yu and Pan, 2003; Zeng, 2004; Wang et al., 2005; Zhu and Starlack, 2007). Some of these variance formulas have been validated for accuracy either experimentally or using Monte Carlo studies, e.g. (Bennett et al., 1984; Bennett and Byer, 1984; Bennett and Byer, 1986; Davis, 1994; Pan and Yu, 2003; Yu and Pan, 2003; Wang et al., 2005; Zhu and Starlack, 2007). On the other hand, image covariance does not seem to have received the same amount of attention. Image covariance for FBP reconstruction has been investigated by (Tanaka and Iinuma, 1975; Tanaka and Iinuma, 1976; Alvarez and Stonestrom, 1979; Duerinckx and Macovski, 1979; Duerinckx and Macovski, 1980; Tanaka and Murayama, 1982; Huesman, 1984). However, to our knowledge, none of these covariance calculation methods have been quantitatively validated for accuracy either experimentally or with Monte Carlo simulation, primarily because their developments aimed at extracting general covariance behavior. Also, we note that none of these methods easily lends itself to accurate and effcient computational evaluation of lesion detectability using some model observers, such as a channelized Hotelling observer, which is discussed next.

In a clinical setting, a common task is the visual detection of lesions from diagnostic CT images. Therefore, an important measure of image quality is lesion detectability. One method to evaluate lesion detectability is through human observer studies (Barrett and Myers, 2004). Because human observer studies can be cumbersome and time-consuming, computerized model observers offer an attractive alternative for the evaluation of lesion detectability. One such model observer is the channelized Hotelling observer (CHO) (Myers and Barrett, 1987; Barrett et al., 1993; Eckstein and Whiting, 1995; Burgess et al., 1997; Abbey et al., 1997; Abbey and Barrett, 2001; Eckstein et al., 2003; Barrett and Myers, 2004; Zhang et al., 2006). It is well known that lesion detectability depends both on noise magnitude and noise correlation (Eckstein and Whiting, 1995; Burgess et al., 1997; Abbey and Barrett, 2001). The CHO incorporates this noise information by using the covariance between all image pixels. The CHO formalism has been applied for purposes of image quality evaluation in many areas of medical imaging. In particular, it has been used to evaluate image reconstruction in nuclear medicine, e.g., (Wollenweber et al., 1999; Bonetto et al., 2000; Zeng and Gullberg, 2002; Kim et al., 2004; Gifford et al., 2005; Gifford et al., 2007), and more recently, in x-ray CT (LaRoque et al., 2007).

In this work, we present a new method for the calculation of covariance for images reconstructed using the classical FBP algorithm for direct reconstruction from fan-beam data. Then we discuss how this method can be used to implement a CHO for the purpose of image quality evaluation. The accuracy of our new covariance method and the resulting CHO are first verified through comparison with results from Monte Carlo noise simulation. Afterwards, to illustrate the utility of our results, we investigate the effects of a bowtie filter and tube current modulation on image noise and lesion detectability.

As a special case, our new covariance method also allows for the calculation of image variance. In this case, we noticed that our method turns out to be similar to equations briefly outlined in a Lawrence Berkeley Laboratory publication (Huesman et al., 1977). Unfortunately, the variance equations in (Huesman et al., 1977) seem to have been forgotten since most of the subsequent research literature on image variance makes no mention of this work.

This paper is organized as follows. In section 2, we review necessary background material. This includes our x-ray data model, our noise model, our models for x-ray tube current modulation and bowtie filters, and also a review of the classical direct fan-beam FBP image reconstruction algorithm. Next, section 3 presents our new covariance method, its simplification for the case of variance calculation, and its application to the evaluation of correlation coeffcient between image pixel values. In section 4, we describe the channelized Hotelling observer methodology and a specific choice for the channel filters. Section 5 contains four sets of experiments using simulated x-ray data. The first set of experiments compare the image covariance and variance estimated from Monte Carlo simulation to the image covariance and variance computed using our new analytical covariance method and its simplification for the case of variance calculation. In the second set of experiments, the performance of a CHO implemented using covariance matrices estimated via Monte Carlo simulation is compared to the performance of a CHO using covariance matrices computed using our new covariance method. The third and fourth sets of experiments evaluate the impacts of a bowtie filter and x-ray tube current modulation on noise and lesion detectability. Finally, section 6 summarizes our results and conclusions.

2. Background

In this section, we describe our notation, necessary background material, and modelling assumptions that are used throughout the paper. In particular, we provide our models for the data acquisition geometry, the x-ray data, the data noise, x-ray tube current modulation, and the bowtie filter. Since this paper studies noise in images reconstructed using the classical direct fan-beam FBP algorithm, we end the section with a description of this algorithm.

2.1. Scanning Geometry and X-ray Data Model

Throughout the paper we represent points in 2-dimensional space using cartesian coordinate vectors. For example, an arbitrary point is given by x = [x, y]T, where the superscript T denotes the transpose operator. The linear attenuation function for the scanned object is designated as μ(x). Technically, μ(x) should be a function of photon energy. However, we neglect this dependence and assume a monochromatic x-ray source. The reconstruction of the linear attenuation function is to be found from x-ray measurements collected in the geometry of 3rd generation CT, which is currently used in most medical CT scanners. This geometry implies fan-beam scanning with a curved detector rotating synchronously with the x-ray source (Kalender, 2006). See Figure 1 for an illustration.

Figure 1
Fan-beam data acquisition geometry

Relative to the scanned object, the source trajectory is a circle given by

a(λ)=[Rocos(λ),Rosin(λ)]Twithλ[λs,λe),
(1)

where Ro is the scanning radius, λ is the polar angle relative to the x-axis, λs is the start angle, and λe is the end angle. By construction, the scanned object is always inside this circle. For any ray emanating from a(λ), let γ be the angle of this ray relative to the line through a(λ) and the origin. Figure 1 depicts the fan-beam geometry with a positive orientation shown for both λ and γ. We denote the ray with angle γ emanating from a(λ) as L(λ, γ).

The number of photons that leave the source assembly and enter into the interrogated object in the direction of L(λ, γ) is represented by Ni(λ, γ). Similarly, No(λ, γ) is the number of photons that exit out of the interrogated object in the direction of L(λ, γ) and hit the detector. In order to model the effects of noise, No(λ, γ) is taken to be a random variable with mean No(λ, γ). The probability model for No is covered in the next section. By definition, the noiseless fan-beam data is

g(λ,γ):=ln(Ni(λ,γ)N¯o(λ,γ)).
(2)

In our x-ray data simulation model, we employ a subsampling of the x-ray tube focal spot and each detector in order to model the blurring that results from the finite size of the focal spot and detector elements, and also to model the shift-variant effect of the x-ray tube anode angle, β, on resolution (Hsieh, 2003, p. 248–252); see Figure 2. Both the x-ray tube focal spot and the detector elements are modeled as planar, rectangular regions in 3-D space. Moreover, a subsampling of each source position is used to model the blurring that results from continuous (uninterrupted) x-ray emission as the x-ray tube rotates (Hsieh, 2003, p. 254).

Figure 2
Subsampling of the x-ray focal spot, the detector element, and the source trajectory; A 2 × 2 subsampling of the focal spot and detector element is shown. The anode angle is β.

Now we are ready to describe more precisely how we calculate No. Represent the x-ray tube focal spot subsamples with indices (i, j), the detector element subsamples with indices (k, l), and the source position subsamples with index m. In addition, let Lijklm be the ray emanating from subsource (i, j) to subdetector (k, l) for subview m. Suppose that there are M total subsample combinations. We compute the mean number of photons out of the interrogated object as

N¯o(λ,γ)=Ni(λ,γ)1Mi,j,k,l,mexp(Lijklmμ(x)dx).
(3)

2.2. Noise Model

For our noise model, we make three simplifying assumptions. First, we assume that the data are statistically uncorrelated for different detectors and for different views. Second, we assume that each detector may be modelled as a photon counting detector. Hence, the number of photons out of the object, No(λ, γ), is assumed to be Poisson distributed with mean No(λ, γ) (Kak and Slaney, 1988). Third, as already mentioned, we assume that the x-ray beam consists of monoenergetic photons.

It is natural to ask how these assumptions relate to the actual physics of the data acquisition and to the data pre-processing. First, the assumption that the data are uncorrelated implies that we ignore crosstalk between detector elements (Hsieh, 2003), as well as the effects of the data correction usually done on CT scanners for scatter (Seibert and Boone, 1988). Second, our assumption of monoenergetic photons means that we ignore the true polyenergetic nature of the x-ray beams in CT scanners. Because beam-hardening is a polyenergetic effect, we are thus neglecting beam-hardening and the beam-hardening correction applied during pre-processing (Hsieh, 2003). Finally, the assumption that the detector may be modelled as a photon counting detector means that we neglect electronic noise, which would contribute a white Gaussian noise component (Whiting et al., 2006). Therefore, the x-ray tube current is assumed to be large enough so that electronic noise may be disregarded.

We use gm(λ, γ) to denote the noisy measured fan-beam data, defined as

gm(λ,γ):=ln(Ni(λ,γ)No(λ,γ)).
(4)

Using equation (2), we may relate the noisy fan-beam data to the noiseless fan-beam data. In particular,

gm(λ,γ)=lnNi(λ,γ)lnNo(λ,γ)
(5)

=lnNi(λ,γ)lnN¯o(λ,γ)+lnN¯o(λ,γ)lnNo(λ,γ)
(6)

=g(λ,γ)+lnN¯o(λ,γ)lnNo(λ,γ)
(7)

=g(λ,γ)+n(λ,γ),
(8)

where n(λ, γ) = ln No(λ, γ) − ln No(λ, γ). Hence, the noisy fan-beam data consists of the noiseless fan-beam data plus non-stationary noise consisting of a deterministic function minus the log of a Poisson random variable. In appendix A, we derive high-order approximations for the mean and the variance of the noisy fan-beam data. These approximations will be used later in the paper.

2.3. X-Ray Tube Current Modulation

Modern CT scanners modulate the x-ray tube current during the scan as a means of reducing x-ray dose to the patient while maintaining image quality (Kalra et al., 2004). The basic idea of tube current modulation (TCM) is to decrease the tube current (and hence, Ni) when the x-ray beam is weakly attenuated and to increase the tube current when the x-ray beam is strongly attenuated. We assume that in the absence of tube current modulation or a bowtie filter (see next section), the number of photons entering the interrogated object is a constant for each ray, i.e. Ni(λ, γ) = Ñi, where Ñi is a constant. The effect of the tube current modulation is to make Ni vary with λ. In other words, Ni(λ,γ)=NiTCM(λ), for some function NiTCM(λ).

The tube current modulation may be modeled by choosing NiTCM(λ) so that the mean number of photons exiting the object along the central ray is constant. However, other (more attractive) options are possible. Following Gies et al. (Gies et al., 1999), we define a modulation parameter α [set membership] [0, 1]. Using the central ray of the λ = 0 view as a reference, the tube current modulation is modeled by

NiTCM(λ)=Nieα(gc(λ)gc(0))wheregc(λ):=g(λ,γ)γ=0,
(9)

with g(λ, γ) given by equations (2) and (3). Observe that if α = 0, then no TCM is applied, and if α = 1.0, then the mean number of photons exiting the object along the central ray is a constant. On the other hand, if α = 0.5 then the mean number of photons exiting the object along the central ray is proportional to the square root of the attenuation factor for this ray, which is generally not constant over all views.

2.4. Bowtie Filter Model

On x-ray CT scanners, the bowtie filter is a shaped piece of material (usually metal) which is placed between the x-ray source and the patient. It is designed to equalize the intensities of the rays hitting the dectector for a given attenuating object (e.g., for a centered water cylinder of given radius) (Hsieh, 2003). Different bowtie filters are required for different scanning field-of-view (FOV) sizes. Typically, CT scanners use one bowtie filter for scans of the head, and another for scans of the body. The primary purpose of the bowtie filter is to decrease the patient radiation dose near the edges of the scanning FOV (Hsieh, 2003). The bowtie filter is also known to help reduce the effects of beam-hardening and to improve detector dynamic range utilization.

The effect of the bowtie filter is to make the number of photons going into the scanned object, Ni, vary with γ. In the presence of tube current modulation and a bowtie filter, we model the number of photons entering the interrogated object as

Ni(λ,γ)=NiTCM(λ)eη(γ)
(10)

where η(γ) is a function which models the effect of the bowtie filter. If no tube current modulation is applied, then we use the constant Ñi in place of NiTCM(λ).

We assume that the bowtie filter is designed for a centered circular water cylinder with linear attenuation function μw(x). Let μb be the linear attenuation coeffcient for the bowtie filter, and let T be its thickness at γ = 0. Then η is defined as

η(γ):=gw(0)gw(γ)+Tμb,
(11)

where gw(γ) is the line integral of μw(x) along the line L(λ, γ) for some fixed λ. Note that because the circular water cylinder is centered at the origin, it is circularly symmetric, and any fixed value of λ will suffice for the definition of gw.

2.5. The Classical Direct Fan-Beam FBP Algorithm

First, we set some notation and definitions. Denote the scanning field-of-view radius by R and define half the fan-angle as γm:= arcsin(R/Ro). Recall that the endpoints of the source trajectory are given by a(λs) and a(λe). For a full-scan, we take [λs, λe) = [0, 2π). In the case of a short-scan, the start angle, λs, is arbitrary and the end angle is λe = λs +π+2γm. Finally, the ray emanating from a(λ) which passes through the point x is given by L(λ, γ*), where

γ(λ,x)=arctan(xsinλ+ycosλRoxcosλysinλ).
(12)

For a circular source trajectory and a curved detector, a classical direct FBP reconstruction from fan-beam data (Kak and Slaney, 1988) estimates μ(x) as

μ(x)=Roλsλe1xa(λ)2γmγmhF(γγ)cosγm(λ,γ)gm(λ,γ)dγdλ.
(13)

In this formula, the outer integral is a backprojection, while the inner integral is a filtering operation defined with

hF(γ):=(γsinγ)2hramp(γ),
(14)

where hramp is an apodized version of the ramp filter kernel, namely,

hramp(γ):=bbφ(σ)σei2πσγdσ.
(15)

In equation (15), b is a cut-off frequency and [var phi](σ) is an apodization window. Both b and [var phi] are used to control noise and finite resolution effects.

Equation (13) also contains a smooth weighting function, m(λ, γ), which handles redundancies in the fan-beam data. This function satisfies

m(λ,γ)+m(λ+π2γ,γ)=1foranyλ[0,2π)andγ[γm,γm].
(16)

For the case of a full-scan, the choice m(λ, γ) = 1/2 is used in this paper. For the case of a short-scan, we use the weighting function

m(λ,γ)=c(λ)c(λ)+c(λ+π2γ),
(17)

where

c(λ)={cos2(π(λλsd)/(2d))ifλsλ<λs+d1ifλs+dλ<λedcos2(π(λλe+d)/(2d))ifλedλ<λe,
(18)

with d being an angular interval over which c(λ) smoothly drops from 1 to 0 (Noo et al., 2002). If d is small, then the weighting function is similar to that of Chen et al. (Chen et al., 2006). On the other hand, if d is large, then the weighting function is similar to Parker weighting (Parker, 1982).

We assume that λ and γ are sampled uniformly over their respective ranges at positions

λn=(n1)Δλ+λsn=1,2,,N
(19)

and

γj=(j+3/4)Δγj=J,J+1,,J1
(20)

where

Δλ=(λeλs)/(N1)andΔγ=2γm/(2J).
(21)

In discrete form, the filtered data may be expressed as

gF(λn,γi)=Δγj=JJ1hF(γiγj)cosγjm(λn,γj)gm(λn,γj),
(22)

with hF defined as in equations (14) and (15) with b = 1/(2Δγ) for the kernel cut-off frequency. Also, the backprojection integral in equation (13) may be written as

μ(x)=RoΔλn=1N1xa(λn)2Ix{gF(λn,γi)}.
(23)

In this equation, x2110x is a linear operator which interpolates the filtered data between the discrete locations γi to find the value at γ*(λn, x). Henceforth, we assume that the operator x2110x is implemented as linear interpolation. Rewriting equation (23) using linear interpolation for x2110x yields

μ(x)=RoΔλn=1N1xa(λn)2{[1un]gF(λn,γk)+ungF(λn,γk+1)},
(24)

with

k:=γ(λn,x)Δγ34andun:=γ(λn,x)Δγ34k,
(25)

where [left floor]z[right floor] is the floor function, which gives the largest integer that is smaller than z.

3. Image Covariance

In this section, we introduce a new fast and accurate method for calculating the covariance between image pixel values reconstructed using the classical direct fan-beam FBP formula. Specifically, we explain how the image covariance relative to a fixed pixel location may be computed in a FBP format. In addition, we examine the special case of this formula for variance calculation and show how it too may be implemented in a FBP framework.

3.1. Covariance Formula

Let r and s be two arbitrary indices, and denote the covariance of two filtered data terms within the projection at λn by

Cn(r,s):=Cov[gF(λn,γr),gF(λn,γs)].
(26)

In appendix B, we show that the covariance between two pixel values, [mu](x) and [mu](x0), can be computed from values of Cn(r, s) as follows

Cov[μ(x),μ(x0)]=(RoΔλ)2n=1N1xa(λn)21x0a(λn)2×{[1un][1vn]Cn(k,l)+un[1vn]Cn(k+1,l)+[1un]vnCn(k,l+1)+unvnCn(k+1,l+1)}.
(27)

In the above equation, the indices are given by

k:=γ(λn,x)Δγ34andl:=γ(λn,x0)Δγ34,
(28)

and the weights are defined as

un:=γ(λn,x)Δγ34kandvn:=γ(λn,x0)Δγ34l.
(29)

The covariance terms in equation (27) may be written in the general form

Cn(r,s)=(Δγ)2j=JJ1hF(γrγj)(hF(γsγj)cos2γj(m(λn,γj))2Var[gm(λn,γj)]),
(30)

where r [set membership] {k, k + 1} and s [set membership] {l, l + 1}. When s [set membership] {l, l + 1} is fixed, equation (30) may be computed effciently as a convolution by pre-weighting Var[gm] and then filtering the outcome with hF. A high-order approximation to the fan-beam data variance, Var[gm(λn, γj)], may be computed using equation (A.7) together with Table A2, which are both provided in appendix A.

Table A2
Coefficients for approximation (A.7) of Var[gm]

Equations (27)(30) constitute a new FBP method to calculate a covariance image relative to a fixed point x0. For each view, we apply equation (30) twice using the fast Fourier transform, first using the weight hF (γlγj) cos2 γj (m(λn, γj))2, and next using the weight hF (γl+1γj) cos2 γj (m(λn, γj))2. Then all values needed for equation (27) are available, and it can be implemented as a backprojection step.

3.2. Variance Formula

We denote the variance of a filtered data term for any index r in the λn projection as

Vn(r):=Var[gF(λn,γr)].
(31)

In the special case of variance calculation, equation (27) reduces to

Var[μ(x)]=(RoΔλ)2n=1N1xa(λn)4×{(1un)2Vn(k)+2(1un)unCn(k,k+1)+un2Vn(k+1)},
(32)

where k and un are as defined in equations (28) and (29).

Equation (30) may be used to express the desired variance and covariance of the filtered data as

Vn(k)=(Δγ)2j=JJ1hF2(γkγj)cos2γj(m(λn,γj))2Var[gm(λn,γj)],
(33)
Cn(k,k+1)=(Δγ)2j=JJ1hF(γkγj)hF(γk+1γj)(cosγjm(λn,γj))2Var[gm(λn,γj)].
(34)

Equation (33) is a convolution which can be computed effciently. As it is written, equation (34) is not a convolution. However, there is a way to express equation (34) as a convolution. First, we define

hG(γ):=hF(γ)hF(γ+Δγ).
(35)

It follows that

hF(γkγj)hF(γk+1γj)=hF(γkγj)hF(γkγj+Δγ)=hG(γkγj).
(36)

Thus, equation (34) may be rewritten in the form

Cn(k,k+1)=(Δγ)2j=JJ1hG(γkγj)cos2γj(m(λn,γj))2Var[gm(λn,γj)],
(37)

which is a convolution.

The image variance may be calculated by first filtering the data variance for each view according to equations (33) and (37) and then backprojecting the results using equation (32). Thus, equations (32), (33), and (37) represent a FBP method for calculating image variance.

As discussed in the introduction, our equations (32)(34) for image variance are similar to equations outlined in (Huesman et al., 1977, p. 29). However, the observation that equation (34) may be implemented effciently as a convolution, as described by equations (35)(37), was not made in (Huesman et al., 1977).

More recently, Pan and Yu suggested another formula for the calculation of image variance for the special case of stationary data noise (Pan and Yu, 2003). Assuming stationary data noise, one may show that their formula is mathematically equivalent to equations (32)(34). However, the Pan and Yu formula does not lend itself to a FBP implementation. A similar formula to that of Pan and Yu was also found by Herman (Herman, 1980).

3.3. Correlation Coeffcient

Each pixel value in the reconstruction, [mu](x), is a linear combination of the fan-beam data, gm(λ, γ), over all measured angles λn and γj. For a typical CT scanner with N = 1160 views per turn and 2J = 672 detectors, there are a total of 1160 × 672 = 779, 520 measurements for a full scan. Therefore, by the multivariate central limit theorem (Anderson, 2003, p. 86), it is reasonable to expect that the image noise is well-approximated by a multivariate normal distribution, irrespective of the underlying probability distribution for the measured fan-beam data.

If the image noise is multivariate normal, then the correlation coeffcient is meaningful as a measure of the statistical dependence between two pixel values (Anderson, 2003, p. 20–22). The correlation coeffcient between the reconstruction values at two locations x and x0 is defined as

ρ(μ(x),μ(x0)):=Cov[μ(x),μ(x0)]Var[μ(x)]Var[μ(x0)].
(38)

Recall that the correlation coeffcient is normalized so that −1 ≤ ρ ≤ 1. Using the formulas from sections 3.1 and 3.2, correlation coeffcient images may be calculated relative to a fixed point x0. Examples of such images are provided in section 5.

4. Lesion Detectability

As mentioned in section 1, the channelized Hotelling observer (CHO) provides a methodology for the evaluation of image quality for a lesion detection task (Barrett et al., 1993; Barrett and Myers, 2004). In this section, we will briefly review the definition of the CHO and describe how the covariance method of section 3.1 may be used for the implementation of a CHO. In addition, we describe a choice of channel filters that is used later for CHO experiments.

4.1. Channelized Hotelling Observers

We consider a signal-known-exactly/background-known-exactly (SKE/BKE) binary classification task in which the goal is to decide whether or not a lesion is present given a noisy image realization. In this task, the only source of noise comes from the measurement of the fan-beam data. The noiseless lesion and background are assumed to be known exactly. The set of noisy images for which no lesion is present is called class 1. Similarly, the set of noisy images for which the lesion is present is class 2.

The CHO produces a statistic, t, which is compared to a threshold tc to decide if the image belongs to class 1 or class 2. In particular, if t > tc, the CHO decides that a lesion is present and that the image belongs to class 2. Otherwise, the CHO decides that no lesion is present and that the image belongs to class 1. In the following, we describe how to find the CHO statistic, t.

Assume that the noisy image realization is represented as a Q × 1 column vector, p. The first step of the CHO computation is to filter the image with P (where P [double less-than sign] Q) channel filters so that the dimensionality of the data set is reduced from Q to P. A choice for these channel filters is described in section 4.2. The impulse response of each channel filter is represented as a column of a Q × P channel matrix, denoted by U. This channel matrix is applied to each image realization to get a P × 1 channel output vector, v = UTp. The CHO statistic is then found by taking the inner product of v with a P × 1 vector wHot called the channelized Hotelling discriminant, i.e., t:=wHotTv.

A number of intermediate quantities are needed to define the channelized Hotelling discriminant. We start by denoting the respective means of p over class 1 and class 2 as [p with macron]1 and [p with macron]2. Also, we write the covariance matrices of p over class 1 and class 2 as Kp1 and Kp2, respectively. The means of the channel outputs are then expressed as v1:= UT[p with macron]1, v2:= UT[p with macron]2. In addition, the channel output covariance matrices for classes 1 and 2 can be obtained as Kv1 = UTKp1U and Kv2 = UTKp2U.

The channelized Hotelling discriminant is defined from two quantities, a P × 1 vector Δv and a P × P matrix Sv, according to

wHot:=Sv1Δv¯.
(39)

The vector Δv is the difference of channel output means, i.e., Δv:= v2v1, and the matrix Sv is the intraclass channel scatter matrix defined as

Sv:=12(Kv1+Kv2).
(40)

A common way to describe the performance of an observer for a binary classification task is through a receiver operating characteristic (ROC) curve, or more specifically, the area under this curve, often written as the AUC (Barrett and Myers, 2004). The ROC curve plots the true positive fraction (TPF) versus the false positive fraction (FPF) for various values of the decision threshold, tc. The AUC is in the range [0.5, 1], where higher values of the AUC indicate better CHO performance for the lesion detection task. If the AUC = 1, then the CHO always detects the lesion, regardless of the decision threshold. On the other hand, if the AUC = 0.5, then for any threshold, TPF = FPF, and the CHO performs no better than arbitrary guessing. If the CHO statistic, t, is normally distributed for each class, then the AUC may be computed as (Barrett and Myers, 2004)

AUC=12[1+erf(12Δv¯TSv1Δv¯)]witherf(z)=2π0zey2dy.
(41)

As discussed in section 3.3, the image noise may be assumed to have a multivariate normal distribution. Because the CHO statistic is calculated as a linear combination of image pixel values, we may thus assume that it is normally distributed for each class. Alternatively, we may invoke the univariate central limit theorem to conclude that the CHO statistic is normally distributed for each class (Barrett and Myers, 2004). Hence, the CHO performance may be accurately estimated with equation (41).

4.2. Channel Selection

In section 4.1, the channelized Hotelling observer methodology was introduced. A key component for the implemention of a CHO is the selection of the channel filters. There are many published choices for the channels, e.g., (Myers and Barrett, 1987; Daly, 1993; Burgess et al., 1997; Abbey and Barrett, 2001; Eckstein et al., 2003; Zhang et al., 2006). One choice for the impulse response of the channels is the Gabor functions, as described in (Eckstein et al., 2003; Zhang et al., 2006); their expression is

Ch(x,y)=exp[4(ln2)((xx0)2+(yy0)2)/ws2]×cos[2πfc((xx0)cosθ+(yy0)sinθ)+ξ].
(42)

These functions consist of a Gaussian modulated by a cosine function. In equation (42), fc denotes the center frequency of the channel, ws is the spatial width of the channel, the point (x0, y0) is the center of the channel, θ is the channel orientation, and ξ is a phase factor. The spatial width of the channel, ws, is related to wf, the full-width-at-half-maximum of the Fourier transform of the channel, through ws = 4(ln 2)/(πwf). The channel center (x0, y0) should always coincide with the lesion center. A choice for the channel parameters will be given in section 5.

5. Computational Evaluations

This section contains several computational examples which illustrate both the accuracy and the utility of our new covariance method and its application to a CHO. For the first set of examples, results found using Monte Carlo simulation are compared to results found using the analytic covariance method. In particular, values of the covariance, the variance, and the AUC for a lesion detection task are compared. The second set of examples illustrates a possible application of our new covariance method. In these examples, we consider how a bowtie filter and x-ray tube current modulation can impact noise and lesion detectability in CT images.

5.1. Computational Parameters

The scanning parameters, which are used for all of the fan-beam reconstructions hereafter, are listed in Table 1. They are representative of current medical CT scanners. In addition, each of the examples in section 5 uses a rectangular window for the ramp filter kernel, i.e. [var phi] [equivalent] 1 in equation (15).

Table 1
Fan-beam scanning parameters

As described in section 2.1, the x-ray data simulation uses a subsampling of each detector, the x-ray tube focal spot, and each source position. All of the examples in section 5 use a 3 × 3 subsampling for both the detector and the x-ray tube focal spot and use 5 subviews for each view angle. Hence, following equation (3), 405 ray integrals are exponentially averaged to compute each value of No.

Numerical evaluation of the analytical expressions for image covariance, image variance, and AUC from sections 3 and 4 requires values for the mean and variance of the fan-beam data, gm. To estimate these values, we use approximations (A.4) and (A.7), with all terms included.

For the CHO implementation, we use the Gabor channels that were decribed earlier by equation (42). We considered four channel passbands, given by [1/64, 1/32], [1/32, 1/16], [1/16, 1/8], and [1/8, 1/4] cycles/pixel, respectively; this choice of passbands is the same as Zeng and Gullberg (Zeng and Gullberg, 2002) and is similar to the passbands used in (Myers and Barrett, 1987; Wollenweber et al., 1999; Kim et al., 2004; Gifford et al., 2007). The center frequencies for the passbands are thus fc = 3/128, 3/64, 3/32, and 3/16 cycles/pixel, respectively, and their spatial channel widths are given by ws = 56.48, 28.24, 14.12, and 7.06 pixels, respectively. In addition, following Eckstein et al. (Eckstein et al., 2003), we use five orientations and two phases for each passband, so that the total number of channels is 40 (4 frequency bands × 5 orientations × 2 phases). The five orientations are θ = 0, 2π/5, 4π/5, 6π/5, and 8π/5 radians and the two phases are ξ = 0 and ξ = π/2 (see equation (42)).

Throughout sections 5.2–5.4, we compute the AUC for a CHO applied to a lesion detection task for a 96 × 96 image with grid spacing Δx = Δy = 0.045 cm, centered on the lesion location. The lesion is 3 mm wide with a relative contrast of 20 HU. For this image size, the spatial impulse response of the lowest frequency CHO channel filter is mostly contained within the image (see equation (42)). More specifically, with fc = 3/128 cycles/pixel and ws = 56.48 pixels, 99.29% of the total energy of the Gabor channel impulse response is contained within the image. Therefore, using an image size bigger than 96 × 96 should not make a large difference in the CHO performance.

5.2. Comparison of Monte Carlo and Analytical Estimates

The examples in this section all use a short-scan source trajectory with λ [set membership] [0, π+2γm), tube current modulation with α = 1.0, and a bowtie filter designed for a circular water cylinder of diameter 24 cm. The scanned object is the FORBILD head phantom (with no ears); see Figure 3. In addition, a smoothing range of d = 8.7267 × 10−2 radians is used in the short-scan weighting function (see equations (17) and (18)).

Figure 3
The FORBILD head phantom (with no ears) used for the examples in section 5.2. For the third example in section 5.2, a lesion is inserted at position (x, y) = (4, −7) cm. The lesion, indicated by the white arrow, is 3 mm wide and has a contrast ...

Each example in this section involves Monte Carlo estimates of either image covariance or image variance. To find a Monte Carlo estimate based on n realizations, a set of n noisy x-ray data realizations is first computed. Next, each noisy data realization is reconstructed to form a noisy image realization. Finally, either the image covariance or image variance is estimated from the set of n noisy image realizations using the standard unbiased estimators (Anderson, 2003).

The first example compares covariance images estimated with Monte Carlo simulation to covariance images computed with the new analytical method of section 3.1. In this example, we compute the covariance image relative to the position (4, −7) cm, in the FORBILD head phantom. The covariance images are computed on a 32 × 32 grid centered at (4, −7) cm with a grid spacing of Δx = Δy = 0.01 cm, and with a photon level of Ñi = 30, 000. Figure 4 (Left) compares the analytical covariance estimate to Monte Carlo estimates found by progressively increasing the number of noise realizations from 50 to 100, 200, 1000, and 10, 000. We see that as the number of noise realizations is increased, the shape of the covariance image becomes more similar to the analytical result. This visual assessment is supported by the plot in Figure 4 (Right), which shows the relative root-mean-square di3erence between the Monte Carlo and analytical results, and by the plot in Figure 5, which shows the profile along the line x = 4 cm for some of the covariance images. For 10, 000 noise realizations, the relative root-mean-square difference is 3.67%. We investigated if this remaining difference is due to the number of noise realizations or to the accuracy of expansion (A.7) that is used for the data variance in the analytical computation of the covariance image. If only the first term of expansion (A.7) is used instead of all terms, we find that the relative root-mean-square difference between the analytical result and the Monte Carlo estimate with 10, 000 noise realizations is 3.74%. Therefore, in this experiment, the last six terms in (A.7) contribute little in covariance accuracy, and the discrepancy of 3.67% that remains between the Monte Carlo result with 10, 000 noise realizations and the analytical result appears to be primarily due to the finite number of noise realizations.

Figure 4
(Left) Covariance images relative to the lesion location, (x, y) = (4, −7) cm, in the FORBILD head phantom of Figure 3 with the lesion removed. The absolute value of the covariance image is displayed instead of the image itself, to enhance the ...
Figure 5
Vertical profile through the center of covariance images shown in Figure 4. The dashed lines with squares, circles, and triangles, are for the Monte Carlo covariance images obtained with the first 100, 200 and 10,000 noise realizations, respectively. ...

The second example compares variance images estimated with Monte Carlo simulation to variance images found using the analytical method of section 3.2. All images were computed on a 128 × 128 grid with Δx = Δy = 0.2 cm and with a photon level of Ñi = 30, 000. Figure 6 (Left) compares the analytical variance estimate to Monte Carlo estimates found by progressively increasing the number of noise realizations from 50 to 100, 200, 1000 and 10 000. As the number of noise realizations increases, the Monte Carlo estimates appear less grainy and become more similar to the analytically computed variance image. The relative root-mean-square difference between the Monte Carlo and analytical results is plotted in figure 6 (right). For 10 000 noise realizations, the relative root-mean-square difference is 1.5%.

Figure 6
Left: variance images for the FORBILD head phantom of figure 3 with the lesion removed. Top row: from left to right, we see the Monte Carlo estimates obtained from 50, 100 and 200 noise realizations. Bottom row: from left to right, the Monte Carlo estimates ...

The third example of this section implements the CHO of section 5.1 for a lesion detection task. The lesion detection task is with a photon level of Ñi = 100 000 and with a lesion placed in the phantom at location (4, −7) cm (see figure 3). In this example, the lesion detectability is measured by computing the AUC (see equation (41)). The AUC value estimated by using Monte Carlo estimates of the necessary covariance matrices is compared to the AUC value found using the analytical covariance method of section 3.1. The Monte Carlo estimates progressively increase the number of noise realizations from 50 to 100, 200, 1000, and 10, 000 noise realizations, respectively. The left plot in Figure 7 shows the AUC estimates as a function of the number of realizations, with a reference line for the analytical estimate of the AUC. The right plot in Figure 7 illustrates the relative difference between the Monte Carlo AUC estimates and the analytical estimate. From both plots, we see that as the number of noise realizations increases, the Monte Carlo AUC estimates become progressively closer to the analytical AUC estimate. For 10, 000 noise realizations, the relative difference between the Monte Carlo and analytical estimates is 0.27%.

Figure 7
Monte Carlo estimates of the AUC value at (x, y) = (4, −7) cm in the Forbild head phantom of Figure 3., with a progressive increase in the number of noise realizations. (Left) The horizontal dashed line gives the analytical estimate of the AUC ...

The new analytical covariance method was found to yield a significant savings in computational time compared to the Monte Carlo simulations. On one hand, for the cases of 200 and 1000 noise realizations, the Monte Carlo AUC estimates required approximately 2 and 11.5 hours of computation on a personal computer to yield estimates differing from the analytical estimate by 10% and 3.7%, respectively; see Figure 7 (Right). On the other hand, the analytical AUC estimate required only 46 minutes of computation on the same personal computer.

It is important to understand that the Monte Carlo estimates in this section are either random variables (in the case of AUC) or spatially varying random processes (in the cases of covariance and variance). Therefore, the Monte Carlo estimates each have some associated standard deviation (which is spatially varying for the cases of covariance and variance). Due to the excessive computational effort required, no attempt was made to estimate the standard deviations for the examples above. Therefore, the Monte Carlo results should only be taken to indicate general trends.

5.3. Impact of a Bowtie Filter

In this section, we evaluate the effects of a bowtie filter on lesion detectability and noise in a circular water cylinder of diameter 38 cm scanned over a full-scan source trajectory, i.e., λ [set membership] [0, 2π). A bowtie filter designed to perfectly equalize the attenuation of the object is considered. When no bowtie filter is present, the photon level is Ñi = 200, 000. In the presence of a bowtie filter, the photon level is increased to Ñi = 200, 000 exp(T μb), to compensate for the finite bowtie filter thickness at γ = 0.

The left plot in Figure 8 shows the AUC for various lesion locations on the x-axis, with and without a bowtie filter. The right plot in Figure 8 illustrates the variance profiles with and without a bowtie filter for points on the x-axis. Without a bowtie filter, the AUC increases with x, but this is not surprising if one remembers that exposure increases with x. A more important observation from Figure 8 is that with a bowtie filter, AUC degrades only slightly with increasing x, while the variance increases rather quickly. We see that the bowtie filter makes the AUC nearly uniform.

Figure 8
AUC value and variance profile for detection of a lesion placed at various locations on the x-axis; the lesion is embedded in a centered, circular water phantom of diameter 38 cm. The dashed curve corresponds to data acquisition with no bowtie filter, ...

We find it interesting to look at correlation coefficient images for this example. Images of correlation coefficient relative to three lesion locations on the x-axis are given in Figure 9. More specifically, correlation coefficient images centered on lesions at x = 0 cm, x = 9 cm and x = 18 cm on the x-axis are shown. Each image is 64 × 64 and has a grid spacing of Δx = Δy = 0.01 cm. The top row of images are for the case of no bowtie filter and the bottom row of images are in the presence of a bowtie filter. See how the presence of a bowtie filter makes the correlation coefficient support at x = 18 cm narrower and more symmetrical, with strongly reduced sidelobes. However, it is not easy to make inferences about the AUC from these images.

Figure 9
Images of correlation coefficient at three lesion locations along the x-axis, for the imaging scenario of Figure 8. The top row is for data acquisition with no bowtie filter, while the bottom row corresponds to using a bowtie-filter. In each row, the ...

5.4. Impact of Tube Current Modulation

The next example explores the influence of x-ray tube current modulation on lesion detectability in an elliptic water cylinder scanned with a full-scan source trajectory. The elliptic water cylinder has a diameter of 38 cm in the x-direction and a diameter of 25 cm in the y-direction. A bowtie filter designed for a circular water cylinder of diameter 38 cm is assumed to be present. The TCM effect is examined for two choices of α, namely, α = 1.0 and α = 0.5. In the presence of TCM, the photon level, Ñi, is modified so that the total exposure for the reconstruction at the origin is the same as without TCM. The total exposure at the origin is approximated as the sum of the number of photons that enter the object along each measured line through the origin. Equalizing the total exposure at the origin leads to photon levels of Ñi= (2.80)(200, 000) exp(T μb) and Ñi = (1.82)(200, 000) exp(T μb) for TCM modulation parameters of α = 1.0 and α = 0.5, respectively.

Figure 10 contains plots of the AUC and image variance along the x-axis. We see that when the phantom becomes elliptical and no TCM is applied, the AUC is not as uniform as it was for the case of a circular cylinder in the presence of a bowtie filter. This is not surprising, since the bowtie filter is designed for a circular water cylinder. However, the AUC non-uniformity is greater than one might have expected.

Figure 10
AUC value and variance profile for detection of a lesion placed at various locations on the x-axis. The lesion is embedded in a centered, ellipse-shaped water phantom (width: 38 cm, height: 25 cm). Data acquisition is performed with a bowtie filter. The ...

Near the origin, TCM with α = 1.0 gives comparable variance to no TCM, as observed by (Gies et al., 1999). This also remains true at large x-values (|x| < 15 cm). On the other hand, lesion detectability is strongly degraded for the case of TCM with α = 1.0. For TCM with α = 0.5, we see that the image variance decreases at the origin, as predicted by (Gies et al., 1999). Moreover, we observe that this remains true away from the origin. The decrease in variance comes with gain in lesion detectability compared to TCM with α = 1.0. However, TCM with α = 0.5 still suffers from lower lesion detectability than the no TCM case, emphasizing that lower image variance does not necessarily guarantee better lesion detectability.

6. Discussion and Conclusions

We presented a new method to analytically determine the covariance between pixel values in reconstructions performed using the classical direct fan-beam FBP algorithm. This new method yields the covariance image relative to any fixed point in an efficient FBP format. Moreover, as a special case, it can deliver variance images. In this case, we noticed that our approach is similar to formulae outlined in (Huesman et al., 1977), which seem to have gone forgotten in the research literature. The utility of our covariance method was demonstrated through the implementation of a channelized Hotelling observer. More specifically, we showed that the results of this paper provide a means to analytically calculate CHO performance metrics, such as the AUC. All of our results were validated by comparison with a Monte Carlo approach. As the number of noise realizations was increased, Monte Carlo estimates of the image covariance, the image variance, and the AUC were found to approach results found using our analytical methods.

As an exciting application of our new covariance method, we investigated the effects of a bowtie filter and x-ray tube current modulation on image noise and lesion detectability for both circular and elliptic water cylinders. These investigations highlighted significant differences between variance and detectability in CT. Specifically that (i) large increases in the variance do not always correspond to large decreases in the AUC (see Figure 8), (ii) methods with a similar variance can have very different detectability properties (see Figure 10), and (iii) compared to no TCM, TCM with a modulation parameter of α = 0.5 results in lower variance as shown in (Gies et al., 1999), but apparently also in lower detectability. The first two observations are not surprising, since in addition to image variance, image correlation is known to have an effect on lesion detectability (Eckstein and Whiting, 1995; Burgess et al., 1997; Abbey and Barrett, 2001). However, we find it illuminating to see these expected observations realized in the context of CT experiments. Regarding the third observation, it is important to note that all of our results are dependent on the reconstruction algorithm that we selected and on the data model. Changes in the reconstruction algorithm and the data model could significantly affect our results, in particular the third observation. To clarify this issue, we are currently investigating refinements of our covariance method using real CT data. In any case, our experiments demonstrate that the ability to accurately evaluate the effects of different dose reduction schemes on noise and lesion detectability is important for the design of CT systems.

Another important question, which was not tackled in this paper, is whether a CHO with our 40 directional Gabor channels reflects human observer performance in x-ray CT. Although a CHO with Gabor channels has been found to indicate trends in human observer performance for image reconstruction in SPECT (Wollenweber et al., 1999), this question does not seem to have been investigated in the context of x-ray CT. Investigation of this issue is another topic of high interest to us.

Acknowledgments

This work was partially supported by NIH grant R01 EB000627 and by a generous grant from the Ben B. and Iris M. Margolis Foundation. Its contents are solely the responsibility of the authors.

Appendix A

In this appendix, we derive high-order approximations for the mean and variance of the noisy fan-beam data, gm. Let mj = E[Xj] denote the jth raw moment of a Poisson random variable X with mean m1. Riordan (Riordan, 1937) showed that the raw moments for the Poisson distribution satisfy

mj+1=m1mj+m1ddm1mjj=1,2,3,
(A.1)

Using the above recurrence relation, we can easily find the raw moments to any order.

In the following, the functional dependence of gm, No, and No on λ and γ is implicit. Recall that No is Poisson distributed with mean No. We assume that No is large enough so that the following manipulations are justified (Rockmore and Macovski, 1977). First, we use equation (6) to rewrite gm as

gm=lnNilnN¯oln(1+(NoN¯o)1)
(A.2)
lnNilnN¯on=1Q(1)n+1n[(NoN¯o)1]n.
(A.3)

We expand this last approximation in powers of No, take the expectation, and then substitute the first Q raw Poisson moments found using equation (A.1). This yields a Qth-order estimate of E[gm], denoted by eQ.

By finding progressively higher-order estimates for E[gm] and collecting the terms in powers of 1/No, one will see that only the first few terms in the expansions stay the same. This indicates that not all of the terms found for a particular estimate are the same as for the true expansion of E[gm]. This issue can be handled by computing two very high-order estimates and then discarding all terms which change between these two estimates. In particular, with the aid of the symbolic algebra capabilities of the Mathematica software package, we obtained mean estimates of orders 19 and 20. By comparing the terms, we found that the first seven terms were the same for the two high-order expansions. Discarding the remaining terms, the mean of the fan-beam data may thus be approximated by

E[gm]ln(NiN¯o)+a1N¯o+a2N¯o2+a3N¯o3+a4N¯o4+a5N¯o5+a6N¯o6+a7N¯o7,
(A.4)

where the values of the coefficients are given in Table A1.

Table A1
Coefficients for approximation (A.4) of E[gm]

Now we focus on the variance of gm. Using a Rth-order Taylor approximation in equation (A.3) and subtracting eQ from both sides gives

gmeQlnNilnN¯on=1R(1)n+1n[(NoN¯o)1]neQ.
(A.5)

We approximate the variance as

Var[gm]E[(gmeq)2].
(A.6)

Next, we substitute (A.5) into (A.6), expand in powers of No, and then use the first R raw Poisson moments to find an approximation to Var[gm]. We call such an expression the variance estimate of order (Q, R).

As we did above for the mean, we found two very high-order estimates of the variance. After expanding these estimates in powers of 1/No, we discarded all terms which changed between the two very high-order estimates. In particular, with the aid of Mathematica, we obtained variance estimates of orders (38, 19) and (40, 20). By comparing the terms, we found that the first seven terms were the same for the two high-order expansions. Discarding the remaining terms, the variance of the fan-beam data may thus be approximated by

Var[gm]b1N¯o+b2N¯o2+b3N¯o3+b4N¯o4+b5N¯o5+b6N¯o6+b7N¯o7,
(A.7)

where the values of the coefficients are given in Table A2.

For both the mean and variance, the first three coefficients found by our procedure agree with the three term expansions found by (Rockmore and Macovski, 1977). Note that in the approximations (A.4) and (A.7), both No and Ni are functions of λ and γ. Consequently, E[gm] and Var[gm] are generally functions of λ and γ.

Appendix B

In this appendix, we derive equations (27) and (30) for the covariance between two image pixel values. As mentioned in section 2.2, we assume that the measured data are uncorrelated for different rays, i.e.

Cov[gm(λm,γi),gm(λn,γj)]=δmnδijVar[gm(λn,γj)],
(B.1)

where the Kronecker delta is defined as usual, namely, δmn = 1 if m = n and δmn = 0 if mn. Since the ramp filtering does not combine data from different views, the filtered data is uncorrelated between views. Hence,

Cov[gF(λm,γi),gF(λn,γj)]=δmnCov[gF(λn,γi),gF(λn,γj)].
(B.2)

Taking the expectation of equation (24) yields

E[μ(x)]=RoΔλn=1N1xa(λn)2{(1un)E[gF(λn,γk)]+unE[gF(λn,γk+1)]}.
(B.3)

Define gF(λn, γk):= gF (λn, γk) − E[gF (λn, γk)]. Also, recall the definitions of the interpolation indices and weights given in equations (28) and (29). The covariance of the image pixel values [mu](x) and [mu](x0) is

Cov[μ(x),μ(x0)]=E[m=1NRoΔλxa(λm)2{(1um)gF(λm,γk)+umgF(λm,γk+1)}×n=1NRoΔλx0a(λn)2{(1vn)gF(λn,γl)+vngF(λn,γl+1)}].
(B.4)

Using the linearity of the Expectation operator and rearranging yields

Cov[μ(x),μ(x0)]=(RoΔλ)2m=1N1xa(λm)2n=1N1x0a(λm)2×{(1um)(1vn)Cov[gF(λm,γk),gF(λn,γl)]+um(1vn)Cov[gF(λm,γk+1),gF(λn,γl)]+(1um)vnCov[gF(λm,γk),gF(λn,γl+1)]+umvnCov[gF(λm,γk+1),gF(λn,γl+1)]}.
(B.5)

Hence, applying equation (B.2) to the above result and using the definition of Cn given in equation (26), we arrive at equation (27).

Next, we derive equation (30). As a special case of equation (B.1), the data are uncorrelated within each view, i.e.

Cov[gm(λn,γi),gm(λn,γj)]=δijVar[gm(λn,γj)].
(B.6)

Taking the expectation of the expression for the filtered data, equation (22), yields

E[gF(λn,γr)]=Δγj=JJ1hF(γrγj)cosγjm(λn,γj)E[gm(λn,γj)].
(B.7)

Define gm(λn, γr):= gm(λn, γr) − E[gm(λn, γr)]. The covariance of two filtered data terms is

Cn(r,s):=Cov[gF(λn,γr),gF(λn,γs)]=E[Δγi=JJ1hF(γrγi)cosγim(λn,γi)gm(λn,γi)×Δγj=JJ1hF(γsγj)cosγjm(λn,γj)gm(λn,γi)].
(B.8)

Using the linearity of the expectation operator, equation (B.6), and then simplifying yields equation (30).

References

  • Abbey CK, Barrett HH. Human- and model-observer performance in ramp-spectrum noise: Effects of regularization and object variability. J Opt Soc Amer A. 2001;18(3):473–488. [PMC free article] [PubMed]
  • Abbey CK, Barrett HH, Eckstein MP. Practical issues and methodology in assessment of image quality using model observers. Proc of SPIE. 1997;3032:182–194.
  • Alvarez RE, Stonestrom JP. Optimal processing of computed tomography images using experimentally measured noise properties. Journ Comp Assist Tom. 1979;3(1):77–84. [PubMed]
  • Anderson T. An Introduction to Multivariate Statistical Analysis. 3. John Wiley & Sons; 2003.
  • Barrett HH, Gordon S, Hershel R. Statistical limitations in transaxial tomography. Computers in Biology and Medicine. 1976;6(4):307–323. [PubMed]
  • Barrett HH, Myers KJ. Foundations of Image Science. John Wiley & Sons, Inc.; 2004.
  • Barrett HH, Yao J, Rolland JP, Myers KJ. Model observers for assessment of image quality. Proc Nat Academy Sci of the USA. 1993;90(21):9758–9765. [PubMed]
  • Bennett K, Byer R. Optical tomography: Experimental verification of noise theory. Optics Letters. 1984;9(7):270–272. [PubMed]
  • Bennett KE, Byer RL. Fan-beam tomography noise theory. J Opt Soc Amer A. 1986;3(5):624–633.
  • Bennett KE, Faris GW, Byer RL. Experimental optical fan beam tomography. Applied Optics. 1984;23(16):2678–2685. [PubMed]
  • Bonetto P, Qi J, Leahy RM. Covariance approximation for fast and accurate computation of channelized hotelling observer statistics. IEEE Trans Nuc Sci. 2000;47(4):1567–1572.
  • Brooks RA, Di Chiro G. Statistical limitations in x-ray reconstructive tomography. Med Phys. 1976;3(4):237–240. [PubMed]
  • Burgess AE, Li X, Abbey CK. Visual signal detectability with two noise components: Anomalous masking effects. J Opt Soc Amer A. 1997;14(9):2420–2442. [PubMed]
  • Chen GH, Tokalkanahalli R, Zhuang T, Nett BE, Hsieh J. Development and evaluation of an exact fan-beam reconstruction algorithm using an equal weighting scheme via locally compensated filtered backprojection (LCFBP) Med Phys. 2006;33(2):475–481. [PubMed]
  • Chesler DA, Riederer SJ, Pelc NJ. Noise due to photon counting statistics in computed x-ray tomography. Journ Comp Assist Tom. 1977;1(1):64–74. [PubMed]
  • Daly S. The visible differences predictor: An algorithm for the assessment of image fidelity. In: Watson AB, editor. Digital Images and Human Vision. MIT Press; Cambridge, MA: 1993. pp. 179–206.
  • Davis GR. The effect of linear interpolation of the filtered projections on image noise in x-ray computed tomography. Journ X-Ray Sci Tech. 1994;4(3):191–199. [PubMed]
  • Duerinckx AJ, Macovski A. Nonlinear polychromatic and noise artifacts in x-ray computed tomography images. Journ Comp Assist Tom. 1979;3(4):519–526. [PubMed]
  • Duerinckx AJ, Macovski A. Information and artifact in computed tomography image statistics. Med Phys. 1980;7(2):127–134. [PubMed]
  • Eckstein MP, Bartroff JL, Abbey CK, Whiting JS, Bochud FO. Automated computer evaluation and optimization of image compression of x-ray coronary angiograms for signal known exactly detection tasks. Optics Express. 2003;11(5):460–475. [PubMed]
  • Eckstein MP, Whiting JS. Lesion detection in structured noise. Academic Radiology. 1995;2:249–253. [PubMed]
  • Faulker K, Moores B. Noise and contrast detection in computed tomography images. Phys Med Biol. 1984;29(4):329–339. [PubMed]
  • Gies M, Kalender WA, Wolf H, Suess C, Madsen MT. Dose reduction in CT by anatomically adapted tube current modulation. I simulation studies, Med Phys. 1999;26(11):2235–2247. [PubMed]
  • Gifford HC, King MA, Pretorius PH, Wells RG. A comparison of human and model observers in multislice LROC studies. IEEE Trans Med Imag. 2005;24(2):160–169. [PubMed]
  • Gifford H, Kinahan PE, Lartizien C, King M. Evaluation of multiclass model observers in PET LROC studies. IEEE Trans Nuc Sci. 2007;54(1):116–123. [PMC free article] [PubMed]
  • Gore J, Tofts P. Statistical limitations in computed tomography. Phys Med Biol. 1978;23(6):1176–1182. [PubMed]
  • Herman GT. On the noise in images produced by computed tomography. Computer Graphics and Image Processing. 1980;12(3):271–285.
  • Hsieh J. Computed Tomography Principles, Design, Artifacts, and Recent Advances. SPIE Press; 2003.
  • Huesman R. A new fast algorithm for the evaluation of regions of interest and statistical uncertainty in computed tomography. Phys Med Biol. 1984;29(5):543–552. [PubMed]
  • Huesman R, Gullberg G, Greenberg W, Budinger T. Technical Report PUB-214. Lawrence Berkeley Laboratory; 1977. Users manual: Donner algorithms for reconstruction tomography.
  • Kak AC, Slaney M. Principles of Computerized Tomographic Imaging. IEEE Press; 1988. republished by SIAM in 2001.
  • Kalender WA. X-ray computed tomography (review) Phys Med Biol. 2006;51(13):R29–R43. [PubMed]
  • Kalra MK, Maher MM, Toth TL, Schmidt B, Westerman BL, Morgan HT, Saini S. Techniques and applications of automatic tube current modulation for CT. Radiology. 2004;233(3):649–657. [PubMed]
  • Kim JS, Kinahan PE, Lartizien C, Comtat C, Lewellen TK. A comparison of planar versus volumetric numerical observers for detection task performance in whole-body pet imaging. IEEE Trans Nuc Sci. 2004;51(1):34–40.
  • LaRoque SJ, Sidky EY, Edwards DC, Pan X. Evaluation of the channelized Hotelling observer for signal detection in 2D tomographic imaging. Medical Imaging 2007, Vol 6515 of Proc of SPIE. 2007:651514.
  • Myers KJ, Barrett HH. Addition of a channel mechanism to the ideal-observer model. J Opt Soc Amer A. 1987;4(12):2447–2457. [PubMed]
  • Noo F, Defrise M, Clackdoyle R, Kudo H. Image reconstruction from fan-beam projections on less than a short scan. Phys Med Biol. 2002;47(14):2525–2546. [PubMed]
  • Pan X. Optimal noise control in and fast reconstruction of fan-beam computed tomography image. Med Phys. 1999;26(5):689–697. [PubMed]
  • Pan X, Yu L. Image reconstruction with shift-variant filtration and its implication for noise and resolution properties in fan-beam computed tomography. Med Phys. 2003;30(4):590–600. [PubMed]
  • Parker DL. Optimal short-scan convolution reconstruction for fanbeam CT. Med Phys. 1982;9(2):254–257. [PubMed]
  • Riederer SJ, Pelc NJ, Chesler DA. The noise power spectrum in computed x-ray tomography. Phys Med Biol. 1978;23(3):446–454. [PubMed]
  • Riordan J. Moment recurrence relations for binomial, poisson, and hypergeometric frequency distributions. The Annals of Mathematical Statistics. 1937;8(2):103–111.
  • Rockmore A, Macovski A. A maximum likelihood approach to transmission image reconstruction from projections. IEEE Trans Nuc Sci. 1977;24(3):1929–1935.
  • Seibert J, Boone J. X-ray scatter removal by deconvolution. Med Phys. 1988;15(4):567–575. [PubMed]
  • Tanaka E, Iinuma T. Correction functions for optimizing the reconstructed image in transverse section scan. Phys Med Biol. 1975;20(5):789–798. [PubMed]
  • Tanaka E, Iinuma TA. Correction functions and statistical noises in transverse section picture reconstuction. Computers in Biology and Medicine. 1976;6(4):295–306. [PubMed]
  • Tanaka E, Murayama H. Proceedings International Workshop on Physics and Engineering in Medical Imaging. IEEE Computer Society; 1982. Properties of statistical noise in positron emisson tomography; pp. 158–164.
  • Wang J, Lu H, Li T, Liang Z. An alternative solution to the nonuniform noise propagation problem in fan-beam FBP image reconstruction. Med Phys. 2005;32(11):3389–3394. [PMC free article] [PubMed]
  • Whiting BR, Massoumzadeh P, Earl OA, O’Sullivan JA, Snyder DL, Williamson JF. Properties of preprocessed sinogram data in x-ray computed tomography. Med Phys. 2006;33(9):3290–3303. [PubMed]
  • Wollenweber S, Tsui B, Lalush D, Frey E, LaCroix K, Gullberg G. Comparison of hotelling observer models and human observers in defect detection from myocardial SPECT imaging. IEEE Trans Nuc Sci. 1999;46(6):2098–2103.
  • Yu L, Pan X. Half-scan fan-beam computed tomography with improved noise and resolution properties. Med Phys. 2003;30(10):2629–2637. [PubMed]
  • Zeng GL. Nonuniform noise propagation by using the ramp filter in fan-beam computed tomography. IEEE Trans Med Imag. 2004;23(6):690–695. [PubMed]
  • Zeng GL, Gullberg GT. A channelized-hotelling-trace collimator design method based on reconstruction rather than projections. IEEE Trans Nuc Sci. 2002;49(5):2155–2158.
  • Zhang Y, Pham BT, Eckstein MP. The effect of nonlinear human visual system components on performance of a channelized hotelling observer in structured backgrounds. IEEE Trans Med Imag. 2006;25(10):1348–1362. [PubMed]
  • Zhu L, Starlack J. A practical reconstruction algorithm for CT noise variance maps using FBP reconstruction. Medical Imaging 2007, Vol 6510 of Proc of SPIE. 2007:651023.