## 1.

## Introduction

3-D diffuse optical tomography (DOT), as an emerging biomedical imaging modality, has made considerable advances in recent years.
^{1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27} So far DOT has been primarily applied to image brain and breast tissues.^{2, 3, 4, 5, 6} Recent phantom and clinical studies show that DOT can also provide quantitative optical images of hand joints for early detection of joint-related diseases.^{7,8,9,10,11,12,13} In DOT, an efficient reconstruction algorithm is critical. Most of the DOT algorithms developed to date make use of so-called model-based image reconstruction schemes. Due to computational complexity, the models used for describing photon migration in biological tissues and image reconstruction have usually been limited to the diffusion approximation (DA) to the Boltzmann radiative transport equation (RTE). However, the DA is not accurate enough to model light transport in some particular tissue regions, including low scattering fluid-filled and high absorption vascular tissues.^{14, 15, 16, 17, 18} In particular, using the DA is very challenging to image small volumes including arthritis in finger joints and body parts in small animals. In such cases, the DA is not able to describe photon migration accurately due to the small optical distance between sources and detectors. To overcome these limitations associated with the DA, several groups have employed the discrete ordinate
$\left({s}_{n}\right)$
approximated RTE to recover the optical parameters of small volume tissues.^{10, 11, 12} However, this type of reconstruction algorithm is time consuming. Typically, the
${s}_{n}$
approximated-RTE-based reconstruction is 40 to 60 times slower than DA based.^{11}

We have developed a 3-D reconstruction method based on simplified spherical harmonics
$\left({p}_{3}\right)$
approximated-RTE or higher order diffusion equations.^{21} Our simulations and phantom tests^{21} have demonstrated that the developed higher order diffusion (HD) equations are able to provide solutions that are significantly more accurate than the DA for simulating photon migration and imaging reconstruction in finger joints. We also found that the HD-model-based reconstruction algorithm is only 1 to 2 times slower than the DA based reconstruction.^{21}

In this study, we focus on 3-D DOT imaging of finger joints for the diagnosis of osteoarthritis (OA) using the developed HD model. One goal of this study is to investigate if the HD model has the potential to improve the reconstruction accuracy for *in-vivo* data as compared to that using the DA.^{8} In addition, we evaluate the performance of the HD model using an x-ray/DOT combined multimodality approach.^{28} The idea of this multimodality approach is to incorporate the high resolution x-ray images of joint tissues into the DOT reconstruction so that both the resolution and accuracy of optical image reconstruction are enhanced. We found from our previous investigations that the improved quantitative capability of imaging absorption and scattering coefficients of the joint tissues allows for more accurate diagnosis of osteoarthritis over the DA-based DOT alone.^{28} When anatomical *a-priori* information becomes available, we would like to study if a more accurate modeling of photon migration in tissue can bring in enhanced accuracy for image reconstruction in the finger joints.

## 2.

## Reconstruction Methods

The migration of photons within the absorbing and scattering media can be described by the time- and energy-independent form of the Boltzmana transport equation as,^{14, 15, 16, 17}

## Eq. 1

$\Omega \cdot \nabla \Phi (r,\Omega )+{\sigma}_{t}\left(r\right)\Phi (r,\Omega )=\int {\sigma}_{s}({\Omega}^{\prime},\Omega ,r)\Phi (r,{\Omega}^{\prime})d{\Omega}^{\prime}+s(r,\Omega ),$## Eq. 2

$\Phi (r,\Omega )=\sum _{l=0}^{N}\sum _{m=0}^{l}(2l+1){p}_{l}^{m}\left(\mathrm{cos}\phantom{\rule{0.2em}{0ex}}\theta \right)[{\psi}_{lm}\left(r\right)\mathrm{cos}\left(m\varphi \right)+{\gamma}_{lm}\left(r\right)\mathrm{sin}\left(m\varphi \right)],$^{21}${p}_{3}$ , ${p}_{5}$ , and ${p}_{7}$ approximations have been performed in various geometries for 2-D problems, and the results validated that the ${p}_{3}$ approximation is sufficiently accurate.

^{20}As such, only the ${p}_{3}$ approximation is used here for an efficient 3-D reconstruction. Based on Eqs. 1, 2, the Eqs. 3, 4, 5, 6, 7, 8 are the derived 3-D ${p}_{3}$ approximated-HD diffusion equations,

^{21}

## Eq. 3

$D(\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}x}+\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}y}+\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}z})-{\mu}_{a}{\psi}_{1}+D(2\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}z}-\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}x}-\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}y})+6D\left(\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}x}\right)+6D(-\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}y})+6D\frac{{\partial}^{2}{\psi}_{5}}{\partial z\partial x}+6D\frac{{\partial}^{2}{\psi}_{6}}{\partial z\partial y}+6D\frac{{\partial}^{2}{\psi}_{4}}{\partial x\partial y}+6D\frac{{\partial}^{2}{\psi}_{4}}{\partial y\partial x}=-S\left(r\right),$## Eq. 4

$\frac{25}{7}D(\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}x}+\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}y})+2D(\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}z}-\frac{1}{2}\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}x}-\frac{1}{2}\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}y})+\frac{60}{7}D(\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}y}-\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}x}-\frac{{\partial}^{2}{\psi}_{4}}{\partial x\partial y}-\frac{{\partial}^{2}{\psi}_{4}}{\partial y\partial x})+\frac{55}{7}D\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}z}+\frac{30}{7}D(\frac{{\partial}^{2}{\psi}_{5}}{\partial z\partial x}+\frac{{\partial}^{2}{\psi}_{6}}{\partial z\partial y})-5{\mu}_{t}^{\prime}{\psi}_{2}=0,$## Eq. 5

$D(\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}x}-\frac{{\partial}^{2}{\psi}_{1}}{{\partial}^{2}y})-\frac{10}{7}D(\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}x}-\frac{{\partial}^{2}{\psi}_{2}}{{\partial}^{2}y})+\frac{30}{7}D(\frac{{\partial}^{2}{\psi}_{5}}{\partial x\partial z}-\frac{{\partial}^{2}{\psi}_{6}}{\partial y\partial z})+\frac{90}{7}D(\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}x}+\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}y}+\frac{1}{3}\frac{{\partial}^{2}{\psi}_{3}}{{\partial}^{2}z})-10{\mu}_{t}^{\prime}{\psi}_{3}=0,$## Eq. 6

$\frac{1}{2}D\frac{{\partial}^{2}{\psi}_{1}}{\partial x\partial y}+\frac{1}{2}D\frac{{\partial}^{2}{\psi}_{1}}{\partial y\partial x}-\frac{5}{7}D\frac{{\partial}^{2}{\psi}_{2}}{\partial x\partial y}-\frac{5}{7}D\frac{{\partial}^{2}{\psi}_{2}}{\partial y\partial x}+\frac{45}{7}D(\frac{{\partial}^{2}{\psi}_{4}}{{\partial}^{2}x}+\frac{{\partial}^{2}{\psi}_{4}}{{\partial}^{2}y}+\frac{1}{3}\frac{{\partial}^{2}{\psi}_{4}}{{\partial}^{2}z})+\frac{15}{7}D(\frac{{\partial}^{2}{\psi}_{5}}{\partial y\partial z}+\frac{{\partial}^{2}{\psi}_{6}}{\partial x\partial z})-5{\mu}_{t}^{\prime}{\psi}_{4}=0,$## Eq. 7

$2D\frac{{\partial}^{2}{\psi}_{1}}{\partial x\partial z}+\frac{10}{7}D\frac{{\partial}^{2}{\psi}_{2}}{\partial x\partial z}+\frac{45}{7}D(\frac{{\partial}^{2}{\psi}_{5}}{{\partial}^{2}z}+\frac{{\partial}^{2}{\psi}_{5}}{{\partial}^{2}x}+\frac{1}{3}\frac{{\partial}^{2}{\psi}_{5}}{{\partial}^{2}y})+\frac{60}{7}D(\frac{{\partial}^{2}{\psi}_{3}}{\partial x\partial z}+\frac{{\partial}^{2}{\psi}_{4}}{\partial y\partial z}+\frac{1}{2}\frac{{\partial}^{2}{\psi}_{6}}{\partial z\partial y})-5{\mu}_{t}^{\prime}{\psi}_{5}=0,$## Eq. 8

$2D\frac{{\partial}^{2}{\psi}_{1}}{\partial y\partial z}+\frac{10}{7}D\frac{{\partial}^{2}{\psi}_{2}}{\partial y\partial z}+\frac{45}{7}D(\frac{{\partial}^{2}{\psi}_{6}}{{\partial}^{2}z}+\frac{1}{3}\frac{{\partial}^{2}{\psi}_{6}}{{\partial}^{2}x}+\frac{{\partial}^{2}{\psi}_{6}}{{\partial}^{2}y})+\frac{60}{7}D(\frac{{\partial}^{2}{\psi}_{4}}{\partial x\partial z}-\frac{{\partial}^{2}{\psi}_{3}}{\partial y\partial z}+\frac{1}{2}\frac{{\partial}^{2}{\psi}_{5}}{\partial x\partial y})-5{\mu}_{t}^{\prime}{\psi}_{6}=0,$## Eq. 9

$\nabla \cdot D\left(r\right)\nabla {\psi}_{1}\left(r\right)-{\mu}_{a}\left(r\right){\psi}_{1}\left(r\right)=-S\left(r\right).$## Eq. 10

$-D\nabla {\psi}_{1}\cdot n=\alpha {\psi}_{1},\phantom{\rule{1em}{0ex}}{\psi}_{2}={\psi}_{3}={\psi}_{4}={\psi}_{5}={\psi}_{6}=0,$^{21, 25}Generally, type 3 BCs are applied in diffuse optical imaging and contain the linear combination of the photon density and the current along the boundary, while type 1 BCs are the commonly used boundary conditions for the RTE model. The primary reason here for using type 3 BCs for the first component and type 1 BCs for other components is for mathematical simplicity. In fact, type 3 BCs can also be applied to components 2 through 6, but we believe this may not add any benefit, as the first component is the measured quantity for image reconstruction. In finite element discretization, the quantities $\psi \left(r\right)$ , ${\mu}_{t}^{\prime}$ , $D$ , and ${\mu}_{a}$ are spatially discretized as,

## Eq. 11

${\psi}_{k}=\sum _{i=1}^{N}{\varphi}_{i}{\left({\psi}_{k}\right)}_{i},\phantom{\rule{1em}{0ex}}{\mu}_{t}^{\prime}=\sum _{i=1}^{N}{\varphi}_{i}{\left({\mu}_{t}^{\prime}\right)}_{i},\phantom{\rule{1em}{0ex}}D=\sum _{i=1}^{N}{D}_{i}{\varphi}_{i},{\mu}_{a}=\sum _{i=1}^{N}{\left({\mu}_{a}\right)}_{i}{\varphi}_{i},$**A**] are integrated over the problem domain (

**V**) and boundary domain (Γ), and {

**b**} is the source vector. The inverse solution can be obtained through the following equations:

## Eq. 13

$\left[\mathbf{A}\right]\{\partial \psi \u2215\partial \chi \}=\{\partial \mathbf{b}\u2215\partial \chi \}-[\partial \mathbf{A}\u2215\partial \chi ]\left\{\psi \right\},$## Eq. 14

$({\mathfrak{I}}^{T}\mathfrak{I}+\lambda \mathbf{I})\Delta \mathbf{\chi}={\mathfrak{I}}^{T}({\mathbf{\psi}}^{o}-{\mathbf{\psi}}^{c}),$**I**, the identity matrix, are used to realize the invertible system [Eq. 14]; $\Delta \chi ={(\Delta {D}_{1},\Delta {D}_{2},\dots ,\Delta {D}_{N},\Delta {\mu}_{a,1},\Delta {\mu}_{a,2},\dots ,\Delta {\mu}_{a,N})}^{T}$ ; ${\psi}^{o}={({\psi}_{1}^{1o},{\psi}_{1}^{2o},\dots ,{\psi}_{1}^{Mo})}^{T}$ , and ${\psi}^{c}={({\psi}_{1}^{1c},{\psi}_{1}^{2c},\dots ,{\psi}_{1}^{Mc})}^{T}$ , where ${\psi}_{1}^{o}$ and ${\psi}_{1}^{c}$ are observed and computed photon density (the first component) for $i=1,2\dots ,M$ boundary locations. It should be noted that only the first component, as the measurements, is used for the inverse computation due to the fact that the other five moment variables are equal to zero at the boundary measurement positions.

For x-ray guided 3-D DOT reconstruction, the following updating equation is used,

## Eq. 15

$\Delta \mathbf{\chi}={({\mathfrak{I}}^{T}\mathfrak{I}+{\mathfrak{I}}^{T}\mathfrak{I}+\lambda \mathbf{I}+{\mathbf{L}}^{T}\mathbf{L})}^{-1}\left[{\mathfrak{I}}^{T}({\mathbf{\psi}}^{o}-{\mathbf{\psi}}^{c})\right],$*a-priori*information is incorporated into the iterative process using the spatially variant filter matrix

**L**. In this study, the Laplacian-type filter matrix was used, and its elements ${\mathbf{L}}_{ij}$ were constructed according to the visible region or tissue type it was associated with, with x-ray derived priors as follows:

^{28, 29, 30}

## Eq. 16

${L}_{ij}=\{\begin{array}{ccc}1& \text{when}\phantom{\rule{0.3em}{0ex}}i=j& \\ -1\u2215nn& \text{when}\phantom{\rule{0.3em}{0ex}}i,j\subset \text{one}\phantom{\rule{0.3em}{0ex}}\text{region}\\ 0& \text{when}\phantom{\rule{0.3em}{0ex}}i,j\subset \text{different}\phantom{\rule{0.3em}{0ex}}\text{region}\end{array}\phantom{\}},$^{18}Thus the final updating equation is modified as follows:

## Eq. 17

${\mathbf{\chi}}_{\text{new}}={\mathbf{\chi}}_{\text{old}}+\zeta \Delta \mathbf{\chi}(0<\zeta \u2a7d1),$## 3.

## Results

In this section, we present reconstruction results using the HD model with phantom and *in-vivo* data. For comparison, reconstructed optical images using the DA-based algorithm are also provided.

Both phantom and *in-vivo* data were collected using our homemade 64 × 64-channel photodiodes-based stand-alone DOT system or the hybrid x-ray/DOT imaging system that integrates a modified mini-C-arm x-ray system with the 64 × 64-channel DOT system, which were all detailed elsewhere.^{8, 22, 28} This DOT system consists of several laser modules, a hybrid light delivery subsystem, a fiber optics/tissue interface, a data acquisition module, and light detection modules.^{22} The cylindrical fiber optics/tissue interface is composed of 64 source and 64 detector fiber bundles that are positioned in four layers along the surface of a plexiglass container and cover a volume of 15 (height) by 30 mm (diameter). In each layer, 16 source and 16 detector fiber bundles are alternatively arranged. The space between the finger and the wall of the plexiglass container is filled with tissue-like phantom materials as coupling media, consisting of distilled water, agar powder, Indian ink, and Intralipid, giving an absorption coefficient of
$0.014\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$
and a reduced scattering coefficient of
$1.0\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$
.^{8, 22} For *in-vivo* study, we focused on imaging the distal interphalangeal (DIP) finger joints from healthy volunteers and OA patients. The study was approved by the Institutional Review Board at the University of Florida. In the optical exam, the index finger of the participant was inserted into the 30-mm-diam cylindrical solid phantom medium.^{8}

For DOT-alone imaging, we just present the *in-vivo* results as we have previously reported the phantom study using both the HD and DA models.^{21} Reconstructions were performed using a 3-D mesh with 3009 nodes and 12,800 tetrahedral elements. The mesh used is uniform in the height direction and the average distance of nodes is actually about 1.25 mm. In addition, in the region around the joint cavity/spacing, the mesh was discretized in terms of number of nodes per
${\mathrm{mm}}^{3}$
to ensure that the average distance of nodes is around 1 mm in the joint spacing. The relatively coarse mesh density generated inside the bone domain has no significant influence on our reconstruction accuracy, which has been demonstrated and validated by our previous simulation tests and phantom experiments.^{8, 18, 21, 28} As such, the mesh size used is good enough for diffuse optical imaging that only has a resolution limit of 2 to 3 mm. For the imaging systems, 64 sources and 64 detectors were distributed uniformly along the surface of the phantom at four planes (
$z=2.5$
,
$z=7.5$
,
$z=12.5$
, and
$z=17.5\phantom{\rule{0.3em}{0ex}}\mathrm{mm}$
; 16 sources and 16 detectors at each plane). The initial guesses used are
${\mu}_{a}=0.03$
and
${\mu}_{s}^{\prime}=1.0$
for all the tests. For the reconstruction using the HD model, the scattering anisotropic factor
$g=0.7$
. Figures 1a and 1b and Figs. 2a and 2b provide reconstructed two selected coronal (
$y\text{-}z$
) slices for a healthy and OA joint, respectively. The recovered quantitative optical properties are also given in Table 1
. We note that the results provided for both DA and HD models required six iterations for the inverse computation. If we use a LU decomposition-based forward solution solver, the DA and HD models require, respectively, 9.0 and 15.0 min per iteration when a FE mesh of 3009 nodes with 12,800 tetrahedral elements is employed in a high performance PC. If a gradient-based solver is used, the computation time can be reduced to 12 min per iteration for the HD model.

## Table 1

Averaged values of absorption and scattering coefficients of bone and joint tissues using the HD and DA models without x-ray guidance.

Test | μa (mm−1) (bone) | μa (mm−1) (joint) | μs′ (mm−1) (bone) | μs′ (mm−1) (joint) | ||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|

HD | DA | Difference | HD | DA | Difference | HD | DA | Difference | HD | DA | Difference | |

OA | 0.049 | 0.053 | 6.3% | 0.041 | 0.0452 | 10% | 2.20 | 2.30 | 4.3% | 1.85 | 1.94 | 4.6% |

Healthy | 0.048 | 0.051 | 6.2% | 0.029 | 0.0337 | 15% | 1.91 | 1.96 | 2.6% | 1.08 | 1.16 | 6.9% |

In the phantom experiments with x-ray guidance, a 30-mm-diam cylindrical solid phantom was used as the background medium. Two 20-mm-diam cylindrical solid objects (3-mm off $z$ -axis) mimicking bones were embedded in the background medium. The spacing (“cartilage”) between the two “bones” was 2.5 mm. The optical properties for the background were ${\mu}_{a}=0.01\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$ and ${\mu}_{s}^{\prime}=1.0\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$ , while the optical properties of the “cartilage” were assumed to be the same as the background medium. The optical properties for the mimicking bones are ${\mu}_{a}=0.07\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$ and ${\mu}_{s}^{\prime}=4.0\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$ . The recovered images using the DA (left column) and HD equations (right column) at selected longitudinal planes are provided in Figs. 3a and 3b .

For *in-vivo* imaging with x-ray guidance, image reconstruction of the DIP finger joints with the optical coupling phantom/media (30 × 20 mm in volume) were performed with a finite element mesh of 2705 nodes and 13,440 tetrahedral elements. This generated mesh is not uniform—the prior x-ray image of bones allows us to use a nodal distance of ∼0.4 mm in the joint spacing, which is accurate enough for joint space imaging. In particular, we are interested in imaging the optical properties of the joint region, and the fine mesh around the cavity is sufficient for our study. Again, the coarse mesh density for the bones has no significant effect on the reconstruction accuracy. The initial optical properties used were optimized based on an x-ray guided forward-fitting algorithm for the approximated joint and bone regions/tissues.^{24} The optical properties for the background were
${\mu}_{a}=0.014\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$
and
${\mu}_{s}^{\prime}=1.0\phantom{\rule{0.3em}{0ex}}{\mathrm{mm}}^{-1}$
, which also served as initial guesses for the background medium. Figures 4a and 4b plot the selected absorption and scattering slices of the recovered 3-D images for an OA joint using the DA and HD models, while Figs. 5a and 5b display the reconstructed slices for a healthy joint. The optical absorption and scattering coefficient values are provided in Table 2 for both the phantom and *in-vivo* data.

## Table 2

Averaged values of absorption and scattering coefficients of bone and joint tissues using the HD and DA models with x-ray guidance. Note that difference = (DA-HD) × 100%/DA.

Test | μa (mm−1) (bone) | μa (mm−1) (joint) | μs′ (mm−1) (bone) | μs′ (mm−1) (joint) | ||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|

HD | DA | Difference | HD | DA | Difference | HD | DA | Difference | HD | DA | Difference | |

Phantom | 0.056 | 0.054 | − 3.8% | 0.0098 | 0.0095 | − 3.2% | 3.51 | 3.53 | 0.6% | 1.003 | 1.02 | 1.3% |

OA | 0.046 | 0.047 | 1.9% | 0.023 | 0.0230 | 2.1% | 2.25 | 2.27 | 0.7% | 1.527 | 1.54 | 0.7% |

Healthy | 0.040 | 0.041 | 2.3% | 0.0098 | 0.0100 | 2.0% | 2.15 | 2.16 | 0.6% | 0.952 | 0.97 | 1.4% |

## 4.

## Discussion

From the reconstructed optical images shown in Figs. 1 and 2, we note that the bones are clearly delineated for both the normal and OA joints using DA and HD models. Importantly, compared with the optical property value of the bones, we observe a large drop in the magnitude of absorption and scattering coefficients within the joint space for the healthy case, while only a small drop is seen for the OA joint. These observations also agree well with the x-ray guided optical findings, as displayed in Figs. 4 and 5. In particular, the image resolution from the hybrid imaging method is significantly improved compared to that from DOT alone. We see that x-ray guided optical images exhibit apparent narrowing of the joint space for the OA joint (Fig. 4) relative to those for the healthy joint (Fig. 5).

It is noted from Figs. 1 and 2 that model errors do exist when the DA is used to handle cases involving low scattering, small source/detector distances or large variations of optical properties. As shown in Table 1, the difference in recovered quantitative optical properties between the two models can be as large as 15% (for the absorption coefficient of joint tissue). This demonstrates that the
$p3$
approximated HD model can produce significantly improved reconstruction of both absorption and scattering images for joint imaging. In addition, it should also be pointed out that while a transport or HD model can certainly reduce a large amount of errors in forward computation,^{21} the model errors existing in the forward calculation do not appear to propagate much into the image reconstruction. This is understandable, because the accuracy of an inverse solution depends not only on the accuracy of the forward model, but also on the quality of the experimental data and the use of robust regularization techniques.

More interestingly, model errors appear to lead to even smaller error in reconstruction when the x-ray *a-priori* structural information is incorporated into DOT reconstruction, as shown in Figs. 3, 4, 5. The difference for the recovered optical properties is less than 4% between the two models, which is evident from Table 2. This is due to the fact that a more accurate modeling of photon migration in tissue can be achieved using the DA when anatomical *a-priori* information becomes available. The use of prior structural information eliminates the need to look for spatial/anatomy information in the optical inversion, which ensures that optical property profiles are the only parameter(s) that need to be recovered.

In summary, we demonstrate in this study that the improvement resulted from HD-based reconstruction is significant over the DA model for DOT-alone imaging of finger joints. The improvement is moderate using the HD model when x-ray guidance is used because of the high quality of reconstruction already available from the DA model under such a situation.

## Acknowledgments

This research was supported in part by the National Institutes of Health (R01 AR048122).

## References

**,” Phys. Today, 48 34 –40 (1995). https://doi.org/10.1063/1.881445 0031-9228 Google Scholar**

*Spectroscopy and imaging with diffusing light***,” Neoplasia, 2 26 –40 (2000). https://doi.org/10.1038/sj.neo.7900082 1522-8002 Google Scholar**

*Noninvasive in vivo characterization of breast cancer using photon migration spectroscopy***,” J. Biomed. Opt., 10 024033 (2005). https://doi.org/10.1117/1.1899183 1083-3668 Google Scholar**

*Coregistered tomography x-ray and optical breast imaging: initial results***,” IS11 87 –120 (1993) Google Scholar**

*A perturbation approach for optical diffusion tomography using continuous-wave and time-resolved data***,” Phys. Med. Biol., 42 825 –840 (1997). https://doi.org/10.1088/0031-9155/42/5/007 0031-9155 Google Scholar**

*Optical imaging in medicine: I: experimental techniques***,” Phys. Med. Biol., 44 2897 –2906 (1999). https://doi.org/10.1088/0031-9155/44/12/303 0031-9155 Google Scholar**

*Photon migration in non-scattering tissue and the effects on image reconstruction***,” J. Biomed. Opt., 9 (3), 474 –480 (2004). https://doi.org/10.1117/1.1691029 1083-3668 Google Scholar**

*Optical biopsy of bone tissue: a step toward the diagnosis of bone pathologies***,” J. Biomed. Opt., 12 034001 (2007). https://doi.org/10.1117/1.2737420 1083-3668 Google Scholar**

*3D diffuse optical tomography imaging of osteoarthritis: Initial results in finger joints***,” J. Biomed. Opt., 7 (1), 88 –92 (2002). https://doi.org/10.1117/1.1427336 1083-3668 Google Scholar**

*Three-dimensional diffuse optical tomography of bones and joints***,” Phys. Med. Biol., 49 1147 –1163 (2004). https://doi.org/10.1088/0031-9155/49/7/005 0031-9155 Google Scholar**

*Sagittal laser optical tomography for imaging of rheumatoid finger joints***,” Appl. Opt., 46 6669 –6679 (2007). https://doi.org/10.1364/AO.46.006669 0003-6935 Google Scholar**

*Transport and diffuse-based optical tomography in small volumes: a comparative study***,” J. Opt. Soc. Am. A, 20 2355 –2364 (2003). https://doi.org/10.1364/JOSAA.20.002355 0740-3232 Google Scholar**

*Generalized diffuse model in optical tomography with clear layer***,” Biomed. Tech., 42 319 –326 (1997). 0013-5585 Google Scholar**

*Development of a finger joint phantom for the optical simulation of early stages of rheumatoid arthritis***,” J. Comp. Physiol., 220 441 –470 (2006). 0373-0859 Google Scholar**

*Light transport in biological tissue based on the simplified spherical harmonics equations***,” Med. Phys., 29 2013 –2023 (2002). https://doi.org/10.1118/1.1500404 0094-2405 Google Scholar**

*A comparison between transport and diffusion calculations using a finite element-spherical harmonics radiation transport method***,” Opt. Lett., 26 1066 –1068 (2001). https://doi.org/10.1364/OL.26.001066 0146-9592 Google Scholar**

*Photon transport forward model for imaging in turbid media***,” Inverse Probl., 14 1107 –1130 (1998). https://doi.org/10.1088/0266-5611/14/5/003 0266-5611 Google Scholar**

*Transport-backtransport method for optical tomography***,” Appl. Opt., 46 2757 –2768 (2007). https://doi.org/10.1364/AO.46.002757 0003-6935 Google Scholar**

*Image reconstruction schemes that combines modified Newton method and efficient initial guess estimate for optical tomography of finger joints***,” J. Opt. Soc. Am. A, 20 92 –98 (2003). https://doi.org/10.1364/JOSAA.20.000092 0740-3232 Google Scholar**

*Light propagation in biological tissues***,” Proc. SPIE, 2389 240 –247 (1995). https://doi.org/10.1117/12.209973 0277-786X Google Scholar**

*Photon migration within the P3 approximation***,” Phys. Med. Biol., 54 65 –88 (2009). https://doi.org/10.1088/0031-9155/54/1/005 0031-9155 Google Scholar**

*A higher-order diffusion model for three-dimensional photon migration and image reconstruction in optical tomography***,” J. Opt. A, Pure Appl. Opt., 7 224 –231 (2005). https://doi.org/10.1088/1464-4258/7/5/003 1464-4258 Google Scholar**

*Three-dimensional diffuse optical tomography of simulated hand joints with a 64 × 64-channel photodiodes-based optical system***,” Appl. Opt., 35 3447 –3458 (1996). https://doi.org/10.1364/AO.35.003447 0003-6935 Google Scholar**

*Enhanced frequency-domain optical image reconstruction in tissues through total-variation minimization***,” Appl. Opt., 39 5256 –5261 (2000). https://doi.org/10.1364/AO.39.005256 0003-6935 Google Scholar**

*Quantitative optical image reconstruction of turbid media by use of direct-current measurements***,” J. Opt. Soc. Am. A, 13 253 –266 (1996). https://doi.org/10.1364/JOSAA.13.000253 0740-3232 Google Scholar**

*Optical image reconstruction using frequency-domain data: simulations and experiments***,” Med. Phys., 22 691 –701 (1995). https://doi.org/10.1118/1.597488 0094-2405 Google Scholar**

*Spatially-varying optical property reconstruction using a finite element diffusion equation approximation***,” Opt. Express, 4 241 –246 (1999). https://doi.org/10.1364/OE.4.000241 1094-4087 Google Scholar**

*Optical image reconstruction based on the third-order diffusion equations***,” J. Biomed. Opt., 13 044006 (2008). https://doi.org/10.1117/1.2965547 1083-3668 Google Scholar**

*Tomographic x-ray-guided three-dimensional diffuse optical tomography of osteoarthritis in the finger joints***,” J. Biomed. Opt., 10 011504 (2005). https://doi.org/10.1117/1.2098627 1083-3668 Google Scholar**

*Combining near-infrared tomography and magnetic resonance imaging to study in vivo breast tissues: implementation of a Laplacian-type regularization to incorporate magnetic resonance structure***,” Appl. Opt., 38 2950 –2961 (1996). https://doi.org/10.1364/AO.38.002950 0003-6935 Google Scholar**

*Spatially variant regularization improves diffuse optical tomography*