Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Wavelength and pulse energy optimization for detecting hypoxia in photoacoustic imaging of the neonatal brain: a simulation study

Open Access Open Access

Abstract

Cerebral hypoxia is a severe injury caused by oxygen deprivation to the brain. Hypoxia in the neonatal period increases the risk for the development of neurological disorders, including hypoxic-ischemic encephalopathy, cerebral palsy, periventricular leukomalacia, and hydrocephalus. It is crucial to recognize hypoxia as soon as possible because early intervention improves outcomes. Photoacoustic imaging, using at least two wavelengths, through a spectroscopic analysis, can measure brain oxygen saturation. Due to the spectral coloring effect arising from the dependency of optical properties of biological tissues to the wavelength of light, choosing the right wavelength-pair for efficient and most accurate oxygen saturation measurement and consequently quantifying hypoxia at a specific depth is critical. Using a realistic neonate head model and Monte Carlo simulations, we found practical wavelength-pairs that quantified regions with hypoxia most accurately at different depths down to 22 mm into the cortex neighboring the lateral ventricle. We also demonstrated, for the first time, that the accuracy of the sO2 measurement can be increased by adjusting the level of light energy for each wavelength-pair. Considering the growing interest in photoacoustic imaging of the brain, this work will assist in a more accurate use of photoacoustic spectroscopy and help in the clinical translation of this promising imaging modality. Please note that explaining the effect of acoustic aberration of the skull is not in the scope of this study.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Tissue oxygen saturation (sO2) level, defined as the ratio of the concentration of oxy-hemoglobin (HbO) to the concentration of the combination of oxy- and deoxy-hemoglobin (HbR) (i.e., total hemoglobin [HbT]), is a critical hemodynamic parameter in clinical practice [111]. Insufficient oxygen supply can lead to transient physiological dysfunction or permanent tissue damage [46]. One such case is cerebral hypoxia, in which the sO2 level of the brain drops below a normal limit, causing damage to the brain and resulting in hypoxic-ischemic encephalopathy associated with variable degrees of sensorimotor and cognitive deficits [4,5]. The severity of brain damage, at least in part, depends on the level of sO2 drop [1214]. It is essential to monitor the sO2 level of the brain parenchyma to detect hypoxia promptly because early intervention is suggested to improve functional outcomes in children with a diagnosis of hypoxic-ischemic encephalopathy [4,5].

There are several diagnostic tools to measure and monitor sO2 in the neonatal brain. Cranial ultrasonography (CUS) has the benefit of being available at the newborn’s bedside; however it has a low sensitivity for detecting the anomalies that are caused by hypoxia due to its inability to measure tissue oxygen saturation. In contrast, second-stage diagnostic tools, magnetic resonance imaging (MRI), computed tomography (CT), and positron emission tomography (PET), have high sensitivity and specificity for the detection of brain injuries. However, they may not be used as frontline modalities for routine screening of at-risk newborns, because MRI takes a long time to scan (∼1 hour), CT uses ionizing radiation, and PET requires a positron-emitting radionuclide. A given neonate has to be transferred to a scanning room to undergo such an imaging test. Near-infrared spectroscopy (NIRS) can potentially assist in measuring sO2 [1518], but it has a poor spatial resolution which limits its use to the mapping at the cortical level [16]. Photoacoustic (PA) imaging works based on the acoustic detection of optical absorption from endogenous tissue chromophores [1923]. Oxy- and deoxy-hemoglobin are the primary chromophores in the tissue. The very high sensitivity of PA to HbO and HbR has been well demonstrated in both animal models [2430] and humans [3134]. Using two wavelengths, the concentrations of HbO and HbR, and consequently, sO2 can be calculated, the results of which can be used for diagnosis, therapy planning, and treatment monitoring [3541].

sO2 can be measured accurately using two wavelengths only if the fluence at a given location is known for both wavelengths. Since optical properties of the tissue (e.g., absorption and scattering coefficients) are wavelength dependent, the fluence for different wavelengths is different; this effect known as spectral coloring [4245], introduces an error in the calculation of sO2. Such error increases in deeper depths due to the nonlinear nature of the light attenuation. Therefore, choosing the right wavelength-pair for efficient and most accurate oxygenation saturation measurement, and consequently accurately quantifying hypoxia at a specific depth is critical. There are several methods proposed to mitigate the spectral coloring effect [4554] including, model-based fluence compensation algorithms [4547], deep learning [49,50], the unmixing method based on fluence eigenspectra [51] and wavelength optimization method [45,5254]. Recently, Holuchi et al. [45] used Monte Carlo (MC) simulations to find optimum wavelength-pairs for accurate sO2 measurements of skin vessels. Since the skin vessels are located in shallow depths, their model does not consider the effect of depth, hence it may not be adequate for quantifying cerebral hypoxic regions in the brain. Further, to our best knowledge, nobody has explored the effect of optimization of the optical energy of wavelength-pairs in the accuracy of oxygen saturation measurements.

Here, we use a neonate head model made from MRI images, and utilize Monte Carlo simulations to find the optimum wavelength combinations and their corresponding optimum energies at different depths, for accurate quantification of the degree of hypoxia in the brain when using photoacoustic spectroscopic imaging. This paper is organized as follows: in Section 2, we explain the principal of sO2 measurement and spectral coloring; in Section 3, we introduce our neonatal head model including geometry, optical properties and specified hypoxic regions; in Section 4, we present the results in the following order: 1) we explain our proposed method using an example wavelength-pair (780 nm & 830 nm) for an example hypoxic region at the depth of 14mm inside the brain tissue, 2) considering equal energies for the example wavelengths, we reported the error in the estimation of hypoxia, 3) we introduced energy optimization using the example wavelength-pair and reported the error in the estimation of hypoxia, 4) we repeated the procedure for all wavelength-pairs between 500-1000 nm in 10 nm steps at 6 different depths, 5) we reported the hypoxia estimation error for all wavelength-pairs with both equal energies and optimized energies at all depths. And finally, in section 5, we discussed several considerations and assumptions that we used in this work in more details.

2. Theory: the principle of spectral coloring

When a tissue is irradiated by a nanosecond pulsed laser, thermoelastic expansion of tissue chromophores occurs, resulting in an initial pressure and emission of acoustic waves. The initial pressure, p0, at the imaging target is proportional to the absorbance, A [J/m3], i.e. the amount of absorbed energy per volume (see Eq. (1)). A is the product of absorption coefficient, μa[1/m], and fluence, F [J/m2].

$${p_0} = \eta \varGamma A = \eta \varGamma {\mathrm{\mu} _a}F$$
where, Γ is Gruneisen parameter, and η is the heat conversion efficiency of the tissue. The acoustic waves generated from the initial pressure sources, travel through the tissue, reaching the ultrasound transducers. Finally, the PA signals acquired from transducers are used to reconstruct a PA image to form a map of relative tissue absorption coefficients. Knowing absorption coefficients at two wavelengths, ${\lambda _1}\; $ and ${\lambda _2}$, one can calculate the concentration of sO2 as the ratio of the concentration of HbO to the concentration of HbT (see Eq. (2)):
$$sO_2^{exact} = \frac{{\; \; [{HbO} ]}}{{\; \; [{HbT} ]}} = \frac{{\; \; [{HbO} ]}}{{\; \; [{HbO} ]+ [{HbR} ]}} = \frac{{{{\varepsilon }_{HbR}}({{\lambda_2}\; } ){{\mu} _a}({{\lambda_1}\; } )- {\mathrm{\varepsilon }_{HbR}}({{\lambda_1}\; } ){{\mu} _a}({{\lambda_2}\; } )}}{{\varDelta {\mathrm{\varepsilon }_{Hb}}\; ({{\lambda_2}} ){{\mu} _a}({{\lambda_1}\; } )- \varDelta {\mathrm{\varepsilon }_{Hb}}\; ({{\lambda_1}} ){{\mu} _a}({{\lambda_2}\; } )}}$$
where ${\mathrm{\varepsilon }_{HbR}}$ and ${\mathrm{\varepsilon }_{HbO}}$ are the extinction coefficients of HbR and HbO, respectively, and $\varDelta {\mathrm{\varepsilon }_{Hb}}$=${\mathrm{\varepsilon }_{HbR}} - {\mathrm{\varepsilon }_{HbO}}$. Having two PA images at two different wavelengths and knowing that ${{\mu} _a}(\lambda )= {p_0}(\lambda )/F(\lambda )$, Eq. (2) can be written as,
$$sO_2^{exact} = \frac{{\; \; [{HbO} ]}}{{\; \; [{HbO} ]+ [{HbR} ]}} = \frac{{{\varepsilon _{HbR}}({{\lambda_2}} ){p_0}({{\lambda_1}} )/F({{\lambda_1}} )- {\varepsilon _{HbR}}({{\lambda_1}} ){p_0}({{\lambda_2}} )/F({{\lambda_2}} )}}{{\varDelta {\varepsilon _{Hb}}({{\lambda_2}} ){p_0}({{\lambda_1}} )/F({{\lambda_1}} )- \varDelta {\varepsilon _{Hb}}({{\lambda_1}} ){p_0}({{\lambda_2}} )/F({{\lambda_2}} )}}$$

Equation (3), suggests that if we know the fluence at a given location for the two wavelengths, ${\lambda _1}$ and ${\lambda _2}$, we can measure the exact value of sO2. In practice, the fluence information is not available, hence p0 is used instead of μa in Eq. (2) (see Eq. (4)).

$$sO_2^{estim} = \frac{{{\varepsilon _{Hb}}({{\lambda_2}} ){p_0}({{\lambda_1}} )- {\varepsilon _{Hb}}({{\lambda_1}} ){p_0}({{\lambda_2}} )}}{{\varDelta {\varepsilon _{Hb}}({{\lambda_2}} ){p_0}({{\lambda_1}} )- \varDelta {\varepsilon _{Hb}}({{\lambda_1}} ){p_0}({{\lambda_2}} )}}$$

It can be seen from Eq. (3) that if F(λ1) = F(λ2) then $sO_2^{estim}$= $sO_2^{exact}$. However, F(λ1) and F(λ2) are wavelength-dependent and generally are not equal. This effect, known as spectral coloring, is the source of error in the estimation of sO2. Although an optimum wavelength pair decreases the fluence difference and improves the accuracy of sO2 measurement, adjusting the energy of each selected wavelength could also decrease the fluence difference, because fluence is directly proportional to the pulse energy. Therefore, we simultaneously optimize the wavelengths and their pulse energy for the most accurate sO2 measurement. A specific wavelength/energy selection that decreases the spectral coloring effect at a given depth, may worsen it at other depths. Therefore, we find and report best combinations of wavelengths at different depths.

3. Materials and methods

3.1 Simulation

We used an MRI atlas of a 37 weeks neonate, made by Cooper et al. [55], as the head model. The size of the model is 136 mm ×185 mm ×146 mm with a resolution of 0.86 mm. The model is composed of an extra-cranial tissue (ECT), cerebrospinal fluid (CSF), gray matter, white matter, cerebellum, and brainstem. ECT includes scalp and skull; the averaged optical property of scalp and skull has been assigned to this layer. Following the work by Steven L. Jacques [56], we combined gray matter, white matter, cerebellum, and brain stem, and named it brain tissue and assigned the “neonatal brain” optical properties used in [56] to it (see Fig. 1(a)). As a result of these approximations, in terms of optical properties, our head model composes of three tissues: ECT, CSF and brain.

 figure: Fig. 1.

Fig. 1. (a) Head model: sagittal view of a 37 weeks neonate atlas used for simulations. ECT: extra-cranial tissue, CSF: cerebrospinal fluid. The areas with different colors show regions of hypoxia. (b) At a given experiment, the selected voxels in normal tissue and a region of hypoxia are used for analysis.

Download Full Size | PDF

To model hypoxia, a pyramid-shaped region is created inside the head just above the ventricle and continues until it reaches the surface of the cortex. The pyramid shape of the region represents the branching of the vessels towards the surface of the cortex [5760]. To assess the degree of hypoxia at different depths, this region is divided into 6 sub-regions, each having a thickness of 4 mm and a center-to-center distance of 4 mm, located at depths 2 mm to 22 mm from the surface of the cortex. We call each sub-region, a region of hypoxia (ROH); they are all shown in Fig. 1(a) with different colors.

As for light illumination, a collimated source with a diameter of 2 cm (area ∼3 cm2) is placed on the surface of the scalp, illuminating the brain tissue (see Fig. 1(b)). Light illumination is performed at 51 wavelengths ranging from 500 nm to 1000 nm in 10 nm steps. For each wavelength, 108 photons were simulated inside the head model using MCX software [61] and resulting fluence distributions were used for initial pressure (Eq. (1)) and $sO_2^{estim}$ (Eq. (4)) calculations. Optical properties used for these simulations are explained in section 3.3.

3.2 Simulation

We defined Hypoxiatrue as

$$Hypoxi{a^{true}} = \left( {1 - \frac{{sO_{2,ROH}^{true}}}{{sO_{2,normal}^{true}}}} \right) \times 100$$
where true superscript denotes the ground truth values used in the simulations. In Eq. (5), $sO_{2,normal}^{true}$ is the oxygen saturation for normal brain tissue and $sO_{2,ROH}^{true}$ is that for regions of hypoxia. We considered four $sO_{2,normal}^{true}$ values including 0.6, 0.7, 0.8 or 0.9 [62,63]. For each $sO_{2,normal}^{true}$, we considered five values of $sO_{2,ROH}^{true}$ such that the Hypoxiatrue in Eq. (5) equals to 10%, 20%, 30% 40% or 50%. In total we simulated 20 different scenarios at each depth, and for wavelengths between 500 nm to 1000 nm.

3.3 Optical properties

We calculated the values of scattering and absorption coefficients of ECT, CSF and brain tissue from [56] for our simulations.

Scattering coefficient: Eq. (6) is used to estimate the reduced scattering coefficient, ${\mu} _s^{\prime}$, of biological tissue for wavelengths 500 nm to 1000 nm.

$${\mu} _s{^{\prime}}(\lambda )= a^{\prime}\left( {{f_{Ray}}{{\left( {\frac{\lambda }{{500({nm} )}}} \right)}^{ - 4}} + ({1 - {f_{Ray}}} ){{\left( {\frac{\lambda }{{500({nm} )}}} \right)}^{ - {b_{Mie}}}}} \right)$$
$a^{\prime}$ is the scaling factor equals to ${\mu} _s{^{\prime}}({\lambda = 500nm} )$, ${f_{Ray}}$ and $1 - {f_{Ray}}$ are contributions of Rayleigh and Mie scattering, respectively. $a^{\prime}$, ${f_{Ray}}$ and ${b_{Mie}}$ are the specific constants to biological tissue. The scattering coefficient, ${{\mu} _s}$, is then obtained from the relationship between ${\mu} _s^{\prime}$ and ${\mu} _s^{}$ (i.e., ${\mu} _s^{\prime} = {{\mu} _s}({1 - g} )$). Here, g is the anisotropy factor. The constants in Eq. (6) used for ECT and brain tissues are presented in Table 1. For ECT, we used averaged values of ${f_{Ray}}$, $a^{\prime}$ and ${b_{Mie}}$ of the skin and bone. The reduced scattering coefficient of CSF was set to 0.3 mm-1 [64].

Tables Icon

Table 1. Constants used to generate ECT and brain tissue scattering coefficient. ECT: extracerebral tissue

Absorption coefficient: Major absorbers in biological tissues are HbR and HbO, water, melanin, and fat. Using Eq. (7), we calculated the tissue absorption coefficient at different wavelengths.

$$\begin{aligned}{{\mu} _a} &= B({sO_2^{true}} ){{\mu} _{a,HbO}} + B({1 - ({sO_2^{true}} )} ){{\mu} _{a,HbR}} + W{{\mu} _{a,water}} + F{{\mu} _{a,fat}}\\&+ M{{\mu} _{a,melanosome}}\end{aligned}$$
where, B indicates the averaged blood volume fraction, W, water content, F, fat content, and M, melanosome volume fraction of a specific tissue. Table 2 presents the values of $sO_2^{true}$, B, W, F, and M that have been found in the literature [6570]. Here, we assumed that the brain tissue consists of water, fat and blood, and that the absorption coefficient of the CSF is same as that of water. Refractive index is set to 1.37 for all tissues.

Tables Icon

Table 2. Constants used in Eq. (7) to generate the absorption coefficient of ECT, CSF, normal brain tissue, and ROH. ECT: extracranial tissue, and CSF: cerebrospinal fluid

4. Results

4.1 Synchronization optimization for efficient averaging

Although all the simulations have been performed in a three-dimensional model, we only presented the coronal cross-section results, where the center of the light illumination source was placed and the fluence was maximum. 120 scenarios were tested: four sO2 levels for normal tissue and five hypoxia levels at six different depths. We also tested these scenarios at 51 wavelengths, totaling 6120 simulations. To better explain the optimization procedure in this study, initially we explained sO2 calculation and hypoxia measurement algorithms for a sample wavelength-pair, i.e., 780 nm and 830 nm; simulations were run for four normal tissue sO2 levels and 5 hypoxia at the depth of 14 mm. We then discussed the results for all 6120 simulations.

Figures 2(a) and 2(b) show initial pressure maps at 780nm and 830nm, respectively. Figure 2(c) shows the estimate sO2 calculated for wavelength-pair 780nm and 830nm using Eq. (4). To quantify hypoxia in Fig. 2(c), we calculated the pixel intensity histogram of hypoxic and surrounding normal tissue region depicted in Fig. 2(c) (see Fig. 3(a)). As seen in Fig. 3(a), there is a distinct separation between two regions in the histogram corresponding to normal and hypoxic tissue; we named the mean value of the normal region, $sO_{2,normal}^{estim}$, and that for the hypoxic region, $sO_{2,ROH}^{estim}$ (indicated by two dashed-red lines in Fig. 3(a)). Similar to Hypoxiatrue in Eq. (5), we defined Hypoxiaestim to quantify the severity of the hypoxic region (see Eq. (8)).

$$Hypoxi{a^{estim}}({{\lambda_1},{\lambda_2}} )= \left( {1 - \frac{{sO_{2,ROH}^{estim}({{\lambda_1},{\lambda_2}} )}}{{sO_{2,normal}^{estim}({{\lambda_1},{\lambda_2}} )}}} \right) \times 100$$

 figure: Fig. 2.

Fig. 2. (a) Map of initial pressure at 780nm, (b) map of initial pressure at 830nm, (c) map of estimated sO2 using the wavelength pair, 780nm and 830nm with $sO_{2,normal}^{true}$=0.7 and $Hypoxi{a^{true}}$=30%. The area inside the shaded region in c shows the hypoxic region and the area between two black borders contains the pixels of surrounding normal tissue which are involved in the calculation of hypoxia.

Download Full Size | PDF

 figure: Fig. 3.

Fig. 3. (a) pixel intensity histogram of estimated sO2 for hypoxic and surrounding normal tissue located at a depth of 14mm with a wavelength-pair 780nm and 830nm. Dashed-red lines indicate the mean value of the normal region, $sO_{2,normal}^{estim}$, and that for the hypoxic region, $sO_{2,ROH}^{estim}$. (b) $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for 20 scenarios (four sO2 levels for normal tissue and five hypoxia) at a depth of 14mm for the wavelength-pair, 780nm, and 830nm. The $Hypoxi{a^{estim}}$ calculated for the case presented in Fig. 3(a) is specified with a red dotted circle in Fig. 3(b).

Download Full Size | PDF

Figure 3(b) shows $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for 20 scenarios at a depth of 14 mm for the wavelength-pair, 780nm and 830nm. In this figure, we studied the accuracy of measuring different levels of $Hypoxi{a^{estim}}$ when the surrounding normal tissue has the $sO_{2,normal}^{true}$ value of 0.6, 0.7, 0.8, or 0.9; the black dashed line shows the ideal $Hypoxi{a^{estim}}$ which is equal to $Hypoxi{a^{true}}$ at any given $sO_{2,normal}^{true}$; the $Hypoxi{a^{estim}}$ calculated for the case presented in Fig. 3(a) is specified with a red dotted circle in Fig. 3(b). It can be seen that there is an underestimation for lower values of $sO_{2,normal}^{true}$ (i.e., 0.6, 0.7 and 0.8), and there is an overestimation for 0.9.

We reported the error for estimation of $Hypoxi{a^{estim}}$ at each wavelength-pair, averaged over 20 scenarios, using Eq. (9). The error for hypoxia measurement for the wavelength-pair 780nm and 830nm at a depth of 14 mm was calculated as 19%.

$$Error({{\lambda_1},{\lambda_2}} )= \frac{1}{{20}}\mathop \sum \nolimits_{all\; sO_{2,normal}^{true}} \mathop \sum \nolimits_{all\; hypoxi{a^{true}}} \left|{\frac{{Hypoxi{a^{estim}}({{\lambda_1},{\lambda_2}} )- Hypoxi{a^{true}}}}{{Hypoxi{a^{true}}}}} \right|\times 100$$

So far, we assumed that the pulse energy illuminated to the defined area, ∼3cm2 (see Section 3.1), is 20 mJ/cm2, which is below the maximum permissible exposure (MPE) dictated by the American National Standards Institute (ANSI) [71]. MPE for wavelengths 500 nm to 700 nm is constant (i.e., 20 mJ/cm2), and increases at wavelengths greater than 700 nm (see the MPE graph versus wavelength in Fig. 4). Additionally, as we explained in Section 2, adjusting the energy of each selected wavelength can also improve the accuracy of hypoxia quantification. Therefore, in our example, we searched for pulse energies for 780 nm and 830 nm that most accurately calculate $Hypoxi{a^{estim}}$. We increased the level of energy in 10% increments from the initial baseline value (see red astros in Fig. 4) until it reached the MPE limit at the given wavelength. MPE was calculated using Eq. (10). For instance, MPE is calculated as 28.9 mJ/cm2 for 780 nm and 36.4 mJ/cm2 for 830 nm.

$$MPE = 2{C_A} \times {10^{ - 2}}({J/c{m^2}} ),\; \; {C_A} = \left\{ \begin{array}{ll} {1.0} &0.400 < \lambda < 0.700\; ({{\mu} m} )\\ {{{10}^{2({\lambda - 0.700} )}}\; \;} &0.700 < \lambda < 1.050\; ({{\mu} m} ) \end{array} \right.$$

 figure: Fig. 4.

Fig. 4. Maximum permissible exposure graph according to ANSI limit (black line), taken from [71]; black arrows show MPE for wavelengths 780nm and 830nm; the red astros show discrete levels of energy that are used in energy optimization.

Download Full Size | PDF

We tested 45 different combinations of energy levels; five energy levels for 780 nm, and nine energy levels for 830 nm, and calculated the Error(780nm, 830nm) (see Eq. (9)) for all energy level combinations. The energy level combination that resulted in the smallest Error was reported as optimum.

Figure 5 (a) shows $Hypoxi{a^{estim}}$=26%, calculated at optimum energy levels (i.e., 28mJ/cm2 at 780nm and 26mJ/cm2 at 830nm), which yielded 3% improvement, compared to when energy levels were equal; for the explanation of how $sO_{2,normal}^{estim}$ and $sO_{2,ROH}^{estim}$ were calculated, see the beginning of this section. Figure 5(b) shows $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for all 20 scenarios at a depth of 14 mm for the wavelength-pair, 780 nm and 830 nm at their optimized level of energy, respectively.

 figure: Fig. 5.

Fig. 5. (a) pixel intensity histogram of estimated sO2 after energy optimization calculated for hypoxic and surrounding normal tissue located at the depth of 14mm with a wavelength-pair, 780nm and 830nm. Dashed red lines indicate the mean value of the normal region, $sO_{2,normal}^{estim}$, and that for the hypoxic region, $sO_{2,ROH}^{estim}$. (b) $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for 20 scenarios (four sO2 levels for normal tissue and five levels for hypoxic regions) at the depth of 14mm for the wavelength-pair, 780nm and 830nm after energy optimization. The $Hypoxi{a^{estim}}$ calculated for the case presented in Fig. 5(a) is specified with a red dotted circle in Fig. 5(b).

Download Full Size | PDF

In Fig. 5 (b), it can be seen that after energy optimization, the error for the estimation of $Hypoxi{a^{estim}}$ when $sO_{2,normal}^{true}$=0.6, 0.7, and 0.8 were decreased, whereas the error was increased when $sO_{2,normal}^{true}$=0.9, resulting in an overall smaller error for hypoxia measurement, i.e., 18%.

4.2 Wavelength optimization

We calculated the error (using Eq. (9)) for all wavelength-pairs in the wavelength range 500 nm to 1000 nm, and depths from 2 mm to 22 mm. For each wavelength-pair we considered all eligible energy-pairs and select the one with the smallest error. The results for depths 2 mm, 6 mm, 10 mm, 14 mm, 18 mm, and 22 mm are plotted in top triangle maps in Fig. 6, where each pixel value indicates the corresponding error for the wavelength-pair λ1 and λ2; those results after energy optimization are plotted in bottom triangle maps. It can be seen that the majority of the wavelength-pairs that produce small error are between 750 nm to 900 nm. At greater depths, the wavelengths near 1000 nm also produce small error.

 figure: Fig. 6.

Fig. 6. Hypoxia estimation error calculated at wavelength-pairs in the wavelength range 500 nm to1000 nm (with steps of 10 nm) without energy optimization (top triangle) and with energy optimization (bottom triangle) for depths (a)2 mm, (b)6 mm, (c)10 mm, (d)14 mm, (e)18 mm and (f)22 mm.

Download Full Size | PDF

Figure 7 shows energy ratio, E(λ2)/E(λ1), at each wavelength-pair for depths 2 mm, 6 mm, 10 mm, 14 mm, 18 mm and 22 mm after energy optimization, where E(λ) represents pulse energy of the selected wavelength, λ.

 figure: Fig. 7.

Fig. 7. Map of energy ratio E(λ2)/E(λ1) at wavelength-pairs (λ1,λ2) in the wavelength range 500 nm to 1000 nm (with steps of 10 nm) for depths (a) 2 mm, (b) 6 mm, (c) 10 mm, (d) 14 mm, (e) 18 mm and (f) 22 mm. E(λ) represents the pulse energy of the selected wavelength, λ.

Download Full Size | PDF

4.3 Transducer noise equivalent pressure and sensitivity

For the practical utility of the proposed method, we should consider the characteristics of real-world transducers. The presented results so far are with the assumptions that (i) any pressure level, however small, can be detected by a US transducer, and that (ii) the PA image reconstruction is perfect. Studying the impact of the choice of an image reconstruction method [72,73] and its performance in hypoxia quantification is out of the scope of this study. Therefore, with the assumption of using a perfect image reconstruction, i.e., the PA image is equal to the initial pressure map, we investigated the detectability of the initial pressure sources.

Wavelength-pairs that fulfill the following two conditions for at least half of the photoacoustic image pixels in Figs.2 (a) and (b) are considered detectable if: (i) the pressure is greater than the noise equivalent pressure (NEP) of the transducer at each wavelength, and (ii) these detected pressures have a difference that is within the sensitivity of the transducer. We, therefore, used the median value of p0 in our computations in this section. Figure 8 shows the boxplot of initial pressure of all pixels in both normal and hypoxic regions for a special case ($sO_{2,normal}^{true}$=0.7 and $Hypoxi{a^{true}}$=30%) at depth of 14 mm for different wavelengths, where both NEP and transducer sensitivity were considered to be 5 Pa. Considering NEP = 5 Pa and sensitivity = 5 Pa, we reproduced the maps in the bottom triangles in Fig. 6, the results are shown in Fig. 9.

 figure: Fig. 8.

Fig. 8. Boxplot of initial pressure in both normal and hypoxic regions for a special case ($sO_{2,normal}^{true}$=0.7 and $Hypoxi{a^{true}}$=30%) at depth 14mm for different wavelengths. Black dots indicate median values of initial pressures. Red dashed line shows the NEP level.

Download Full Size | PDF

 figure: Fig. 9.

Fig. 9. Consideration of transducer NEP and sensitivity in hypoxia estimation errors calculated at wavelength-pairs in the wavelength range 500 nm to1000 nm (with steps of 10 nm) with energy optimization for depths (a) 2 mm, (b)6 mm, (c)10 mm, (d)14 mm, (e)18 mm and (f) 22 mm. Grey areas in all maps are the areas where at least one of the detectability conditions has not been fulfilled.

Download Full Size | PDF

To find the best wavelength-pairs we sorted the results of the 1275 wavelength-pairs based on their error percentage in an ascending order. Table 3 shows ten best wavelength-pairs for the depth of 14 mm. Similar tables have been generated for other depths, see Supplement 1, Tables S1 to S6; these tables can be used as a guideline for wavelength-pair selection.

Tables Icon

Table 3. Optimized wavelength-pairs at depth=14 mm with corresponding Error (λ1, λ2) sorted in ascending order for equal and optimized energies without and with considering NEP

4.4 Linearity

For a given value of normal tissue oxygen saturation, $sO_{2,normal}^{true}$, we performed a linear regression on the estimated hypoxia and calculated the linearity by the slope, β and coefficient of determination, R2. β<1 and β>1 indicates underestimation and overestimation in the estimated hypoxia, respectively. R2 closer to 1 means a better fit and a more accurate measurement of relative changes in hypoxia. The linearity metrics provide a measure of how accurately we detect changes in the hypoxia for a given $sO_{2,normal}^{true}$; this means even if the hypoxia result was not accurate, it can indicate the changes in hypoxia. Such analysis can be used in longitudinal studies and treatment monitoring events.

Figure 10 shows the results of the linear regression at depth 14 mm. Here we used the best wavelength-pair obtained after wavelength and energy optimization, i.e. 800 nm and 830 nm (See Table 3). Table 4 presents the results for different depths from 2 mm to 22 mm with the corresponding optimized wavelength pairs (See Supplement 1, Tables S1-S6); as can be seen there is a good linearity between $Hypoxi{a^{estim}}$ and $Hypoxi{a^{true}}$ at all these wavelength-pairs.

 figure: Fig. 10.

Fig. 10. Linear regressions calculated for estimated hypoxia at the depth 14 mm with optimized wavelength and energy, for the wavelength-pair 800 nm and 830 nm. Different markers correspond to different $sO_{2,normal}^{true}$ and red dashed lines are linear fits where β is the slope and R2 is the coefficient of determination.

Download Full Size | PDF

Tables Icon

Table 4. Slope and R2 of linear regression of estimated hypoxia for given normal tissue oxygen saturations at different depths with optimized wavelength and energy.

5. Discussion

Cerebral hypoxia in newborns can lead to severe neurological disorders if it remains undetected early [4,5]; thus, an accurate measurement of the brain tissue oxygen saturation is highly warranted. Photoacoustic imaging, using at least two wavelengths, through spectroscopic analysis, has the potential to measure brain oxygen saturation [23,24,26]. The present study showed that the spectral coloring effect, arising from the dependency of optical properties of biological tissue to the wavelength of light, can be greatly reduced by choosing an optimum wavelength-pair with their optimized energy level. Our method allowed an efficient and most accurate hypoxia quantification at different depths from cortex surface down to lateral ventricle, where usually hemorrhage and consequently hypoxia occur [74]. Without wavelength/energy optimization, the error in the quantification of hypoxia becomes significantly large which does not justify the use of photoacoustic spectroscopy for the measurement of SO2 and the detection of hypoxia. Table 5 lists the number of wavelength-pairs with hypoxia estimation errors ≤ 10%, ≤ 20% and ≤ 30% at different depths from 2 mm to 22 mm. It is evident that there are only a few suitable wavelength-pairs at deeper levels (see numbers at depths 18 mm and 22 mm in Table 5). At depth 18 mm, before energy optimization there are 2 and 24 wavelength-pairs with errors ≤ 20% and ≤ 30%, respectively; there is 4 wavelength-pairs with error ≤ 30% at depth 22 mm; and there is no wavelength-pair with hypoxia estimation errors ≤ 20% (see upper triangles in Figs.6(e), 6(f)).

Tables Icon

Table 5. Slope and R2 of linear regression of estimated hypoxia for given normal tissue oxygen saturations at different depths with optimized wavelength and energy.

It can be seen that energy optimization becomes more important at deeper levels due to the stronger effect of the spectral coloring as a result of the nonlinear increase of the attenuation of light in the tissue with depth. After energy optimization, the number of wavelength-pairs with errors ≤ 20% and ≤ 30% increases to 25 and 158, respectively, at the depth of 18 mm. At the depth 22 mm, the number of wavelength-pairs with error ≤ 30% increases to 39;also, one wavelength-pair becomes available with error ≤ 20% (see lower triangles in Figs.6(e), 6(f)).

According to the results presented in Fig. 6, the majority of optimal wavelength-pairs are within 750 nm to 900 nm for all depths. After all optimizations and taking into account practical transducer parameters i.e., NEP and sensitivity, the number of wavelength-pairs with error ≤30% or less, at depths 2 mm, 6 mm, 10 mm, 14 mm, 18 mm and 22 mm are calculated as: 105, 257, 215, 97, 45 and 2, respectively. Figure 9(f) shows the importance of our energy adjustment protocol; one can see that the eligible wavelength-pairs at a depth of 22 mm are near 1000 nm; these are wavelength-pairs that were not available before energy optimization (see upper triangles in Fig. 6(f) and Supplement 1, Table S6). Our results show that the wavelengths in the range of 500 nm to 650 nm were not selected, probably due to the large optical attenuation in this spectral window which led to initial pressures lower than NEP of the transducer (see condition (i) in section 4.3). Further, there are wavelengths in the range of 700 nm to 900 nm that were not selected, because the difference in the initial pressures generated by these two wavelengths were below the sensitivity of the transducer (see condition (ii) in section 4.3). The smaller number of wavelength pairs at depth 2 mm is probably due to the location of ROI which begins at the surface of the cortex and located below the CSF and scalp layer. Deeper ROIs are surrounded by the brain tissue.

We have taken into account several considerations/ assumptions in our proposed methodology: (i) we considered the initial pressure maps as PA images, and did not perform PA image reconstruction (i.e., we assumed that the PA reconstruction was perfect), thus we were not concerned about acoustic properties of the tissue and PA reconstruction complexities including skull aberration [22,75] mainly longitudinal wave and shear wave conversion, because they were completely independent of the spectral coloring effect, and therefore would have not changed our findings about the selected wavelength-pairs; (ii) although the NEP and sensitivity thresholds we used in section 4.3, i.e., 5 Pa, is lower than those of commercially available US transducers, by increasing the area of illumination, more light can be delivered and higher pressures can be produced (for instance by increasing the diameter of the area of illumination from 2 cm to 4 cm, the optical energy will be tripled at the depth of 14 mm). It is also worth noting that we include the noise equivalent pressure and sensitivity analyses assuming that if the signal is detectable then it can perfectly be reconstructed. We set the value of initial pressure as detection criterion meaning that if the signal is not detectable at the origin then it will not be detectable after acoustic propagation attenuation/aberration through the brain tissue and skull/scalp. However, it does not rule out the possibility that if the signal is detectable at the origin it still might not be detectable at the location of the transducers due to the acoustic distortion. This simplistic analysis was to show that our obtained suitable wavelength-pairs become unavailable when some practical acoustic considerations are taken into account, and also to show the importance of energy optimization in making new wavelength-pairs available in the presence of the practical acoustic considerations (see Fig. 9 and Fig. 6); (iii) our available digital head phantom has the resolution of 0.86 mm, which is coarser than the resolution of current photoacoustic tomography systems, therefore, in practice higher resolution SO2 maps of capillaries and arteries instead of tissue can be obtained using the proposed algorithm; (iv) the energy spectrum of the lasers we will use in practice may not necessarily reach the level of MPE provided by ANSI (see the laser profile of Opotek HE Mobile in Supplement 1, Fig. S1 [76]); (vi) dimension of the hypoxic ROI has been chosen according to the clinical information in the literature [4,5]; (vii) although more number of wavelengths can improve the accuracy of sO2 calculations and is in our plan to study in the future, we used the two-wavelength approach for simplicity., it also greatly complicates the problem and it should be investigated separately; (viii) assigning most accurate optical properties to different regions of the brain is one critical step in this study, we used the average of the reported values in the literature for the best approximation. We used an averaged optical property of scalp and skull for ECT layer (skull + scalp), considering the effect of melanin. We did not account for the variation of the concentration of melanin for different skin types. We plan to investigate the effect of skin type in our findings in a future study.

Figure 11 shows how the error in the estimation of hypoxia at different depths changes if true absorption coefficients for both wavelengths vary by +1% (blue), +2% (red), +5(cyan), -1% (green), -2% (pink) and -5 (orange) relative to the original values that we used in this study for the ECT tissue (Fig. 11(a)) as well as the brain tissue (Fig. 11(b)). Black dashed line shows the original error for the best wavelength-pair (with equal energies) at different depths (see Supplement 1, Tables S1-S6). It can be seen that changes in hypoxia estimation are non-trivial and do not follow changes in absorption coefficients which is likely due to the nonlinear nature of the problem Changes in the absorption coefficients lead to a change in fluence distributions and consequently in initial pressure distributions which in turn affects the estimation of sO2 value in the region of hypoxia and its surrounding normal tissue, and finally in the value of hypoxia. The many stages of calculation makes interpreting individual cases difficult. Such complexity though in the analysis of spectral coloring event well justifies the use of simulations. The low percentage of error at 10 mm in both graphs is probably due to the particular values of absorption coefficients of layers above the 10 mm. These results do not necessarily lead to a universal conclusion for other cases with different set of parameters. Please note that we assumed that the tissue optical properties at both wavelengths were equally modified.

 figure: Fig. 11.

Fig. 11. Changes in the error of estimated hypoxia at different depths when absorption coefficient is changed by +1%(blue), +2%(red), +5(cyan), -1%(green), -20%(pink) and -5(orange) relative to the original values of absorption coefficient considered for (a) ECT layer and (b) brain tissue. Dashed black line shows the original error. The errors are computed for the best wavelength pairs with equal energy at different depths (see Supplement 1, Tables S1-S6).

Download Full Size | PDF

In addition to optical properties, variations in the region of interest or geometry can change our current results: (i) if a region of interest with a different size, shape or at a different location is chosen, for which our methodology will be performed and new results will be obtained; (ii) if a subject with a different head geometry than our head phantom is utilized. The latter is a great topic for a comprehensive future study in which subject dependency of our methodology/results to different head geometries and age will be evaluated. In the same study, comparing the results obtained from an averaged head model from several same-age newborns to those obtained from one of the newborns’ head model will be investigated.

In clinical applications, where imaging a region is important instead of imaging a structure at a specific depth, a superposed image of several wavelength pairs will be used. We plan to investigate this in detail in a future study.

Our presented concept and methodology can certainly be extended to adult human brain where we expect being able to measure the hypoxia at much less depths. Indeed, the results will significantly be impacted by the attenuating and aberrating effects of the thick scalp and bony very thick skull (compared to those in newborn) [7781].

6. Conclusion

Cerebral hypoxia is a severe injury that increases the risk of development of neurological disorders in the neonatal period. It is crucial to recognize hypoxia as soon as possible because early intervention improves outcomes. With spectroscopic photoacoustic imaging, using at least two wavelengths, one can measure sO2 and use it for the evaluation of hypoxia. Due to the spectral coloring effect, choosing the right wavelength-pair for most accurate sO2 measurement is essential. Using a realistic neonate head model and Monte Carlo simulations, we studied the error in sO2 measurement due to the spectral coloring and identified optimum wavelength-pairs at different depths from the cortex down to the lateral ventricle (these optimum wavelength-pairs are reported in Table 3 as well as Supplement 1, Tables S1 to S6). We also demonstrated, for the first time, that the accuracy of sO2 measurements can be improved by adjusting the level of light energy for each wavelength-pair, with which a greater number of wavelength-pairs become available. Although PA spectroscopy using more wavelengths improves the accuracy of oxygenation measurement, implementation of near real-time multi-wavelength spectroscopy, using the current high energy pulsed lasers, required to study changes in hemodynamic parameters is not easily feasible. Considering the growing interest in photoacoustic brain imaging, this work could assist in a more efficient and accurate use of photoacoustic spectroscopy and help in the clinical translation of this promising imaging modality.

Funding

National Institutes of Health (R01EB027769-01, R01EB028661-01, R21 DA052657-01).

Acknowledgements

We would like to thank Robert J. Cooper from University College London for providing neonatal head model.

Disclosures

The authors have no conflicts of interest to report.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Supplemental document

See Supplement 1 for supporting content.

References

1. J. A. Howlett, F. J. Northington, M. M. Gilmore, A. Tekes, T. A. Huisman, C. Parkinson, S.-E. Chung, J. M. Jennings, J. J. Jamrogowicz, and A. C. Larson, “Cerebrovascular autoregulation and neurologic injury in neonatal hypoxic–ischemic encephalopathy,” Pediatr. Res. 74(5), 525–535 (2013). [CrossRef]  

2. P. Bickler, J. Feiner, M. Rollins, and L. Meng, “Tissue oximetry and clinical outcomes,” Anesth. Analg. 124(1), 72–82 (2017). [CrossRef]  

3. G. M. Hoffman, “Neurologic monitoring on cardiopulmonary bypass: what are we obligated to do?” The Annals of Thoracic Surgery 81(6), S2373–S2380 (2006). [CrossRef]  

4. S. Shankaran, A. R. Laptook, R. A. Ehrenkranz, J. E. Tyson, S. A. McDonald, E. F. Donovan, A. A. Fanaroff, W. K. Poole, L. L. Wright, and R. D. Higgins, “Whole-body hypothermia for neonates with hypoxic–ischemic encephalopathy,” N. Engl. J. Med. 353(15), 1574–1584 (2005). [CrossRef]  

5. G. Natarajan, A. Pappas, and S. Shankaran, “Outcomes in childhood following therapeutic hypothermia for neonatal hypoxic-ischemic encephalopathy (HIE),” In Proceedings of Seminars in perinatology; pp. 549-555.

6. C. E. Batista, H. T. Chugani, C. Juhász, M. E. Behen, and S. Shankaran, “Transient hypermetabolism of the basal ganglia following perinatal hypoxia,” Pediatric neurology 36(5), 330–333 (2007). [CrossRef]  

7. J.K.G Avanaki and Gelovani, “Ultrasound and multispectral photoacoustic systems and methods for brain and spinal cord imaging through acoustic windows,” Google Patents: 2020.

8. K. Kratkiewicz, R. Manwar, A. Rajabi-Estarabadi, J. Fakhoury, J. Meiliute, S. Daveluy, D. Mehregan, and K. M. Avanaki, “Photoacoustic/ultrasound/optical coherence tomography evaluation of melanoma lesion and healthy skin in a swine model,” Sensors 19, 2815 (2019). [CrossRef]  

9. K. Kratkiewicz, R. Manwar, Y. Zhou, M. Mozaffarzadeh, and K. Avanaki, “Technical considerations in the Verasonics research ultrasound platform for developing a photoacoustic imaging system,” Biomed. Opt. Express 12(2), 1050–1084 (2021). [CrossRef]  

10. J. I. Matchynski, R. Manwar, K. J. Kratkiewicz, R. Madangopal, V. A. Lennon, K. M. Makki, A. L. Reppen, A. R. Woznicki, B. T. Hope, and S. A. Perrine, “Direct measurement of neuronal ensemble activity using photoacoustic imaging in the stimulated Fos-LacZ transgenic rat brain: a proof-of-principle study,” Photoacoustics 24, 100297 (2021). [CrossRef]  

11. M. Nasiriavanaki, J. Xia, H. Wan, A. Q. Bauer, J. P. Culver, and L. V. Wang, “High-resolution photoacoustic tomography of resting-state functional connectivity in the mouse brain,” Proc. Natl. Acad. Sci. 111(1), 21–26 (2014). [CrossRef]  

12. H. C. Glass, K. B. Nash, S. L. Bonifacio, A. J. Barkovich, D. M. Ferriero, J. E. Sullivan, and M. R. Cilio, “Seizures and magnetic resonance imaging–detected brain injury in newborns cooled for hypoxic-ischemic encephalopathy,” The Journal of pediatrics 159(5), 731–735.e1 (2011). [CrossRef]  

13. W. L. Lanier, “Cerebral metabolic rate and hypothermia: their relationship with ischemic neurologic injury,” Journal of neurosurgical anesthesiology 7(3), 216–221 (1995). [CrossRef]  

14. B. E. Lingwood, K. R. Dunster, G. N. Healy, L. C. Ward, and P. B. Colditz, “Cerebral impedance and neurological outcome following a mild or severe hypoxic/ischemic episode in neonatal piglets,” Brain Res. 969(1-2), 160–167 (2003). [CrossRef]  

15. T. Durduran, R. Choe, W. B. Baker, and A. G. Yodh, “Diffuse optics for tissue monitoring and tomography,” Rep. Prog. Phys. 73(7), 076701 (2010). [CrossRef]  

16. A. T. Eggebrecht, S. L. Ferradal, A. Robichaux-Viehoever, M. S. Hassanpour, H. Dehghani, A. Z. Snyder, T. Hershey, and J. P. Culver, “Mapping distributed brain function and networks with diffuse optical tomography,” Nat. Photonics 8(6), 448–454 (2014). [CrossRef]  

17. A. P. Gibson, T. Austin, N. Everdell, M. Schweiger, S. R. Arridge, J. Meek, J. Wyatt, D. T. Delpy, and J. C. Hebden, “Three-dimensional whole-head optical tomography of passive motor evoked responses in the neonate,” NeuroImage 30(2), 521–528 (2006). [CrossRef]  

18. A. Kleinschmidt, H. Obrig, M. Requardt, K.-D. Merboldt, U. Dirnagl, A. Villringer, and J. Frahm, “Simultaneous recording of cerebral blood oxygenation changes during human brain activation by magnetic resonance imaging and near-infrared spectroscopy,” J. Cereb. Blood Flow Metab. 16(5), 817–826 (1996). [CrossRef]  

19. K. J. Cash, C. Li, J. Xia, L. V. Wang, and H. A. Clark, “Optical drug monitoring: photoacoustic imaging of nanosensors to monitor therapeutic lithium in vivo,” ACS Nano 9(2), 1692–1698 (2015). [CrossRef]  

20. D. Wang, Y. Wu, and J. Xia, “Review on photoacoustic imaging of the brain using nanoprobes,” Neurophotonics 3(1), 010901 (2016). [CrossRef]  

21. L. V. Wang and J. Yao, “A practical guide to photoacoustic tomography in the life sciences,” Nat. Methods 13(8), 627–638 (2016). [CrossRef]  

22. S. Na and L. V. Wang, “Photoacoustic computed tomography for functional human brain imaging,” Biomed. Opt. Express 12(7), 4056–4083 (2021). [CrossRef]  

23. S. Na, J. J. Russin, L. Lin, X. Yuan, P. Hu, K. B. Jann, L. Yan, K. Maslov, J. Shi, and D. J. Wang, “Massively parallel functional photoacoustic computed tomography of the human brain,” Nat. Biomed. Eng. 453, 1–9 (2021). [CrossRef]  

24. B. Chen, J. Yang, D. Wu, D. Zeng, Y. Yi, N. Yang, and H. Jiang, “Photoacoustic imaging of cerebral hypoperfusion during acupuncture,” Biomed. Opt. Express 6(9), 3225–3234 (2015). [CrossRef]  

25. M. Jeon, J. Kim, and C. Kim, “Multiplane spectroscopic whole-body photoacoustic imaging of small animals in vivo,” Med. Biol. Eng. Comput. 54(2-3), 283–294 (2016). [CrossRef]  

26. B. Ning, N. Sun, R. Cao, R. Chen, K. K. Shung, J. A. Hossack, J.-M. Lee, Q. Zhou, and S. Hu, “Ultrasound-aided multi-parametric photoacoustic microscopy of the mouse brain,” Sci. Rep. 5, 18775 (2015). [CrossRef]  

27. K. Sivasubramanian, V. Periyasamy, and M. Pramanik, “Non-invasive sentinel lymph node mapping and needle guidance using clinical handheld photoacoustic imaging system in small animal,” J. Biophotonics 11(1), e201700061 (2018). [CrossRef]  

28. C. Wood, K. Harutyunyan, J. De La Cerda, C. Kaffes, N. Z. Millward, S. Shanmugavelandy, M. Konopleva, and R. Bouchard, “Assessment of blood oxygen saturation using spectroscopic photoacoustic imaging as a biomarker for disease progression in a small-animal leukemia model,” In Proceedings of Medical Imaging 2018: Ultrasonic Imaging and Tomography; p. 105800W.

29. G. Xu, L. A. Johnson, J. Hu, J. R. Dillman, P. D. Higgins, and X. Wang, “Detecting inflammation and fibrosis in bowel wall with photoacoustic imaging in a Crohn's disease animal model,” In Proceedings of Photons Plus Ultrasound: Imaging and Sensing (2015), p. 932347.

30. L.-D. Liao, M.-L. Li, H.-Y. Lai, Y.-Y. I. Shih, Y.-C. Lo, S. Tsang, P. C.-P. Chao, C.-T. Lin, F.-S. Jaw, and Y.-Y. Chen, “Imaging brain hemodynamic changes during rat forepaw electrical stimulation using functional photoacoustic microscopy,” NeuroImage 52(2), 562–570 (2010). [CrossRef]  

31. M. Petri, I. Stoffels, J. Jose, J. Leyh, A. Schulz, J. Dissemond, D. Schadendorf, and J. Klode, “Photoacoustic imaging of real-time oxygen changes in chronic leg ulcers after topical application of a haemoglobin spray: a pilot study,” Journal of Wound Care 25(2), 87–91 (2016). [CrossRef]  

32. J. Reber, M. Willershäuser, A. Karlas, K. Paul-Yuan, G. Diot, D. Franz, T. Fromme, S. V. Ovsepian, N. Bézière, and E. Dubikovskaya, “Non-invasive measurement of brown fat metabolism based on optoacoustic imaging of hemoglobin gradients,” Cell Metab. 27(3), 689–701.e4 (2018). [CrossRef]  

33. M. Toi, Y. Asao, Y. Matsumoto, H. Sekiguchi, A. Yoshikawa, M. Takada, M. Kataoka, T. Endo, N. Kawaguchi-Sakita, and M. Kawashima, “Visualization of tumor-related blood vessels in human breast by photoacoustic imaging system with a hemispherical detector array,” Sci. Rep. 7, 41970 (2017). [CrossRef]  

34. J. Weber, P. C. Beard, and S. E. Bohndiek, “Contrast agents for molecular photoacoustic imaging,” Nat. Methods 13, 639 (2016). [CrossRef]  

35. N. Chang, 3rd W. Goodson, F. Gottrup, and T. K. C. Hunt, “Direct measurement of wound and tissue oxygen tension in postoperative patients,” Annals of surgery 197, 470 (1983). [CrossRef]  

36. C. P. Chao, C. G. Zaleski, and A. C. Patton, “Neonatal hypoxic-ischemic encephalopathy: multimodality imaging findings,” RadioGraphics 26(suppl_1), S159–S172 (2006). [CrossRef]  

37. H. W. Hopf, T. K. Hunt, J. M. West, P. Blomquist, W. H. Goodson, J. A. Jensen, K. Jonsson, P. B. Paty, J. M. Rabkin, and R. A. Upton, “Wound tissue oxygen tension predicts the risk of wound infection in surgical patients,” Arch. Surg. 132(9), 997–1004 (1997). [CrossRef]  

38. D. S. Vikram, J. L. Zweier, and P. Kuppusamy, “Methods for noninvasive imaging of tissue hypoxia,” Antioxid. Redox Signaling 9(10), 1745–1756 (2007). [CrossRef]  

39. J. Yao, L. Wang, J.-M. Yang, K. I. Maslov, T. T. Wong, L. Li, C.-H. Huang, J. Zou, and L. V. Wang, “High-speed label-free functional photoacoustic microscopy of mouse brain in action,” Nat. Methods 12, 407 (2015). [CrossRef]  

40. A. Hariri, J. Wang, Y. Kim, A. Jhunjhunwala, D. L. Chao, and J. V. Jokerst, “In vivo photoacoustic imaging of chorioretinal oxygen gradients,” J. Biomed. Opt. 23(03), 1 (2018). [CrossRef]  

41. A.-R. Mohammadi-Nejad, M. Mahmoudzadeh, M. S. Hassanpour, F. Wallois, O. Muzik, C. Papadelis, A. Hansen, H. Soltanian-Zadeh, J. Gelovani, and M. Nasiriavanaki, “Neonatal brain resting-state functional connectivity imaging modalities,” Photoacoustics 10, 1–19 (2018). [CrossRef]  

42. J. Laufer, D. Delpy, C. Elwell, and P. Beard, “Quantitative spatially resolved measurement of tissue chromophore concentrations using photoacoustic spectroscopy: application to the measurement of blood oxygenation and haemoglobin concentration,” Phys. Med. Biol. 52(1), 141–168 (2007). [CrossRef]  

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

44. K. Maslov, H. F. Zhang, and L. V. Wang, “Effects of wavelength-dependent fluence attenuation on the noninvasive photoacoustic imaging of hemoglobin oxygen saturation in subcutaneous vasculature in vivo,” Inverse Problems 23(6), S113–S122 (2007). [CrossRef]  

45. R. Hochuli, L. An, P. C. Beard, and B. T. Cox, “Estimating blood oxygenation from photoacoustic images: can a simple linear spectroscopic inversion ever work?” J. Biomed. Opt. 24(12), 1 (2019). [CrossRef]  

46. M. A. Naser, D. R. Sampaio, N. M. Muñoz, C. A. Wood, T. M. Mitcham, W. Stefan, K. V. Sokolov, T. Z. Pavan, R. Avritscher, and R. R. Bouchard, “Improved photoacoustic-based oxygen saturation estimation with SNR-regularized local fluence correction,” IEEE Trans. Med. Imaging 38(2), 561–571 (2019). [CrossRef]  

47. W. C. Vogt, X. Zhou, R. Andriani, K. A. Wear, T. J. Pfefer, and B. S. Garra, “Photoacoustic oximetry imaging performance evaluation using dynamic blood flow phantoms with tunable oxygen saturation,” Biomed. Opt. Express 10(2), 449–464 (2019). [CrossRef]  

48. A. Hussain, W. Petersen, J. Staley, E. Hondebrink, and W. Steenbergen, “Quantitative blood oxygen saturation imaging using combined photoacoustics and acousto-optics,” Opt. Lett. 41(8), 1720–1723 (2016). [CrossRef]  

49. J. Gröhl, T. Kirchner, T. Adler, and L. Maier-Hein, “Estimation of blood oxygenation with learned spectral decoloring for quantitative photoacoustic imaging (LSD-qPAI),” arXiv preprint arXiv:1902.05839 2019.

50. C. Cai, K. Deng, C. Ma, and J. Luo, “End-to-end deep neural network for optical inversion in quantitative photoacoustic imaging,” Opt. Lett. 43(12), 2752–2755 (2018). [CrossRef]  

51. S. Tzoumas, A. Nunes, I. Olefir, S. Stangl, P. Symvoulidis, S. Glasl, C. Bayer, G. Multhoff, and V. Ntziachristos, “Eigenspectra optoacoustic tomography achieves quantitative blood oxygenation imaging deep in tissues,” Nat. Commun. 7(1), 12121 (2016). [CrossRef]  

52. Z. Yuan and H. Jiang, “Simultaneous recovery of tissue physiological and acoustic properties and the criteria for wavelength selection in multispectral photoacoustic tomography,” Opt. Lett. 34(11), 1714–1716 (2009). [CrossRef]  

53. M. Sivaramakrishnan, K. Maslov, H. F. Zhang, G. Stoica, and L. V. Wang, “Limitations of quantitative photoacoustic measurements of blood oxygenation in small vessels,” Phys. Med. Biol. 52(5), 1349–1361 (2007). [CrossRef]  

54. J. Xiao, Z. Yuan, J. He, and H. Jiang, “Quantitative multispectral photoacoustic tomography and wavelength optimization,” J. X-Ray Sci. Technol. 18(4), 415–427 (2010). [CrossRef]  

55. S. Brigadoi, P. Aljabar, M. Kuklisova-Murgasova, S. R. Arridge, and R. J. Cooper, “A 4D neonatal head model for diffuse optical imaging of pre-term to term infants,” NeuroImage 100, 385–394 (2014). [CrossRef]  

56. S. L. Jacques, “Optical properties of biological tissues: a review,” Phys. Med. Biol. , 58(11), R37–R61 (2013). [CrossRef]  

57. B. K. Hendricks, J. Hartman, and A. A. Cohen-Gadol, “Cerebrovascular operative anatomy: an immersive 3D and virtual reality description,” Operative Neurosurgery 15(6), 613–623 (2018). [CrossRef]  

58. J. H. Helthuis, T. P. Van Doormaal, B. Hillen, R. L. Bleys, A. A. Harteveld, J. Hendrikse, A. Van der Toorn, M. Brozici, J. J. Zwanenburg, and A. Van der Zwan, “Branching Pattern of the Cerebral Arterial Tree,” Anat Rec 302(8), 1434–1446 (2019). [CrossRef]  

59. D. Purves, G. Augustine, D. Fitzpatrick, L. Katz, A. LaMantia, J. McNamara, and S. Williams, “The blood supply of the brain and spinal cord,” Neuroscience 2, (2001).

60. W. B. DeLONG, “Anatomy of the middle cerebral artery: the temporal branches,” Stroke 4(3), 412–418 (1973). [CrossRef]  

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

62. N. Baik, B. Urlesberger, B. Schwaberger, G. M. Schmölzer, L. Mileder, A. Avian, and G. Pichler, “Reference ranges for cerebral tissue oxygen saturation index in term neonates during immediate neonatal transition after birth,” Neonatology 108(4), 283–286 (2015). [CrossRef]  

63. M. M. Tisdall, C. Taylor, I. Ilias Tachtsidis, T. S. Leung, C. Pritchard, C. E. Elwell, and M. Smith, “The effect on cerebral tissue oxygenation index of changes in the concentrations of inspired oxygen and end-tidal carbon dioxide in healthy adult volunteers,” Anesthesia and Analgesia 109, 906 (2009). [CrossRef]  

64. A. Custo, W. M. Wells Iii, A. H. Barnett, E. M. Hillman, and D. A. Boas, “Effective scattering coefficient of the cerebral spinal fluid in adult head models for diffuse optical imaging,” Appl. Opt. 45(19), 4747–4755 (2006). [CrossRef]  

65. H. Mitchell, T. Hamilton, F. Steggerda, and H. Bean, “The chemical composition of the adult human body and its bearing on the biochemistry of growth,” J. Biol. Chem. 158(3), 625–637 (1945). [CrossRef]  

66. S. Ijichi, T. Kusaka, K. Isobe, K. Okubo, K. Kawada, M. Namba, H. Okada, T. Nishida, T. Imai, and S. Itoh, “Developmental changes of optical properties in neonates determined by near-infrared time-resolved spectroscopy,” Pediatr. Res. 58(3), 568–573 (2005). [CrossRef]  

67. J. Zhao, H. S. Ding, X. Hou, C. Zhou, and B. Chance, “In vivo determination of the optical properties of infant brain using frequency-domain near-infrared spectroscopy,” J. Biomed. Opt. 10(2), 024028 (2005). [CrossRef]  

68. N. Choudhury, R. Samatham, and S. L. Jacques, “Linking visual appearance of skin to the underlying optical properties using multispectral imaging,” In Proceedings of Photonic Therapeutics and Diagnostics VI, p. 75480G.

69. B. Weaver, G. Staddon, and M. Pearson, “Tissue blood content in anaesthetised sheep and horses,” Comparative Biochemistry and Physiology Part A: Physiology 94(3), 401–404 (1989). [CrossRef]  

70. J. S. Harrison, P. Rameshwar, V. Chang, and P. Bandari, “Oxygen saturation in the bone marrow of healthy volunteers,” Blood, The Journal of the American Society of Hematology 99, 394 (2002). [CrossRef]  

71. ANSI, “American National Standard for safe use of lasers (ANSI 136.1),” ANSI 136.1-2007 (The Laser Institute of America, 2007).

72. M. Mozaffarzadeh, A. Mahloojifar, M. Orooji, S. Adabi, and M. Nasiriavanaki, “Double-stage delay multiply and sum beamforming algorithm: Application to linear-array photoacoustic imaging,” IEEE Trans. Biomed. Eng. 65(1), 31–42 (2018). [CrossRef]  

73. P. Omidi, M. Zafar, M. Mozaffarzadeh, A. Hariri, X. Haung, M. Orooji, and M. Nasiriavanaki, “A novel dictionary-based image reconstruction for photoacoustic computed tomography,” Appl. Sci. 8(9), 1570 (2018). [CrossRef]  

74. M. J. Siegel, G. D. Shackelford, J. M. Perlman, and K. H. Fulling, “Hypoxic-ischemic encephalopathy in term infants: diagnosis and prognosis evaluated by ultrasound,” Radiology 152(2), 395–399 (1984). [CrossRef]  

75. B. Liang, S. Wang, F. Shen, Q. H. Liu, Y. Gong, and J. Yao, “Acoustic impact of the human skull on transcranial photoacoustic imaging,” Biomed. Opt. Express 12(3), 1512–1528 (2021). [CrossRef]  

76. N. Meimani, N. Abani, J. Gelovani, and M. R. Avanaki, “A numerical analysis of a semi-dry coupling configuration in photoacoustic computed tomography for infant brain imaging,” Photoacoustics 7, 27–35 (2017). [CrossRef]  

77. R. Manwar, K. Kratkiewicz, and K. Avanaki, “Overview of ultrasound detection technologies for photoacoustic imaging,” Micromachines 11(7), 692 (2020). [CrossRef]  

78. R. Manwar, K. Kratkiewicz, and K. Avanaki, “Investigation of the effect of the skull in transcranial photoacoustic imaging: a preliminary ex vivo study,” Sensors 20(15), 4189 (2020). [CrossRef]  

79. R. Manwar, X. Li, S. Mahmoodkalayeh, E. Asano, D. Zhu, and K. Avanaki, “Deep learning protocol for improved photoacoustic brain imaging,” J. Biophotonics 13(10), e202000212 (2020). [CrossRef]  

80. L. Mohammadi, H. Behnam, J. Tavakkoli, and K. Avanaki, “Skull acoustic aberration correction in photoacoustic microscopy using a vector space similarity model: a proof-of-concept simulation study,” Biomed. Opt. Express 11(10), 5542–5556 (2020). [CrossRef]  

81. L. Mohammadi, H. Behnam, J. Tavakkoli, and M. Avanaki, “Skull’s photoacoustic attenuation and dispersion modeling with deterministic ray-tracing: Towards real-time aberration correction,” Sensors 19(2), 345 (2019). [CrossRef]  

Supplementary Material (1)

NameDescription
Supplement 1       Supplement 1

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (11)

Fig. 1.
Fig. 1. (a) Head model: sagittal view of a 37 weeks neonate atlas used for simulations. ECT: extra-cranial tissue, CSF: cerebrospinal fluid. The areas with different colors show regions of hypoxia. (b) At a given experiment, the selected voxels in normal tissue and a region of hypoxia are used for analysis.
Fig. 2.
Fig. 2. (a) Map of initial pressure at 780nm, (b) map of initial pressure at 830nm, (c) map of estimated sO2 using the wavelength pair, 780nm and 830nm with $sO_{2,normal}^{true}$=0.7 and $Hypoxi{a^{true}}$=30%. The area inside the shaded region in c shows the hypoxic region and the area between two black borders contains the pixels of surrounding normal tissue which are involved in the calculation of hypoxia.
Fig. 3.
Fig. 3. (a) pixel intensity histogram of estimated sO2 for hypoxic and surrounding normal tissue located at a depth of 14mm with a wavelength-pair 780nm and 830nm. Dashed-red lines indicate the mean value of the normal region, $sO_{2,normal}^{estim}$, and that for the hypoxic region, $sO_{2,ROH}^{estim}$. (b) $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for 20 scenarios (four sO2 levels for normal tissue and five hypoxia) at a depth of 14mm for the wavelength-pair, 780nm, and 830nm. The $Hypoxi{a^{estim}}$ calculated for the case presented in Fig. 3(a) is specified with a red dotted circle in Fig. 3(b).
Fig. 4.
Fig. 4. Maximum permissible exposure graph according to ANSI limit (black line), taken from [71]; black arrows show MPE for wavelengths 780nm and 830nm; the red astros show discrete levels of energy that are used in energy optimization.
Fig. 5.
Fig. 5. (a) pixel intensity histogram of estimated sO2 after energy optimization calculated for hypoxic and surrounding normal tissue located at the depth of 14mm with a wavelength-pair, 780nm and 830nm. Dashed red lines indicate the mean value of the normal region, $sO_{2,normal}^{estim}$, and that for the hypoxic region, $sO_{2,ROH}^{estim}$. (b) $Hypoxi{a^{estim}}$ versus $Hypoxi{a^{true}}$ for 20 scenarios (four sO2 levels for normal tissue and five levels for hypoxic regions) at the depth of 14mm for the wavelength-pair, 780nm and 830nm after energy optimization. The $Hypoxi{a^{estim}}$ calculated for the case presented in Fig. 5(a) is specified with a red dotted circle in Fig. 5(b).
Fig. 6.
Fig. 6. Hypoxia estimation error calculated at wavelength-pairs in the wavelength range 500 nm to1000 nm (with steps of 10 nm) without energy optimization (top triangle) and with energy optimization (bottom triangle) for depths (a)2 mm, (b)6 mm, (c)10 mm, (d)14 mm, (e)18 mm and (f)22 mm.
Fig. 7.
Fig. 7. Map of energy ratio E(λ2)/E(λ1) at wavelength-pairs (λ1,λ2) in the wavelength range 500 nm to 1000 nm (with steps of 10 nm) for depths (a) 2 mm, (b) 6 mm, (c) 10 mm, (d) 14 mm, (e) 18 mm and (f) 22 mm. E(λ) represents the pulse energy of the selected wavelength, λ.
Fig. 8.
Fig. 8. Boxplot of initial pressure in both normal and hypoxic regions for a special case ($sO_{2,normal}^{true}$=0.7 and $Hypoxi{a^{true}}$=30%) at depth 14mm for different wavelengths. Black dots indicate median values of initial pressures. Red dashed line shows the NEP level.
Fig. 9.
Fig. 9. Consideration of transducer NEP and sensitivity in hypoxia estimation errors calculated at wavelength-pairs in the wavelength range 500 nm to1000 nm (with steps of 10 nm) with energy optimization for depths (a) 2 mm, (b)6 mm, (c)10 mm, (d)14 mm, (e)18 mm and (f) 22 mm. Grey areas in all maps are the areas where at least one of the detectability conditions has not been fulfilled.
Fig. 10.
Fig. 10. Linear regressions calculated for estimated hypoxia at the depth 14 mm with optimized wavelength and energy, for the wavelength-pair 800 nm and 830 nm. Different markers correspond to different $sO_{2,normal}^{true}$ and red dashed lines are linear fits where β is the slope and R2 is the coefficient of determination.
Fig. 11.
Fig. 11. Changes in the error of estimated hypoxia at different depths when absorption coefficient is changed by +1%(blue), +2%(red), +5(cyan), -1%(green), -20%(pink) and -5(orange) relative to the original values of absorption coefficient considered for (a) ECT layer and (b) brain tissue. Dashed black line shows the original error. The errors are computed for the best wavelength pairs with equal energy at different depths (see Supplement 1, Tables S1-S6).

Tables (5)

Tables Icon

Table 1. Constants used to generate ECT and brain tissue scattering coefficient. ECT: extracerebral tissue

Tables Icon

Table 2. Constants used in Eq. (7) to generate the absorption coefficient of ECT, CSF, normal brain tissue, and ROH. ECT: extracranial tissue, and CSF: cerebrospinal fluid

Tables Icon

Table 3. Optimized wavelength-pairs at depth=14 mm with corresponding Error (λ1, λ2) sorted in ascending order for equal and optimized energies without and with considering NEP

Tables Icon

Table 4. Slope and R2 of linear regression of estimated hypoxia for given normal tissue oxygen saturations at different depths with optimized wavelength and energy.

Tables Icon

Table 5. Slope and R2 of linear regression of estimated hypoxia for given normal tissue oxygen saturations at different depths with optimized wavelength and energy.

Equations (10)

Equations on this page are rendered with MathJax. Learn more.

p 0 = η Γ A = η Γ μ a F
s O 2 e x a c t = [ H b O ] [ H b T ] = [ H b O ] [ H b O ] + [ H b R ] = ε H b R ( λ 2 ) μ a ( λ 1 ) ε H b R ( λ 1 ) μ a ( λ 2 ) Δ ε H b ( λ 2 ) μ a ( λ 1 ) Δ ε H b ( λ 1 ) μ a ( λ 2 )
s O 2 e x a c t = [ H b O ] [ H b O ] + [ H b R ] = ε H b R ( λ 2 ) p 0 ( λ 1 ) / F ( λ 1 ) ε H b R ( λ 1 ) p 0 ( λ 2 ) / F ( λ 2 ) Δ ε H b ( λ 2 ) p 0 ( λ 1 ) / F ( λ 1 ) Δ ε H b ( λ 1 ) p 0 ( λ 2 ) / F ( λ 2 )
s O 2 e s t i m = ε H b ( λ 2 ) p 0 ( λ 1 ) ε H b ( λ 1 ) p 0 ( λ 2 ) Δ ε H b ( λ 2 ) p 0 ( λ 1 ) Δ ε H b ( λ 1 ) p 0 ( λ 2 )
H y p o x i a t r u e = ( 1 s O 2 , R O H t r u e s O 2 , n o r m a l t r u e ) × 100
μ s ( λ ) = a ( f R a y ( λ 500 ( n m ) ) 4 + ( 1 f R a y ) ( λ 500 ( n m ) ) b M i e )
μ a = B ( s O 2 t r u e ) μ a , H b O + B ( 1 ( s O 2 t r u e ) ) μ a , H b R + W μ a , w a t e r + F μ a , f a t + M μ a , m e l a n o s o m e
H y p o x i a e s t i m ( λ 1 , λ 2 ) = ( 1 s O 2 , R O H e s t i m ( λ 1 , λ 2 ) s O 2 , n o r m a l e s t i m ( λ 1 , λ 2 ) ) × 100
E r r o r ( λ 1 , λ 2 ) = 1 20 a l l s O 2 , n o r m a l t r u e a l l h y p o x i a t r u e | H y p o x i a e s t i m ( λ 1 , λ 2 ) H y p o x i a t r u e H y p o x i a t r u e | × 100
M P E = 2 C A × 10 2 ( J / c m 2 ) , C A = { 1.0 0.400 < λ < 0.700 ( μ m ) 10 2 ( λ 0.700 ) 0.700 < λ < 1.050 ( μ m )
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.