An ideal-observer framework to investigate signal detectability in diffuse optical imaging

Abhinav K. Jha; Eric Clarkson; Matthew A. Kupinski

doi:10.1364/BOE.4.002107

1. Introduction

Diffuse optical tomography (DOT) is emerging as a novel non-invasive medical imaging modality. The objective in DOT is to probe the tissue with near-infrared (NIR) light to determine the absorption and scattering properties of the tissue, which can then be used by the physicians for patient diagnosis [1–4]. DOT has found applications in breast-cancer detection and characterization [5, 6], in functional brain imaging [7, 8], in imaging of small joints for early diagnosis of rheumatoid arthritis [9], and in small-animal imaging for studying physiological processes and pathologies [10]. Many research groups are engaged in designing novel DOT imaging systems [11–14], and thus there is a requirement for methods to assess the performance of these systems in clinically relevant tasks [15, 16]. One such important task in DOT is the detection of absorptive heterogeneous regions, which we refer to as the signal, within the tissue. In DOT, often the anomalous regions are distinguished by a higher absorption coefficient compared to the rest of the tissue [17–22]. However, the sensitivity of diffuse-light measurements drops off quickly with penetration depth due to the high scattering in the tissue, which leads to poor depth resolution in DOT [23,24]. Thus, an important design parameter for a DOT system is the signal detectability provided by the system at larger penetration depth. Similarly, the variation in signal detectability as a function of the signal size and contrast are also important design parameters.

To investigate signal detectability in DOT for system optimization, an appropriate method is to compute the performance of the Bayesian ideal observer on the task of detecting the signal in the tissue [25, 26]. The Bayesian ideal observer uses the likelihood ratio to compute a test statistic from the image data and compares it against a threshold. This observer uses all information available to make decisions, and is optimized by several criteria. For example, the ideal observer maximizes the area under the receiver operating characteristic (ROC) curve, a quantity also known as the area under curve (AUC). Moreover, the ideal observer uses the raw data, and does not suffer from information loss due to the reconstruction process. However, determining the AUC of an ideal observer is often computationally complex [19, 26]. Thus, an approximation to the AUC of the ideal observer that is faster to compute can be very useful for image quality studies.

The use of Fisher information to determine the information content of measured data when performing an estimation task is quite well established [25]. For a given estimation task, the inverse of the Fisher information matrix (FIM) gives the Cramér-Rao bound (CRB), which is a lower bound on the variance of any unbiased estimator that processes the measured data. More recently, it has been also shown that when we are trying to detect a small change in the vector parameter of a parametrized family of probability distribution functions (PDFs), there is a connection between the Fisher information and the performance of the ideal observer on a detection task, as measured by the AUC [26, 27]. Therefore, the Fisher information can serve as a useful surrogate figure of merit (SFoM) that approximates ideal observer performance in detection tasks. In this paper, our objective is to demonstrate the usage of Fisher information as a SFoM to evaluate signal detectability in DOT. An important quantity required to evaluate the Fisher information in DOT is the gradient of the photon distribution function with respect to the absorption coefficient. We will derive an analytic expression for this purpose, along with an analytic expression for the gradient of the photon distribution function with respect to the scattering coefficient. These derived analytic expressions can also be used to perform image reconstruction in DOT using gradient-based reconstruction schemes.

Currently most methods to evaluate DOT systems or investigate signal detectability in DOT systems do not use observer studies [16]. Metrics such as contrast-to-noise ratio (CNR) and positional errors (PE) [17, 18, 24, 28] have been used to evaluate the performance of depth resolution provided by DOT reconstruction algorithms. The PE metric quantifies the error in locating a signal, while the CNR gives a measure of the whether the object can be detected from the background. These metrics are not ideal to investigate signal detectability for DOT system optimization since they require a reconstruction step that often leads to information loss from the acquired data. Ziegler et al.[29, 30] have devised a statistical method to use raw data to evaluate DOT systems, but their method essentially just gives us one point on the ROC curve, and therefore is not comprehensive. A probabilistic detection theory has been devised by Morgan et al.[20], which evaluates the detection performance of a system based on the area of overlap of the signal present and signal absent PDFs of the data. However, this approach is not suitable for multidimensional data, since it is very complicated to compute area under multidimensional PDFs, and often, different components of the data vector contribute different weights to the final test statistic that is used to make a decision about the presence of the signal. There have been previous studies on using Hotelling observers to evaluate various DOT system configurations [19, 31].However, a framework to evaluate signal detectability in DOT systems using ideal observers has not been proposed. It is the objective of this paper to present a framework that uses the ideal observer, using a easier-to-compute SFoM, to evaluate signal detectability in DOT. We first derive the framework to evaluate detectability for a general DOT setup. Using software to simulate photon propagation developed by our research group [32,33], we present the implementation of the proposed method for a specific test DOT system. We then demonstrate the application of the derived framework to investigate signal detectability in this DOT system for some simple signal known exactly/background known exactly (SKE/BKE) tasks.

2. Methods

2.1. Relation between signal detectability and Fisher information in DOT

In this section, we will analyze the relation between signal detectability and Fisher information for a general diffuse optical imaging setup. Let us consider a DOT setup with M detector elements. Let the flux measured by the m^th detector element be denoted by g_m, and the image acquired using the M detector elements by the M-dimensional vector g. Let us also denote the signal-absent and the signal-present hypothesis by H₀ and H₁, respectively. The objective in the detection task is to determine if g is a sample from the signal-absent PDF, denoted by pr(g|H₀), or the signal-present PDF, denoted by pr(g|H₁). The Bayesian ideal observer performs this task by evaluating the likelihood ratio t(g) and comparing it to a threshold. The likelihood ratio t(g) is given by

t (g) = \frac{pr (g | H_{1})}{pr (g | H_{0})} .

To use this observer in the DOT context, let us denote the absorption and scattering coefficients at location r by μ_a(r) and μ_s(r), respectively. We discreteize this functions using a certain spatial basis with basis functions given by ϕ_n(r). The representation of the absorption and scattering functions in this basis is given by

μ_{a} (r) = \sum_{n = 1}^{N} μ_{a, n} ϕ_{n} (r)

μ_{s} (r) = \sum_{n = 1}^{N} μ_{s, n} ϕ_{n} (r)

where μ_a,n and μ_s,n denote the absorption and scattering coefficients in the considered spatial basis. These basis functions could be the commonly used voxel basis functions. However, we are not placing any restriction on the spatial support of these basis functions. The basis function could in fact represent the support of different anatomical structures in the tissue, where a given anatomical structure is characterized by the same absorption and scattering coefficient. It must be pointed out that in general, the equality in Eqs. (2) and (3) holds in the limit of having an infinite number of basis functions, since we are approximating a continuous function by a finite basis set. Let us denote the N-dimensional vector of the absorption coefficients by μ_a and the vector of scattering coefficients by μ_s. For notational simplicity, let us denote μ = {μ_a, μ_s}. Also, again for notational simplicity, we define a general coefficient function μ(r) that denotes the absorption/scattering coefficient at location r. The representation of this general coefficient function in the spatial basis is given by

μ (r) = \sum_{n = 1}^{N} μ_{n} ϕ_{n} (r),

where μ_n denotes the absorption/scattering coefficients in the considered spatial basis.

In the detection task we consider, the signal-present hypothesis is distinguished from the signal-absent hypothesis by a small change in the parameter vector μ. For the signal-present hypothesis H₁ and signal-absent hypothesis H₀, the parameter vectors have value μ and μ₀, respectively. Thus, using Eq. (1), the likelihood ratio is given by

t_{μ} (g) = \frac{pr (g | μ)}{pr (g | μ_{0})} .

To relate the ideal-observer AUC with the Fisher information, we apply the general framework proposed in Clarkson et al.[27]. We first define a detectability parameter, which is related to the ideal-observer AUC as [26]

AUC (μ) = \frac{1}{2} [1 + erf (\frac{d (μ)}{2})],

where the ideal-observer AUC and detectability when the parameter vector is μ are denoted by AUC(μ) and d(μ), respectively. When the change between μ and μ₀ is small, we obtain a second-order approximation to the square of the detectability [27]:

d^{2} (μ) ≅ {(μ - μ_{0})}^{†} F (μ_{0}) (μ - μ_{0})

where F(μ₀) denotes the FIM evaluated at μ₀. Therefore, the detectability, and consequently the AUC for an ideal observer trying to detect a small change in μ, is directly related to the Fisher information evaluated at μ₀. We describe the procedure to compute the elements of the FIM in the next section.

We would like to mention here that in DOT, often the relative change in absorption coefficient between signal and background could be relatively large, due to the very low value of the absorption coefficient in the background. However, this change in absorption coefficient still translates to a very small change in the final image since the absorption coefficient is generally much lower than the scattering coefficient, and always appears in conjunction with the scattering coefficient in the equation that describes light propagation in tissue. Thus, even when the relative change between the background and signal absorption coefficient values is high, when added with the scattering coefficient, this change is relatively low. As we have also verified, even a large relative change in the value of the absorption coefficient does not cause a significant relative change in the final image. Therefore, the framework that we suggest is valid even for a relatively large change in the absorption coefficient between the signal and the rest of the tissue.

2.2. Computing the Fisher information

To compute the elements of the FIM, we must first obtain the likelihood of the image data as a function of the scattering and absorption coefficient vectors. The major noise source in DOT is often the Poisson-distributed shot noise [19, 34]. Let ḡ_m(μ) denote the noiseless mean image data as a function of the scattering and absorption coefficient vectors μ = {μ_s, μ_a}. Since the Poisson noise in individual detector elements is independent of each other, the PDF for the acquired image data g in our DOT setup is given by

pr (g | μ) = \prod_{m = 1}^{M} \frac{{[{\bar{g}}_{m} (μ)]}^{g_{m}} exp [- {\bar{g}}_{m} (μ)]}{g_{m}!}

Taking the logarithm on both sides, we get the log-likelihood of the image data, which we denote by ℒ(ḡ|μ), as

ℒ (\bar{g} | μ) = log pr (g | μ) = \sum_{m = 1}^{M} g_{m} log {\bar{g}}_{m} (μ) - {\bar{g}}_{m} (μ) - log (g_{m}!) .

The (i, j)^th element of the FIM F(μ) is given by [25]

F_{i j} (μ) = - {〈 \frac{\partial^{2} ℒ}{\partial μ_{j} \partial μ_{i}} 〉}_{g | μ}

where angled brackets denote statistical expectations, and where μ_i denotes the i^th absorption/scattering coefficient as defined in Eq. (4). To obtain the elements of the FIM, we take the derivative of the log-likelihood (Eq. (9)) with respect to the coefficient μ_i. This yields

\frac{\partial ℒ}{\partial μ_{i}} = \sum_{m = 1}^{M} (\frac{g_{m}}{{\bar{g}}_{m} (μ)} - 1) \frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{i}} .

Further taking the derivative of the above expression with respect to μ_j, we get

\frac{\partial^{2} ℒ}{\partial μ_{j} \partial μ_{i}} = \sum_{m = 1}^{M} (- \frac{g_{m}}{{\bar{g}}_{m} {(μ)}^{2}} \frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{j}}) \frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{i}} + (\frac{g_{m}}{{\bar{g}}_{m} (μ)} - 1) \frac{\partial^{2} {\bar{g}}_{m} (μ)}{\partial μ_{j} \partial μ_{i}} .

Averaging the above expression over g and using Eq. (10), we obtain

F_{i j} (μ) = - {〈 \frac{\partial^{2} ℒ}{\partial μ_{j} \partial μ_{i}} 〉}_{g | μ} = \sum_{m = 1}^{M} \frac{1}{{\bar{g}}_{m} (μ)} \frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{j}} \frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{i}},

since 〈g_m〉_g|_μ = ḡ_m(μ). We observe that to evaluate the elements of the FIM, we have to evaluate the gradient of the mean image data with respect to the absorption and scattering coefficients. This is the topic of the next section.

2.3. Evaluating the gradient of the mean image data

To evaluate the gradient of the mean image data, we have developed an analytic expression based on the integral form of the radiative transport equation (RTE), which is the equation used to describe photon propagation through any media. The fundamental radiometric quantity that we describe using the RTE is the photon distribution function w(r, ŝ, ℰ, t). In terms of photons, w(r, ŝ, ℰ, t)ΔVΔΩΔℰ can be interpreted as the number of photons contained in volume ΔV centered on the 3D position vector r = (x,y,z), traveling in a solid angle ΔΩ about direction ŝ = (θ, ϕ), and having energies between ℰ and ℰ + Δℰ at time t. Another radiometric quantity required to describe the RTE is the emission function Ξ(r, ŝ, ℰ, t). Like the distribution function, Ξ(r, ŝ, ℰ, t)ΔVΔΩΔℰ can be interpreted as the number of photons injected per second into volume ΔV in energy range Δℰ, over solid angle ΔΩ, and at time t. In the DOT implementation, we assume a mono-energetic time-independent emission source, so that the emission function can be written as Ξ(r, ŝ). Also, the dominant scattering mechanism in DOT is elastic scattering and the scattered photon does not lose any energy. Since there is no other energy-loss mechanism for the photons, we can completely drop the dependence of the distribution function on energy. Also, since we consider a time-independent emission source, the dependence of the distribution function on time is also dropped, and we write the distribution function as w(r, ŝ).

Let 𝒦 and 𝒳 denote the scattering and attenuation operators in integral form, which represent the effect of scattering, and the effect of attenuation and propagation of photons, respectively. The effect of the scattering operator 𝒦 on the distribution function is given by [25]

[𝒦 w] (r, \hat{s}) = μ_{s} (r) c_{m} \int_{4 π} d Ω^{'} p (\hat{s}, {\hat{s}}^{'} | r) w (r, {\hat{s}}^{'}),

where c_m denote the speed of light in the medium, ŝ and ŝ′ denote the direction of the outgoing and incoming photons, respectively, and p(ŝ, ŝ′|r) denotes the scattering phase function, which in biological tissue is typically given by the Henyey-Greenstein function [35, 36]:

p (\hat{s}, {\hat{s}}^{'} | r) = \frac{1}{4 π} {\frac{1 - α^{2}}{{[1 + α^{2} - 2 α cos (\hat{s} \cdot {\hat{s}}^{'})]}^{3 / 2}}}

where the anisotropy factor α characterizes the angular distribution of scattering in the tissue. The attenuation operator 𝒳 is the standard attenuated X-ray transform, and its effect on the distribution function is given by [25]

[𝒳 w] (r, \hat{s}) = \frac{1}{c_{m}} \int_{0}^{\infty} d l w (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})],

where μ_tot(r) = μ_a(r) + μ_s(r) is the total attenuation coefficient. In terms of the defined attenuation and scattering operators, for a mono-energetic time-independent source Ξ(r, ŝ), the RTE can be written as [25, 32]

w = 𝒳 Ξ + 𝒳 𝒦 w

To evaluate the gradient of the mean image with respect to the n^th absorption/scattering coefficient μ_n, we note that ḡ_m(μ) is computed as

{\bar{g}}_{m} (μ) = \int d^{3} r \int d Ω h_{m} (r, \hat{s}) w (r, \hat{s}),

where h_m(r, ŝ) denotes the sensitivity of the m^th detector pixel to the distribution function. We can write the above equation in inner product notation as

{\bar{g}}_{m} (μ) = (h_{m}, w) .

Taking the gradient of ḡ_m(μ) with respect to μ_n, we get

\frac{\partial {\bar{g}}_{m} (μ)}{\partial μ_{n}} = (h_{m}, \frac{\partial w}{\partial μ_{n}}) .

Thus, to evaluate the gradient of the mean image data, we have to evaluate the gradient of the distribution function with respect to μ_n. Using Eq. (17), we can write the expression for this gradient as

\frac{\partial w}{\partial μ_{n}} = \frac{\partial 𝒳}{\partial μ_{n}} Ξ + \frac{\partial 𝒳}{\partial μ_{n}} 𝒦 w + 𝒳 \frac{\partial 𝒦}{\partial μ_{n}} w + 𝒳 𝒦 \frac{\partial w}{\partial μ_{n}} .

Of the four terms in the above equation, the first two terms can be considered as a new operator

\frac{\partial 𝒳}{\partial μ_{n}}

acting on another distribution function. Let us evaluate the first term. We know from Eq. (16) that the effect of the attenuated X-ray transform on the source function is given by

[𝒳 Ξ] (r, \hat{s}) = \frac{1}{c_{m}} \int_{0}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] .

Note that

\begin{array}{l} μ_{tot} (r) & = μ_{s} (r) + μ_{a} (r) \\ = \sum_{n = 1}^{N} (μ_{s, n} + μ_{a, n}) ϕ_{n} (r) \end{array}

where in the last step, we have used Eqs. (2) and (3). Using the above representation and taking the derivative of Eq. (22) with respect to μ_n, where as mentioned earlier, μ_n could denote the absorption/scattering coefficient, we get

\frac{\partial 𝒳}{\partial μ_{n}} Ξ (r, \hat{s}) = - \frac{1}{c_{m}} \int_{0}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] [\int_{0}^{l} d l^{″} ϕ_{n} (r - \hat{s} l^{″})] .

Defining the step function

step (l) = {\begin{cases} 0, l < 0 \\ 1, otherwise \end{cases}

we can rewrite the integral over l″ in Eq. (24) as

\begin{array}{l} \frac{\partial 𝒳}{\partial μ_{n}} Ξ (r, \hat{s}) = - \frac{1}{c_{m}} \times \\ \int_{0}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] [\int_{0}^{\infty} d l^{″} [step (l^{″}) - step (l^{″} - l)] ϕ_{n} (r - \hat{s} l^{″})] \\ = - \frac{1}{c_{m}} \int_{0}^{\infty} d l \int_{0}^{\infty} d l^{″} [step (l^{″}) - step (l^{″} - l)] ϕ_{n} (r - \hat{s} l^{″}) Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] \end{array}

where in the second step, we have done a simple rearrangement of terms. Interchanging the order of integration, we obtain

\begin{array}{l} \frac{\partial 𝒳}{\partial μ_{n}} Ξ (r, \hat{s}) \\ = - \frac{1}{c_{m}} \int_{0}^{\infty} d l^{″} ϕ_{n} (r - \hat{s} l^{″}) \int_{0}^{\infty} d l [step (l^{″}) - step (l^{″} - l)] Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] \\ = - \frac{1}{c_{m}} \int_{0}^{\infty} d l^{″} ϕ_{n} (r - \hat{s} l^{″}) \int_{l^{″}}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{0}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})], \end{array}

where in the second step, we use the fact that as l″ varies from 0 to ∞, step(l″) = 1, and that 1 − step(l″ − 1) exists only when l > l″. Splitting the exponential integral over l′ into two parts, we get

\begin{array}{l} \frac{\partial 𝒳}{\partial μ_{n}} Ξ (r, \hat{s}) = - \frac{1}{c_{m}} \times \\ \int_{0}^{\infty} d l^{″} ϕ_{n} (r - \hat{s} l^{″}) exp [- \int_{0}^{l^{″}} d l^{'} μ_{tot} (r - \hat{s} l^{'})] \int_{l^{″}}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{l^{″}}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'}) d l^{'}] . \end{array}

This equation is now starting to resemble another attenuation transform, except for the integral over l. Simplifying this integral further by replacing l − l″ by l̃, we get

\begin{array}{l} \int_{l^{″}}^{\infty} d l Ξ (r - \hat{s} l, \hat{s}) exp [- \int_{l^{″}}^{l} d l^{'} μ_{tot} (r - \hat{s} l^{'})] \\ = \int_{0}^{\infty} d \tilde{l} Ξ (r - \hat{s} l^{″} - \hat{s} \tilde{l}, \hat{s}) exp [- \int_{l^{″}}^{l^{″} + \tilde{l}} d l^{'} μ_{tot} (r - \hat{s} l^{'})] \\ = \int_{0}^{\infty} d \tilde{l} Ξ (r - \hat{s} l^{″} - \hat{s} \tilde{l}, \hat{s}) exp [- \int_{0}^{\tilde{l}} d \hat{l} μ_{tot} (r - \hat{s} l^{″} - \hat{s} \hat{l})], \end{array}

where in the last step, we do another change of variables in the integral over l′ by replacing l′ − l″ by l̂. We recognize that the last term is very similar to the attenuated X-ray transform. More precisely, using Eq. (16)

\int_{0}^{\infty} d \tilde{l} Ξ (r - \hat{s} l^{″} - \hat{s} \tilde{l}, \hat{s}) exp [- \int_{0}^{\tilde{l}} d \hat{l} μ_{tot} (r - \hat{s} l^{″} - \hat{s} \hat{l})] = c_{m} [𝒳 Ξ] (r - \hat{s} l^{″}, \hat{s}) .

Substituting this in Eq. (28), we get

\frac{\partial 𝒳}{\partial μ_{n}} Ξ (r, \hat{s}) = - \frac{1}{c_{m}} \int_{0}^{\infty} d l^{″} ϕ_{n} (r - \hat{s} l^{″}) c_{m} {[𝒳 Ξ] (r - \hat{s} l^{″}, \hat{s})} exp [- \int_{0}^{l^{″}} d l^{'} μ_{tot} (r - \hat{s} l^{'})],

which, using Eq. (16), can be written in operator notation as

\frac{\partial 𝒳}{\partial μ_{n}} Ξ = - c_{m} 𝒳 (ϕ_{n} 𝒳 Ξ)

This gives the first term in Eq. (21). The second term in Eq. (21) is very similar to the first term, and can be written as

\frac{\partial 𝒳}{\partial μ_{n}} 𝒦 w = - c_{m} 𝒳 (ϕ_{n} 𝒳 𝒦 w)

Using Eqs. (32) and (33) we get

\begin{array}{l} \frac{\partial 𝒳}{\partial μ_{n}} Ξ + \frac{\partial 𝒳}{\partial μ_{n}} 𝒦 w & = - c_{m} 𝒳 [ϕ_{n} (𝒳 Ξ + 𝒳 𝒦 w)] \\ = - c_{m} 𝒳 (ϕ_{n} w), \end{array}

where in the second step, we have used Eq. (17).

The third term in Eq. (21) involves a gradient of the scattering operator. Evaluating this gradient is easy, since the scattering operator, given by Eq. (14), is linear in μ_s. Using Eq. (4), the gradient of the scattering operator is obtained as

\begin{array}{l} \frac{\partial 𝒦}{\partial μ_{n}} w & = \frac{\partial}{\partial μ_{n}} c_{m} μ_{s} (r) \int d Ω^{'} p (\hat{s}, {\hat{s}}^{'}) w (r, {\hat{s}}^{'}) \\ = β c_{m} ϕ_{n} (r) \int d Ω^{'} p (\hat{s}, {\hat{s}}^{'}) w (r, {\hat{s}}^{'}), \end{array}

where β = 0 when μ_n denotes the absorption coefficient and 1 when μ_n denotes the scattering coefficient. Defining another operator 𝒦₁, whose effect on the distribution function is given by

[𝒦_{1} w] (r, \hat{s}) = c_{m} \int d Ω^{'} p (\hat{s}, {\hat{s}}^{'}) w (r, {\hat{s}}^{'}),

we can rewrite Eq. (35) in operator notation as

\frac{\partial 𝒦}{\partial μ_{n}} w = β ϕ_{n} 𝒦_{1} w .

Thus, the third term in Eq. (21) can be written as

𝒳 \frac{\partial 𝒦}{\partial μ_{n}} w = 𝒳 (β ϕ_{n} 𝒦_{1} w) .

Substituting Eq. (34) and (38) in Eq. (21), we get

\frac{\partial w}{\partial μ_{n}} = 𝒳 (- c_{m} ϕ_{n} w + β ϕ_{n} 𝒦_{1} w) + 𝒳 𝒦 \frac{\partial w}{\partial μ_{n}} .

Comparing with Eq. (17), we realize that the above equation is just the RTE, with the source term given by

s_{n} = - c_{m} ϕ_{n} w + β ϕ_{n} 𝒦_{1} w .

Thus, to evaluate the gradient of the photon distribution function with respect to μ_n, we just have to compute each term of the RTE with the source term as −ϕ_nc_mw +βϕ_n𝒦₁w, i.e. the distribution function (−c_mw + β𝒦₁w) that exists over the spatial support of the n^th basis function.

To summarize, to evaluate the detectability of the signal at a certain location, we use Eq. (7). The Fisher information in this equation is computed using Eq. (13). Computing the Fisher information requires determining the gradient of the photon distribution function, which is evaluated by executing the RTE with a source term that is dependent on the spatial support of the signal.

2.4. Implementation

We implement the proposed framework to investigate signal detectability in DOT using the C programming language on a computing system with a 2.27 GHz Intel Xeon quad core E5520 processor as the central processing unit (CPU) running a 64-bit Linux operating system, and consisting of four NVIDIA Tesla C2050 graphics processing units (GPUs). To evaluate the detectability, we must compute the mean image data and its gradient. Computing these quantities requires solving the RTE. We have developed mathematical methods to solve the RTE using a Neumann-series form [32, 33]. We have also developed software to model light propagation in completely non-uniform three-dimensional small-geometry anisotropic scattering media for a DOT setup [33]. The software has been developed on the NVIDIA Tesla C2050 graphics processing units (GPUs) using the compute unified device architecture (CUDA) programming framework for computational efficiency, and provides up to two orders of magnitude speedup compared to a non-GPU implementation. The Neumann-series RTE software is developed to solve the RTE with a laser beam as the source term, and execution of this software computes the mean image data for the DOT setup. Using this software, the distribution function in all the voxels is also computed. This quantity is stored, since it is required to compute the source term (Eq. (40)) for the RTE to compute the gradient of the mean image data.

To evaluate the gradient of the mean image data, we first note that the RTE for the gradient of the photon distribution function given by Eq. (39) can be written in the following form as a Neumann series [32]:

\frac{\partial w}{\partial μ_{n}} = 𝒳 s_{n} + 𝒳 𝒦 𝒳 s_{n} + 𝒳 𝒦 𝒳 s_{n} + \dots

where s_n is the source term as given by Eq. (40). To implement this equation, we further rewrite in a spherical harmonic and voxel basis as below [32, 33]:

{W_{d}}^{'} = {AS}_{n} + {ADAS}_{n} + {ADADAS}_{n} + \dots

where A and D denote the attenuation and scatter operators in voxel and spherical harmonic basis respectively, and S_n and W_d′ are the source term and the derivative of the photon distribution function in the spherical harmonic and voxel basis, respectively. We have derived the expressions and implemented the attenuation and scattering operators on the GPU for a nonuniform anisotropic scattering media [33]. To evaluate the gradient, we need to substitute the source term by the expression given by Eq. (40) in the spherical harmonic and voxel basis. This quantity can be easily computed from the distribution function that we had computed and stored earlier. Following this, an iterative application of the attenuation and scattering operators until convergence, and a final application of Eq. (20) yields the gradient of the mean image data [33]. Using the mean image data and its gradient, the Fisher information, and therefore the detectability, can be evaluated for any spatial support of the scattering or absorptive signal.

2.5. Methods for the specific SKE/BKE task

To demonstrate the application of the proposed scheme, we consider the task of detecting a high-absorption region, i.e. the signal, in an otherwise homogeneous background. Let the spatial support of the background and the signal be denoted by ϕ₀(r) and ϕ₁(r), respectively. Also, let the absorption coefficient of the background and signal be denoted by μ_a,0 and μ_a,1, respectively. The signal and the background are assumed to have the same scattering coefficient. Let us also denote the difference between μ_a,0 and μ_a,1 by Δμ_a. Using Eq. (7), the expression for detectability in this case will be given by

d^{2} (μ) = Δ μ_{a}^{2} F_{a 1, a 1},

where, using Eq. (13),F_a1,a1 is given by

F_{a 1, a 1} = \sum_{m = 1}^{M} \frac{1}{{\bar{g}}_{m}} {(\frac{\partial {\bar{g}}_{m}}{\partial μ_{a, 1}})}^{2} .

Using Eq. (20), the gradient of the image data with respect to μ_a,1 is given by

\frac{\partial {\bar{g}}_{m}}{\partial μ_{a, 1}} = (h_{m}, \frac{\partial w}{\partial μ_{a, 1}}) .

To evaluate the gradient of the photon distribution with respect to μ_a,1, in accordance with Eq. (39) we execute the RTE with the source term as

s_{n} = - c_{m} ϕ_{1} w .

Therefore, the mean output image ḡ and the gradient of the mean output image are evaluated using the RTE. To vary the signal depth, we vary the spatial support of the signal, given by ϕ₁(r). The variation of detectability with other parameters, such as the size and contrast of the signal, and the scattering coefficient of the background can also be studied using this framework.

Before presenting the simulation studies, we would like to mention that the Neumann-series method is accurate in many scenarios where the conventionally used diffusion approximation breaks down, such as DOT setups with collimated light sources [32,37,38], optically-thin media [32, 39], media with low-scattering void-like regions [33, 40] and media where the absorption coefficient is similar to the scattering coefficient [33,41,42]. Investigating signal detectability in media where the diffusion approximation breaks down is an important research problem [16]. In our experiments, we simulate DOT setups in which the diffusion approximation breaks down, and study the variation in signal detectability as we vary the signal and phantom properties. Thus these studies present some interesting insights on signal detectability in DOT, and also evaluates phantom configurations where the diffusion approximation breaks down. However, our main intent behind presenting these experimental results is to demonstrate the validity of the proposed framework.

3. Experiments and results

We now demonstrate the application of the proposed framework to investigate signal detectability for a particular test DOT system, as the various properties of the phantom and the signal are varied. As mentioned earlier, our intent in performing these observer studies is to validate the proposed framework. Thus, we will be simulating DOT setups where we can predict the expected observer output, and validate whether the proposed SFoM also computes a similar result. Our simulated DOT setup is as shown in Fig. 1, where a collimated 10 mW laser source emits a beam of NIR light with a transverse circular profile with radius 0.5 mm. The beam is incident approximately on the center of the entrance face of the scattering medium, where the center of the beam has x − y coordinates as (0.5 mm, 0.5 mm). The scattering medium is a 3-D phantom of dimensions 2 × 2 × 2 cm³, characterized by a reduced scattering coefficient μ′_s = μ_s(1 − g). The background of the phantom is homogeneous with absorption and scattering coefficients given by 0.01 cm⁻¹ and 1 cm⁻¹, respectively. This is a small-geometry media where the light source is collimated, which is a case where the diffusion approximation breaks down [37, 39]. The experiments are carried out for non-reentry boundary conditions, so the refractive index of the medium is fixed to unity. The value of the anisotropic coefficient α of the media is set to 0, partly because it is easier to predict the expected observer output with such media. The signal is spherical in shape with a radius of 1.5 mm and has the same scattering coefficient as the background, but a higher absorption coefficient of 0.011 cm⁻¹. The DOT setup consists of a pixellated detector in the x−y plane, consisting of 20×20 pixels, and of dimensions 2×2 cm², which is in contact with the entrance face of the tissue and measures the reflected intensity. In the Neumann-series method, for all the experiments, the highest degree in the spherical harmonic expansion is 3. For this DOT setup and the specified signal and phantom properties, the computation time to determine the SFoM was close to a minute.

Fig. 1 The simulated DOT setup for our experiments.

Download Full Size | PDF

3.1. Signal detectability as a function of depth

In the first experiment, the signal is placed approximately at the center of the x − y plane. The x − y coordinates of the center of the signal are (0.5 mm, 0.5 mm), and its location is varied along the z dimension. We would expect that since we are using only the reflected intensity to make inferences about the presence of the signal, the detectability would decrease as the signal depth increases. The detectability is computed using the developed software and plotted as a function of the signal location in Fig. 2, and we observe the expected trend.

Fig. 2 The signal detectability as a function of the depth of the signal for on-axis and off-axis signal locations.

Download Full Size | PDF

We then shift the location of the signal to an off-axis position (x,y) = (5.5 mm, 5.5 mm), and repeat the above experiment. The plot of detectability vs. signal depth is plotted in Fig. 2. As we would expect, when the depth increases, the signal detectability reduces. Also, we expect that the detectability should be lower for an off-axis signal compared to an on-axis signal, and the plot shows this trend.

3.2. Signal detectability as a function of the scattering coefficient of the tissue

In the next experiment, we vary the scattering coefficient of the tissue. The signal has the same scattering coefficient as the background. The object is placed at four different on-axis locations at 2.5 mm, 7.5 mm, 12.5 mm, and 17.5 mm from the reflected face, respectively. It would be expected that as the scattering coefficient increases, due to increased scattering, the detectability of the signal will decrease. However, the process of scattering also contributes to attenuation, and thus, an increase in the scattering coefficient also causes an increase in the value of the attenuation coefficient of the signal, which should lead to increased detectability.

We plot the signal detectability as a function of the scattering coefficient for these four locations of the signal in Fig. 3. We find that the plot reflects the discussed trade-off. When the signal depth is small, scattering has very little effect on detectability. Thus, increasing the scattering coefficient up to a certain value leads to improved detectability. However, as the signal depth increases, the effect of increase in scattering in the media tends to outweigh the increase in the attenuation coefficient. At moderate signal-depth values, the detectability is maximized at an optimal scattering coefficient value. For higher signal-depth values, the scattering effect dominates and an increase in scattering coefficient leads to a decrease in detectability.

Fig. 3 The signal detectability as a function of the scattering coefficient for different signal depths

Download Full Size | PDF

3.3. Signal detectability as a function of signal size and signal depth

Due to the poor spatial resolution in DOT, it is important to investigate the variation in detectability as a function of the signal size. To study this variation, in the next experiment, we vary the radius of the spherical signal from 0.5 to 2 mm, and for each size, vary the depth from 2.5 mm to 16.5 mm. The result of this experiment is plotted in Fig. 4. As expected, an increase in the size of the signal leads to an increase in detectability, for all values of signal depth.

Fig. 4 The signal detectability as a function of the signal size for different signal depths

Download Full Size | PDF

3.4. Signal detectability as a function of signal size and signal contrast

Another important study in DOT is the variation in detectability as a function of the signal size and signal contrast [16], where signal contrast refers to the difference between the absorption coefficients of the signal and the background. We can perform this observer study using the proposed framework,. In our study, we vary the absorption coefficient of a spherical signal from 0.011 cm⁻¹ to 0.05 cm⁻¹. Simultaneously, the radius of the signal is also increased from 0.5 mm to 2.5 mm. The signal is placed at a depth of 2.5 mm from the exit face of the tissue. From the results plotted in Fig. 5, we observe that as the signal contrast increases, the detectability also increases, which is an expected result.

Fig. 5 The signal detectability as a function of the signal contrast for different signal sizes

Download Full Size | PDF

4. Discussions and conclusions

In this paper, we have suggested a framework to evaluate signal detectability in DOT using a SFoM that approximates the performance of the Bayesian ideal observer and is based on the Fisher information. We have implemented this framework using the Neumann-series RTE to simulate light propagation through the tissue. The software has been developed on NVIDIA GPUs for computational efficiency. This framework can be used for various tasks such as evaluating the detectability of an absorptive or scattering signal in any general DOT setup, quantifying the performance of various DOT systems on different detection tasks, and comparing different DOT systems for different detection tasks.

We have also derived an expression to compute the gradient of the mean image data as a function of the scattering and absorption coefficients. This expression is very useful since it just requires executing the RTE with a different source term. Moreover, while we have used the Neumann-series form of the RTE to evaluate the gradient, it can also be computed using other methods to solve the RTE, such as the differential methods [41, 43, 44] and the Monte Carlo methods [45]. The different photon-transport-simulation methods have trade-offs in terms of accuracy and computational efficiency for different kind of phantoms, and thus, depending on the phantom type, the appropriate method can be used to compute the gradient, and implement the proposed SFoM-based scheme. The expression to evaluate the gradient also has other applications, such as in image reconstruction or in studying information content of DOT data with respect to estimation tasks [47].

We have demonstrated the application of the developed SFoM to study signal detectability as a function of signal depth, signal size, signal contrast, and the scattering coefficient of the tissue. The SFoM predicts that while signal detectability decreases with an increase in depth, and increases with an increase in signal contrast and signal size, the relation is not as straightforward when we vary the scattering coefficient of the tissue. While varying the scattering coefficient of the tissue, it is observed that for a given signal depth, there are optimum values of scattering coefficients at which the detectability is maximized. The SFoM predicts expected trends, and thus these studies help to validate the proposed framework.

In this paper, we have investigated the signal detectability for a simple SKE/BKE task. However, the developed framework is general and can be used for studying more rigorous and realistic detection tasks. To apply this scheme in more realistic media, we would have to account for various factors, of which two important factors are medium inhomogeneity and signal and background randomness. Our scheme can be easily applied to heterogeneous media as long as we simulate photon transport using a method that accounts for tissue non-uniformity, such as the Neumann-series-based method that we have developed [33]. To investigate signal detectability with random backgrounds or random signals, we would have to determine the Fisher information for these scenarios. The general form of the Fisher information when the signal and background are random has been derived in Kupinski et al.[46], and can be determined for our case using the methods suggested in this paper. Since the formalism that relates detectability to Fisher information (Eq. (7)) is completely general, observer studies with random signals or backgrounds can be performed using Fisher information as the SFoM.

In the experiments we have performed, we often obtain high detectability values. Detectability values greater than 5 yield very similar AUC values (≈ 1). In such scenarios detectability cannot be used as an indicator to compare the performance of the two systems. However, we obtain such high detectability values in our experiments because we have chosen a very simple SKE/BKE task, the signal is not too deep into the tissue, the scattering coefficient is low, and often the signal is on-axis. For more realistic tasks, the detectability values would be much lower. For example, we observe in Fig. 2 that when the signal is more than 1.5 cm deep into the tissue and off axis, the detectability is less than 0.5. Similarly, when the scattering coefficient is high and the signal is deep, as in the lower two plots in Fig. 3, the detectability values are very low, even when the signal is on axis. When the detectability values lie in the range of 0 to 5, the suggested framework yields different AUC values and can thus be used to compare systems.

The simulation study in this paper has been performed for media with low values of scattering coefficients. This is because the developed Neumann-series RTE approach is computationally intensive for media that have very high values of μ_sH, where μ_s denotes the scattering coefficient of the phantom and H denotes the length of the tissue [32]. With improvement in computational capacity, it should be possible to perform this study using the Neumann-series RTE with higher values of scattering coefficients [33]. Alternatively, we can use the differential or Monte Carlo methods, which are computationally more efficient when the scattering coefficient is high. Also, while in our simulation studies, we have chosen our imaging domain to be a non-diffusion regime, the SFoM-based technique can be definitely used in media where the diffusion approximation is valid, and in those cases, the diffusion-approximation-based methods can be used to evaluate the mean image data.

In our detection task, the signal is a high-absorption region, but there are instances in DOT where we are interested in detecting scattering heterogeneities [19]. While we have not demonstrated this through experiments, the method can easily take care of scattering heterogeneities. We can evaluate the gradient of the photon distribution function with respect to the scattering coefficient by executing the RTE with the source term given by Eq. (40), where we set β = 1. The computed gradient value can then be used to evaluate the Fisher information with respect to the scattering coefficients using Eq. (13). Following this, the detectability for a scattering inhomogeneity can be computed using Eq. (7).

We have performed our analysis with a Poisson noise model, since this is often the major noise source in DOT systems [19, 34]. Our analysis can be extended to detectors that have other noise sources. For example, in the DOT setup being considered, apart from the Poisson noise, the detector might have significant thermal and readout noise, which can be modeled using a Gaussian distribution. Under some assumptions, the combined effect of these Poisson and Gaussian noise sources can be approximated by another normal PDF [48], the first and second-order statistics of which can be determined [47]. The Fisher information for a Gaussian noise model is easy to determine. For detectors with other noise sources, we would need to derive the PDF of the noise distribution, and thus the expression for the Fisher information, following which the proposed scheme can be applied.

Acknowledgments

This work was supported by Canon U. S. A. Inc, and National Institute of Biomedical Imaging and Bioengineering of National Institute of Health under grant number RC1-EB010974, R01-EB000803 and P41-EB002035, and a William Wolfe Family Endowed Scholarship. The authors would also like to thank Dr. Harrison Barrett and the reviewers for their comments.

References and links

1. A. P. Gibson, J. C. Hebden, and S. R. Arridge, “Recent advances in diffuse optical imaging,” Phys. Med. Biol. 50, 1–43 (2005). [CrossRef]

2. D. A. Boas, D. H. Brooks, E. L. Miller, C. A. DiMarzio, M. Kilmer, R. J. Gaudette, and Q. Zhang, “Imaging the body with diffuse optical tomography,” IEEE Signal Process. Mag. 18, 57–75 (2001). [CrossRef]

3. A. Gibson and H. Dehghani, “Diffuse optical imaging,” Phil. Tran. A. Math. Phys. Eng. Sci. 367, 3055–3072 (2009). [CrossRef]

4. H. Dehghani, S. Srinivasan, B. W. Pogue, and A. Gibson, “Numerical modelling and image reconstruction in diffuse optical tomography,” Phil. Trans. Royal Soc. A 367, 3073–3093 (2009). [CrossRef]

5. H. Dehghani, B. W. Pogue, S. P. Poplack, and K. D. Paulsen, “Multiwavelength three-dimensional near-infrared tomography of the breast: initial simulation, phantom, and clinical results,” App. Optics 42, 135–146 (2003). [CrossRef]

6. S. Srinivasan, B. W. Pogue, S. Jiang, H. Dehghani, C. Kogel, S. Soho, J. J. Gibson, T. D. Tosteson, S. P. Poplack, and K. D. Paulsen, “In vivo hemoglobin and water concentrations, oxygen saturation, and scattering estimates from near-infrared breast tomography using spectral reconstruction,” Acad. Radiol. 13, 195–202 (2006). [CrossRef] [PubMed]

7. T. Austin, A. P. Gibson, G. Branco, R. M. Yusof, S. R. Arridge, J. H. Meek, J. S. Wyatt, D. T. Delpy, and J. C. Hebden, “Three dimensional optical imaging of blood volume and oxygenation in the neonatal brain,” Neuroimage 31, 1426–1433 (2006). [CrossRef] [PubMed]

8. B. W. Zeff, B. R. White, H. Dehghani, B. L. Schlaggar, and J. P. Culver, “Retinotopic mapping of adult human visual cortex with high-density diffuse optical tomography,” Proc. Nat. Acad. Sciences 104, 12169–12174 (2007). [CrossRef]

9. A. H. Hielscher, A. D. Klose, A. K. Scheel, B. Moa-Anderson, M. Backhaus, U. Netz, and J. Beuthan, “Sagittal laser optical tomography for imaging of rheumatoid finger joints,” Phys. Med. Biol. 49, 1147–1163 (2004). [CrossRef] [PubMed]

10. A. H. Hielscher, “Optical tomographic imaging of small animals,” Curr. Opinion in Biotech. 16, 79–88 (2005). [CrossRef]

11. A. Li, E. L. Miller, M. E. Kilmer, T. J. Brukilacchio, T. Chaves, J. Stott, Q. Zhang, T. Wu, M. Chorlton, R. H. Moore, D. B. Kopans, and D. A. Boas, “Tomographic optical breast imaging guided by three-dimensional mammography,” Appl. Opt. 42, 5181–5190 (2003). [CrossRef] [PubMed]

12. X. Intes, J. Yu, A. Yodh, and B. Chance, “Development and evaluation of a multi-wavelength multi-channel time resolved optical instrument for NIR/MRI mammography co-registration,” in “Proceedings of the IEEE 28th Annual Northeast Bioengineering Conference,” (2002), pp. 91–92.

13. G. Gulsen, O. Birgul, M. B. Unlu, R. Shafiiha, and O. Nalcioglu, “Combined diffuse optical tomography (DOT) and MRI system for cancer imaging in small animals,” Tech. Cancer Res. Treatment 5, 351–363 (2006).

14. N. Biswal, Y. Xu, and Q. Zhu, “Imaging tumor oxyhemoglobin and deoxyhemoglobin concentrations with ultrasound-guided diffuse optical tomography.” Tech. Cancer Res. Treatment 10, 417 (2011).

15. S. van de Ven, S. Elias, A. Wiethoff, M. van der Voort, A. Leproux, T. Nielsen, B. Brendel, L. Bakker, M. van der Mark, W. Mali, and P. Luijten, “Diffuse optical tomography of the breast: initial validation in benign cysts,” Mol. Imaging Biol. 11, 64–70 (2009). [CrossRef]

16. B. W. Pogue, S. C. Davis, X. Song, B. A. Brooksby, H. Dehghani, and K. D. Paulsen, “Image analysis methods for diffuse optical tomography,” J. Biomed. Opt. 11, 33001 (2006). [CrossRef] [PubMed]

17. V. C. Kavuri, Z. J. Lin, F. Tian, and H. Liu, “Sparsity enhanced spatial resolution and depth localization in diffuse optical tomography,” Biomed. Opt. Express 3, 943–957 (2012). [CrossRef] [PubMed]

18. H. Niu, Z. J. Lin, F. Tian, S. Dhamne, and H. Liu, “Comprehensive investigation of three-dimensional diffuse optical tomography with depth compensation algorithm,” J. Biomed Opt. 15, 046005 (2010). [CrossRef] [PubMed]

19. D. Kang and M. A. Kupinski, “Signal detectability in diffusive media using phased arrays in conjunction with detector arrays,” Opt. Express 19, 12261–12274 (2011). [CrossRef] [PubMed]

20. S. P. Morgan, “Detection performance of a diffusive wave phased array,” Appl. Opt. 43, 2071–2078 (2004). [CrossRef] [PubMed]

21. Y. Chen, C. Mu, X. Intes, and B. Chance, “Signal-to-noise analysis for detection sensitivity of small absorbing heterogeneity in turbid media with single-source and dual-interfering-source,” Opt. Express 9, 212–224 (2001). [CrossRef] [PubMed]

22. S. Morgan and K. Yong, “Controlling the phase response of a diffusive wave phased array system,” Opt. Express 7, 540–546 (2000). [CrossRef] [PubMed]

23. J. P. Culver, A. M. Siegel, J. J. Stott, and D. A. Boas, “Volumetric diffuse optical tomography of brain activity,” Opt. Lett. 28, 2061–2063 (2003). [CrossRef] [PubMed]

24. H. Niu, P. Guo, X. Song, and T. Jiang, “Improving depth resolution of diffuse optical tomography with an exponential adjustment method based on maximum singular value of layered sensitivity,” Chin. Opt. Lett. 6, 886–888 (2008). [CrossRef]

25. H. H. Barrett and K. J. Myers, Foundations of Image Science (Wiley, 2004), 1st ed.

26. F. Shen and E. Clarkson, “Using Fisher information to approximate ideal-observer performance on detection tasks for lumpy-background images,” J. Opt. Soc. Am. A 23, 2406–2414 (2006). [CrossRef]

27. E. Clarkson and F. Shen, “Fisher information and surrogate figures of merit for the task-based assessment of image quality,” J. Opt. Soc. Am. A 27, 2313–2326 (2010). [CrossRef]

28. Y. Zhan, A. T. Eggebrecht, J. P. Culver, and H. Dehghani, “Image quality analysis of high-density diffuse optical tomography incorporating a subject-specific head model,” Front Neuroenergetics 4, 6 (2012). [CrossRef] [PubMed]

29. R. Ziegler, B. Brendel, A. Schipper, R. Harbers, M. v. Beek, H. Rinneberg, and T. Nielsen, “Investigation of detection limits for diffuse optical tomography systems: I. Theory and experiment,” Phys. Med. Biol. 54, 399–412 (2009). [CrossRef]

30. R. Ziegler, B. Brendel, H. Rinneberg, and T. Nielsen, “Investigation of detection limits for diffuse optical tomography systems: II. Analysis of slab and cup geometry for breast imaging,” Phys. Med. Biol. 54, 413–431 (2009). [CrossRef]

31. S. Young, M. A. Kupinski, and A. K. Jha, “Estimating signal detectability in a model diffuse optical imaging system,” in “Biomedical Optics,” (Optical Society of America, 2010), p. BSuD26.

32. A. K. Jha, M. A. Kupinski, T. Masumura, E. Clarkson, A. A. Maslov, and H. H. Barrett, “Simulating photon-transport in uniform media using the radiative transfer equation: A study using the Neumann-series approach,” J. Opt. Soc. Amer. A 29, 1741–1757 (2012). [CrossRef]

33. A. K. Jha, M. A. Kupinski, H. H. Barrett, E. Clarkson, and J. H. Hartman, “Three-dimensional Neumann-series approach to model light transport in nonuniform media,” J. Opt. Soc. Am. A 29, 1885–1899 (2012). [CrossRef]

34. V. Toronov, E. D’Amico, D. Hueber, E. Gratton, B. Barbieri, and A. Webb, “Optimization of the signal-to-noise ratio of frequency-domain instrumentation for near-infrared spectro-imaging of the human brain,” Opt. Express 11, 2717–2729 (2003). [CrossRef] [PubMed]

35. L. G. Henyey and J. L. Greenstein, “Diffuse radiation in the galaxy,” Astrophys. J. 93, 70–83 (1941). [CrossRef]

36. M. Chu, K. Vishwanath, A. D. Klose, and H. Dehghani, “Light transport in biological tissue using three-dimensional frequency-domain simplified spherical harmonics equations,” Phys. Med. Biol. 54, 2493–2509 (2009). [CrossRef] [PubMed]

37. T. Tarvainen, M. Vauhkonen, V. Kolehmainen, and J. P. Kaipio, “Hybrid radiative-transfer-diffusion model for optical tomography,” Appl. Opt. 44, 876–886 (2005). [CrossRef] [PubMed]

38. T. Spott and L. O. Svaasand, “Collimated light sources in the diffusion approximation,” Appl. Opt. 39, 6453–6465 (2000). [CrossRef]

39. Z. Q. Zhang, I. P. Jones, H. P. Schriemer, J. H. Page, D. A. Weitz, and P. Sheng, “Wave transport in random media: the ballistic to diffusive transition,” Phys. Rev. E 60, 4843–4850 (1999). [CrossRef]

40. E. Aydin, C. de Oliveira, and A. Goddard, “A finite element-spherical harmonics radiation transport model for photon migration in turbid media,” J. Quant. Spectr. Rad. Trans. 84, 247–260 (2004). [CrossRef]

41. A. Klose and E. Larsen, “Light transport in biological tissue based on the simplified spherical harmonics equations,” J. Comput. Phys. 220, 441–470 (2006). [CrossRef]

42. A. H. Hielscher and R. E. Alcouffe, “Discrete-ordinate transport simulations of light propagation in highly forward scattering heterogeneous media,” in “Advances in Optical Imaging and Photon Migration,” (Optical Society of America, 1998), p. ATuC2.

43. S. R. Arridge and W. R. B. Lionheart, “Nonuniqueness in diffusion-based optical tomography,” Opt. Lett. 23, 882–884 (1998). [CrossRef]

44. S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20, 299–309 (1993). [CrossRef] [PubMed]

45. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17, 20178–20190 (2009). [CrossRef] [PubMed]

46. M. A. Kupinski, E. Clarkson, K. Gross, and J. W. Hoppin, “Optimizing imaging hardware for estimation tasks,” in “Proc. SPIE Medical Imaging ,” (2003), pp. 309–313. [CrossRef]

47. A. K. Jha, “Retrieving Information from Scattered Photons in Medical Imaging,” Ph.D. thesis, College of Optical Sciences, University of Arizona, Tucson, AZ, USA (2013).

48. B. W. Miller, “High-Resolution Gamma-Ray Imaging with Columnar Scintillators,” Ph.D. thesis, College of Optical Sciences, University of Arizona, Tucson, AZ, USA (2011).

An ideal-observer framework to investigate signal detectability in diffuse optical imaging

Abstract

1. Introduction

2. Methods

2.1. Relation between signal detectability and Fisher information in DOT

2.2. Computing the Fisher information

2.3. Evaluating the gradient of the mean image data

2.4. Implementation

2.5. Methods for the specific SKE/BKE task

3. Experiments and results

3.1. Signal detectability as a function of depth

3.2. Signal detectability as a function of the scattering coefficient of the tissue

3.3. Signal detectability as a function of signal size and signal depth

3.4. Signal detectability as a function of signal size and signal contrast

4. Discussions and conclusions

Acknowledgments

References and links

Cited By

Figures (5)

Equations (46)

Biomedical Optics Express