## A least-squares fixed-point iterative algorithm for multiple illumination photoacoustic tomography |

Biomedical Optics Express, Vol. 4, Issue 10, pp. 2224-2230 (2013)

http://dx.doi.org/10.1364/BOE.4.002224

Acrobat PDF (1012 KB)

### Abstract

The optical absorption of tissues provides important information for clinical and pre-clinical studies. The challenge in recovering optical absorption from photoacoustic images is that the measured pressure depends on absorption and local fluence. One reconstruction approach uses a fixed-point iterative technique based on minimizing the mean-squared error combined with modeling of the light source to determine optical absorption. With this technique, convergence is not guaranteed even with an accurate measure of optical scattering. In this work we demonstrate using simulations that a new multiple illumination least squares fixed-point iteration algorithm improves convergence - even with poor estimates of optical scattering.

© 2013 OSA

## 1. Introduction

6. J. Ripoll and V. Ntziachristos, “Quantitative point source photoacoustic inversion formulas for scattering and absorbing media,” Phys. Rev. E **71**, 031912 (2005). [CrossRef]

8. B. Banerjee, S. Bagchi, R. M. Vasu, and D. Roy, “Quantitative photoacoustic tomography from boundary pressure measurements: noniterative recovery of optical absorption coefficient from the reconstructed absorbed energy map,” J. Opt. Soc. Am. A **25**, 2347–2356 (2008). [CrossRef]

9. B. T. Cox, S. R. Arridge, K. P. Köstli, and P. C. Beard, “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt. **45**, 1866–1875 (2006). [CrossRef] [PubMed]

10. L. Yin, Q. Wang, Q. Zhang, and H. Jiang, “Tomographic imaging of absolute optical absorption coefficient in turbid media using combined photoacoustic and diffusing light measurements,” Opt. Lett. **32**, 2556–2558 (2007). [CrossRef] [PubMed]

11. T. Jetzfellner, D. Razansky, A. Rosenthal, R. Schulz, K. H. Englmeier, and V. Ntziachristos, “Performance of iterative optoacoustic tomography with experimental data,” Appl. Phys. Lett. **95**, 013703 (2009). [CrossRef]

12. B. T. Cox, S. R. Arridge, and P. C. Beard, “Estimating chromophore distributions from multiwavelength photoacoustic images,” J. Opt. Soc. Am. A **26**, 443–455 (2009). [CrossRef]

13. G. Bal and K. Ren, “On multi-spectral quantitative photoacoustic tomography in diffusive regime,” Inverse Probl. **28**, 025010 (2012). [CrossRef]

14. G. Bal and G. Uhlmann, “Inverse diffusion theory of photoacoustics,” Inverse Probl. **26**, 085010 (2010). [CrossRef]

15. R. J. Zemp, “Quantitative photoacoustic tomography with multiple optical sources,” Appl. Opt. **49**, 3566–3572 (2010). [CrossRef] [PubMed]

16. P. Shao, B. Cox, and R. J. Zemp, “Estimating optical absorption, scattering, and grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt. **50**, 3145–3154 (2011). [CrossRef] [PubMed]

15. R. J. Zemp, “Quantitative photoacoustic tomography with multiple optical sources,” Appl. Opt. **49**, 3566–3572 (2010). [CrossRef] [PubMed]

16. P. Shao, B. Cox, and R. J. Zemp, “Estimating optical absorption, scattering, and grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt. **50**, 3145–3154 (2011). [CrossRef] [PubMed]

17. P. Shao, T. Harrison, and R. J. Zemp, “Iterative algorithm for multiple illumination photoacoustic tomography (mipat) using ultrasound channel data,” Biomed. Opt. Express **3**, 3240–3249 (2012). [CrossRef] [PubMed]

4. B. Cox, J. G. Laufer, S. R. Arridge, and P. C. Beard, “Quantitative spectroscopic photoacoustic imaging: a review,” J. Biomed. Opt. **17**, 061202 (2012). [CrossRef] [PubMed]

18. H. Gao, S. Osher, and H. Zhao, “Quantitative photoacoustic tomography,” in “*Mathematical Modeling in Biomedical Imaging II: Optical, Ultrasound, and Opto-Acoustic Tomographiess*,” vol. 2035 of *Lecture Notes in Mathematics: Mathematical Biosciences Subseries*, H. Ammari, ed. (Springer-Verlag, Berlin, 2011), pp. 131–158.

19. K. Ren, H. Gao, and H. Zhao, “A Hybrid Reconstruction Method for Quantitative PAT,” SIAM J. Imaging Sci. **6**, 32–55 (2013). [CrossRef]

9. B. T. Cox, S. R. Arridge, K. P. Köstli, and P. C. Beard, “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt. **45**, 1866–1875 (2006). [CrossRef] [PubMed]

11. T. Jetzfellner, D. Razansky, A. Rosenthal, R. Schulz, K. H. Englmeier, and V. Ntziachristos, “Performance of iterative optoacoustic tomography with experimental data,” Appl. Phys. Lett. **95**, 013703 (2009). [CrossRef]

## 2. Theory

### 2.1. Single illumination

*p*

_{0}has already been reconstructed using any of the various techniques available [1

1. M. Xu and L. V. Wang, “Analytic explanation of spatial resolution related to bandwidth and detector aperture size in thermoacoustic or photoacoustic reconstruction,” Phys. Rev. E **67**, 056605 (2003). [CrossRef]

3. Z. Guo, C. Li, L. Song, and L. V. Wang, “Compressed sensing in photoacoustic tomography in vivo,” J. Biomed. Opt. **15**, 021311 (2010). [CrossRef] [PubMed]

*μ*by modeling the fluence distribution Φ. The initial pressure distribution is actually a combination of

_{a}*μ*, Φ and the Grüneisen parameter Γ which for simplicity we assume is uniform and constant.

_{a}*p*

_{0}takes the form in equation 1.

*r*,

*μ*,

_{a}*μ′*): a function varying over spatial position

_{s}**, absorption**

*r**μ*(

_{a}**), and the reduced scattering coefficient**

*r**μ′*(

_{s}**). For the simple iterative technique we are interested in [11**

*r*11. T. Jetzfellner, D. Razansky, A. Rosenthal, R. Schulz, K. H. Englmeier, and V. Ntziachristos, “Performance of iterative optoacoustic tomography with experimental data,” Appl. Phys. Lett. **95**, 013703 (2009). [CrossRef]

*μ′*is taken to be uniform. As in that work, we assume the diffusion equation is applicable, and use the finite element method solver available in the TOAST toolkit [20

_{s}20. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. **15**, R41 (1999). [CrossRef]

*μ̂*= 0, the iterative method of [9

_{a}9. B. T. Cox, S. R. Arridge, K. P. Köstli, and P. C. Beard, “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt. **45**, 1866–1875 (2006). [CrossRef] [PubMed]

*i*: use the forward model to estimate Φ

*; calculate an error parameter*

^{i}*σ*is a regularization parameter); and repeat with

*i*=

*i*+ 1 until the error is within acceptable bounds, or the solution for

*μ̂*has converged.

_{a}### 2.2. Extension to multiple illuminations

*N*×

*N*image due to illuminations

*k*= 1,...,

*S*. The reconstructed initial pressures can be modeled as

*(*

_{k}**) is the estimated fluence due to illumination**

*r**k*, and

*μ̂*(

_{a}**) is the estimated optical absorption coefficient distribution. For simplicity, Γ is considered to be spatially constant. From the**

*r**M*=

*N*

^{2}pixels of each of the

*S*reconstructed images, one can form a vector of the observation data

*p̂*_{0}=

**, where**

*Âμ̂*_{a}**= [**

*μ̂*_{a}*μ̂*(

_{a}

*r*_{1}),...,

*μ̂*(

_{a}

*r**)]*

_{M}*is an*

^{T}*M*× 1 column vector of estimated optical absorption coefficients, and where

**= [**

*Â*

*Â*_{1},...,

*Â**]*

_{S}*with*

^{T}

*Â**= Γ × diag(*

_{k}**Φ**̂

*), where*

_{k}**Φ**̂

*= [Φ̂*

_{k}*(*

_{k}

*r*_{1}),..., Φ̂

*(*

_{k}

*r**)]*

_{M}*.*

^{T}**such that the error between the observations and the computed images**

*μ̂*_{a}*ε*(

**) = ||**

*μ̂*_{a}

*p̂*_{0}−

**||**

*Âμ̂*_{a}^{2}is a minimum. The least squares solution to this problem comes from solving (

*Â*

^{T}**)**

*Â*

*μ̂**=*

_{a}

*Â*

^{T}

*p̂*_{0}for the vector

*μ̂**. Here*

_{a}**= (**

*μ̂*_{a}

*Â*

^{T}**)**

*Â*^{−1}

*Â*

^{T}

*p̂*_{0}, we have

_{k}^{(0)}(

**,**

*r**β. β*is intended only to ensure long-term convergence of the algorithm which may otherwise diverge due to noise effects, or in the case of experimental work, non-idealities. This results in Eq. (3) for iteration

*i*+ 1. For this equation, the absorption coefficients are guaranteed to be non-negative given that reconstructed initial pressures are themselves non-negative.

## 3. Simulations

20. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. **15**, R41 (1999). [CrossRef]

**95**, 013703 (2009). [CrossRef]

*from the forward simulation), and*

_{k}*μ̂*

_{a}^{(i)}, the current quantitative image. This is the same measure as the second quality measure from the work of Jetzfellner et al. [11

**95**, 013703 (2009). [CrossRef]

**95**, 013703 (2009). [CrossRef]

*μ′*= 20cm

_{s}^{−1},

*μ*= 1.5cm

_{a}^{−1}for the inclusion, and

*μ*= 0.2cm

_{a}^{−1}for the main body of the phantom. With

*σ*= 0.001, we get very similar results as that work in terms of convergence with the varying estimates of

*μ′*used for reconstruction. The most concerning aspect of the previous results is that the reconstruction does not converge even with the correct, uniform

_{s}*μ′*= 20cm

_{s}^{−1}. Figure 1(b) shows the initial pressure estimate and true fluence for a single uniform illumination. The first iteration of the algorithm gives a reasonable but inaccurate result seen in Fig. 1(c), but after 30 iterations, the solution is clearly diverging as in Fig. 1(d). We can instead use the four images with resulting fluences given in 1(e). By applying the MIPAT technique detailed above, setting

*β*

^{2}=

*σ*, and keeping all other parameters the same, we obtain very similar results for the first iteration in Fig. 1(f) (this is expected since it is the case with an initial guess of a uniform

*μ*), and a much improved estimate of

_{a}*μ*after 30 iterations in Fig. 1(g).

_{a}*β*plays an important role in ensuring convergence in this case, so we explored the problem space in two dimensions: both over signal-to-noise ratio (SNR) and

*β*. To make a fair comparison between different MIPAT illumination patterns, for each of

*S*illuminations the total optical power delivered was

*β*and

*σ*. Figure 3 shows our results with 4 and 16 illuminations.

## 4. Discussion and conclusions

*μ′*, a simple iterative technique does not necessarily result in a convergent solution. Our simulated results are in very good agreement with experiments from Jetzfellner et al. [11

_{s}**95**, 013703 (2009). [CrossRef]

16. P. Shao, B. Cox, and R. J. Zemp, “Estimating optical absorption, scattering, and grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt. **50**, 3145–3154 (2011). [CrossRef] [PubMed]

*μ′*is not well-known. While the speed of convergence is similar for the cases presented here, this figure shows that in all cases, the root-mean-squared error is minimized with the correct estimate for

_{s}*μ′*as one would expect. Figure 2 illustrates an interesting behavior of the numerical method: for all values of

_{s}*μ′*used in the reconstruction, a local minimum appears to be reached after two or three iterations before the error measure increases to eventually either diverge or converge to a different value. Indeed, fixed-point iteration numerical methods are not guaranteed to converge at a global minimum, and in fact are not guaranteed to converge at all.

_{s}*β*regularization parameter. Figure 3 demonstrates the effects of SNR and

*β*for 4 and 16 illuminations. In general, the trend that can be seen from these simulations is that increasing

*β*can improve convergence in noisy situations, but that comes at the expense of reconstruction quality. The results seen in the figure are fairly intuitive, showing that low SNR hurts reconstruction quality, as does increasing

*β*. However, there seems to be an ideal

*β*for this set of simulations around 0.01 that provides reasonably accurate convergent solutions in low SNR conditions, while not significantly impacting reconstruction in high SNR conditions. While Fig. 3 demonstrates a range of SNR values, low SNR will be more typical of practical systems. In this case, more illuminations and larger regularization parameters are shown to be helpful. The choice of the regularization level using experimental data should be a topic of future work. A value of

*β*that is too high will tend to cause the algorithm to converge in as little as a single iteration while providing inaccurate reconstruction. The choice of these regularization parameters (

*β*for MIPAT and

*σ*for the single illumination case) is certainly non-trivial, and in this work we were guided by the experimental work by Jetzfellner et al. [11

**95**, 013703 (2009). [CrossRef]

*β*that produced stable, yet accurate results is on the same order of magnitude as the

*σ*= 0.001 used in that work if we use

*σ*=

*β*

^{2}.

*S*sources each having

21. B. Cox, T. Tarvainen, and S. Arridge, “Multiple illumination quantitative photoacoustic tomography using transport and diffusion models,” in “*Tomography and Inverse Transport Theory*,” G. Bal, D. Finch, J. Schotland, P. Kuchment, and P. Stefanov, eds. (American Mathematical Society, Providence, RI, USA, 2012), pp. 1–12.

**95**, 013703 (2009). [CrossRef]

*μ′*is not well known. Not only that, but our investigation has shown that MIPAT can be used even in realistic noisy images by appropriately selecting the regularization parameter. While this technique does require a somewhat more complicated experimental setup than a single illumination method, the benefits of applying multiple illuminations may outweigh the drawbacks. Compared to previous multiple illumination algorithms, the fixed-point method discussed here does not require inversion of large Jacobian matrices [17

_{s}17. P. Shao, T. Harrison, and R. J. Zemp, “Iterative algorithm for multiple illumination photoacoustic tomography (mipat) using ultrasound channel data,” Biomed. Opt. Express **3**, 3240–3249 (2012). [CrossRef] [PubMed]

21. B. Cox, T. Tarvainen, and S. Arridge, “Multiple illumination quantitative photoacoustic tomography using transport and diffusion models,” in “*Tomography and Inverse Transport Theory*,” G. Bal, D. Finch, J. Schotland, P. Kuchment, and P. Stefanov, eds. (American Mathematical Society, Providence, RI, USA, 2012), pp. 1–12.

## Acknowledgments

## References and links

1. | M. Xu and L. V. Wang, “Analytic explanation of spatial resolution related to bandwidth and detector aperture size in thermoacoustic or photoacoustic reconstruction,” Phys. Rev. E |

2. | L. Wang, “Tutorial on photoacoustic microscopy and computed tomography,” IEEE J. Sel. Top. Quant. |

3. | Z. Guo, C. Li, L. Song, and L. V. Wang, “Compressed sensing in photoacoustic tomography in vivo,” J. Biomed. Opt. |

4. | B. Cox, J. G. Laufer, S. R. Arridge, and P. C. Beard, “Quantitative spectroscopic photoacoustic imaging: a review,” J. Biomed. Opt. |

5. | G. Bal and K. Ren, “Multi-source quantitative photoacoustic tomography in a diffusive regime,” Inverse Probl. |

6. | J. Ripoll and V. Ntziachristos, “Quantitative point source photoacoustic inversion formulas for scattering and absorbing media,” Phys. Rev. E |

7. | Z. Yuan and H. Jiang, “Quantitative photoacoustic tomography: Recovery of optical absorption coefficient maps of heterogeneous media,” Appl. Phys. Lett. |

8. | B. Banerjee, S. Bagchi, R. M. Vasu, and D. Roy, “Quantitative photoacoustic tomography from boundary pressure measurements: noniterative recovery of optical absorption coefficient from the reconstructed absorbed energy map,” J. Opt. Soc. Am. A |

9. | B. T. Cox, S. R. Arridge, K. P. Köstli, and P. C. Beard, “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt. |

10. | L. Yin, Q. Wang, Q. Zhang, and H. Jiang, “Tomographic imaging of absolute optical absorption coefficient in turbid media using combined photoacoustic and diffusing light measurements,” Opt. Lett. |

11. | T. Jetzfellner, D. Razansky, A. Rosenthal, R. Schulz, K. H. Englmeier, and V. Ntziachristos, “Performance of iterative optoacoustic tomography with experimental data,” Appl. Phys. Lett. |

12. | B. T. Cox, S. R. Arridge, and P. C. Beard, “Estimating chromophore distributions from multiwavelength photoacoustic images,” J. Opt. Soc. Am. A |

13. | G. Bal and K. Ren, “On multi-spectral quantitative photoacoustic tomography in diffusive regime,” Inverse Probl. |

14. | G. Bal and G. Uhlmann, “Inverse diffusion theory of photoacoustics,” Inverse Probl. |

15. | R. J. Zemp, “Quantitative photoacoustic tomography with multiple optical sources,” Appl. Opt. |

16. | P. Shao, B. Cox, and R. J. Zemp, “Estimating optical absorption, scattering, and grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt. |

17. | P. Shao, T. Harrison, and R. J. Zemp, “Iterative algorithm for multiple illumination photoacoustic tomography (mipat) using ultrasound channel data,” Biomed. Opt. Express |

18. | H. Gao, S. Osher, and H. Zhao, “Quantitative photoacoustic tomography,” in “ |

19. | K. Ren, H. Gao, and H. Zhao, “A Hybrid Reconstruction Method for Quantitative PAT,” SIAM J. Imaging Sci. |

20. | S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. |

21. | B. Cox, T. Tarvainen, and S. Arridge, “Multiple illumination quantitative photoacoustic tomography using transport and diffusion models,” in “ |

**OCIS Codes**

(100.0100) Image processing : Image processing

(100.3010) Image processing : Image reconstruction techniques

(110.5120) Imaging systems : Photoacoustic imaging

(110.6960) Imaging systems : Tomography

**ToC Category:**

Image Reconstruction and Inverse Problems

**History**

Original Manuscript: July 3, 2013

Revised Manuscript: September 5, 2013

Manuscript Accepted: September 11, 2013

Published: September 24, 2013

**Citation**

Tyler Harrison, Peng Shao, and Roger J. Zemp, "A least-squares fixed-point iterative algorithm for multiple illumination photoacoustic tomography," Biomed. Opt. Express **4**, 2224-2230 (2013)

http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-4-10-2224

Sort: Year | Journal | Reset

### References

- M. Xu and L. V. Wang, “Analytic explanation of spatial resolution related to bandwidth and detector aperture size in thermoacoustic or photoacoustic reconstruction,” Phys. Rev. E67, 056605 (2003). [CrossRef]
- L. Wang, “Tutorial on photoacoustic microscopy and computed tomography,” IEEE J. Sel. Top. Quant.14, 171–179 (2008). [CrossRef]
- Z. Guo, C. Li, L. Song, and L. V. Wang, “Compressed sensing in photoacoustic tomography in vivo,” J. Biomed. Opt.15, 021311 (2010). [CrossRef] [PubMed]
- B. Cox, J. G. Laufer, S. R. Arridge, and P. C. Beard, “Quantitative spectroscopic photoacoustic imaging: a review,” J. Biomed. Opt.17, 061202 (2012). [CrossRef] [PubMed]
- G. Bal and K. Ren, “Multi-source quantitative photoacoustic tomography in a diffusive regime,” Inverse Probl.27, 075003 (2011). [CrossRef]
- J. Ripoll and V. Ntziachristos, “Quantitative point source photoacoustic inversion formulas for scattering and absorbing media,” Phys. Rev. E71, 031912 (2005). [CrossRef]
- Z. Yuan and H. Jiang, “Quantitative photoacoustic tomography: Recovery of optical absorption coefficient maps of heterogeneous media,” Appl. Phys. Lett.88, 231101 (2006). [CrossRef]
- B. Banerjee, S. Bagchi, R. M. Vasu, and D. Roy, “Quantitative photoacoustic tomography from boundary pressure measurements: noniterative recovery of optical absorption coefficient from the reconstructed absorbed energy map,” J. Opt. Soc. Am. A25, 2347–2356 (2008). [CrossRef]
- B. T. Cox, S. R. Arridge, K. P. Köstli, and P. C. Beard, “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt.45, 1866–1875 (2006). [CrossRef] [PubMed]
- L. Yin, Q. Wang, Q. Zhang, and H. Jiang, “Tomographic imaging of absolute optical absorption coefficient in turbid media using combined photoacoustic and diffusing light measurements,” Opt. Lett.32, 2556–2558 (2007). [CrossRef] [PubMed]
- T. Jetzfellner, D. Razansky, A. Rosenthal, R. Schulz, K. H. Englmeier, and V. Ntziachristos, “Performance of iterative optoacoustic tomography with experimental data,” Appl. Phys. Lett.95, 013703 (2009). [CrossRef]
- B. T. Cox, S. R. Arridge, and P. C. Beard, “Estimating chromophore distributions from multiwavelength photoacoustic images,” J. Opt. Soc. Am. A26, 443–455 (2009). [CrossRef]
- G. Bal and K. Ren, “On multi-spectral quantitative photoacoustic tomography in diffusive regime,” Inverse Probl.28, 025010 (2012). [CrossRef]
- G. Bal and G. Uhlmann, “Inverse diffusion theory of photoacoustics,” Inverse Probl.26, 085010 (2010). [CrossRef]
- R. J. Zemp, “Quantitative photoacoustic tomography with multiple optical sources,” Appl. Opt.49, 3566–3572 (2010). [CrossRef] [PubMed]
- P. Shao, B. Cox, and R. J. Zemp, “Estimating optical absorption, scattering, and grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt.50, 3145–3154 (2011). [CrossRef] [PubMed]
- P. Shao, T. Harrison, and R. J. Zemp, “Iterative algorithm for multiple illumination photoacoustic tomography (mipat) using ultrasound channel data,” Biomed. Opt. Express3, 3240–3249 (2012). [CrossRef] [PubMed]
- H. Gao, S. Osher, and H. Zhao, “Quantitative photoacoustic tomography,” in “Mathematical Modeling in Biomedical Imaging II: Optical, Ultrasound, and Opto-Acoustic Tomographiess,” vol. 2035 of Lecture Notes in Mathematics: Mathematical Biosciences Subseries, H. Ammari, ed. (Springer-Verlag, Berlin, 2011), pp. 131–158.
- K. Ren, H. Gao, and H. Zhao, “A Hybrid Reconstruction Method for Quantitative PAT,” SIAM J. Imaging Sci.6, 32–55 (2013). [CrossRef]
- S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl.15, R41 (1999). [CrossRef]
- B. Cox, T. Tarvainen, and S. Arridge, “Multiple illumination quantitative photoacoustic tomography using transport and diffusion models,” in “Tomography and Inverse Transport Theory,” G. Bal, D. Finch, J. Schotland, P. Kuchment, and P. Stefanov, eds. (American Mathematical Society, Providence, RI, USA, 2012), pp. 1–12.

## Cited By |
Alert me when this paper is cited |

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

« Previous Article | Next Article »

OSA is a member of CrossRef.