1.

Introduction

The endonasal transsphenoidal approach to pituitary tumor resection is a minimally invasive technique that requires insertion of surgical instruments through the nose to remove sphenoid bone and underlying pituitary tumors.¹ Although the procedure is generally safe,² morbidity and mortality rates rise to 14% to 23% and 24% to 26%, respectively, if iatrogenic injury to the internal carotid arteries (ICAs) occurs.^3–5 Current intraoperative guidance techniques, such as stereotactic guidance and endoscopy, enable monitoring of the ICAs in close proximity to the surgical site; however, they suffer from two primary limitations. First, stereotactic guidance is subject to registration errors, which can become increasingly large as patient anatomy is disrupted during surgery and deviates from the anatomy in preoperative x-ray computed tomography (CT) or magnetic resonance images. Second, endoscopy is unable to identify the ICAs when they are obscured by bone or other tissues in the operative path. Our group is investigating the use of transcranial photoacoustic imaging as an intraoperative imaging technique for real-time visualization of the ICAs to address these two limitations, as detailed in our original research papers ^6–11 and summarized in the literature surveys on this topic.^12,13

To achieve photoacoustic imaging for guidance of endonasal transsphenoidal surgery, we propose the insertion of light-transmitting fiber optic devices in the nasal cavity, similar to other surgical tools.^6–12 This light source will then excite the hemoglobin within the ICAs, converting the absorbed optical energy to acoustic energy that is received by an externally placed ultrasound receiver. The external ultrasound receiver placement results in a transcranial photoacoustic imaging scenario, which is challenged by acoustic interactions with bone and is known to degrade image quality.^14–16 Previous work from our group developed and demonstrated a simulation method to identify naturally occurring acoustic windows in the adult human skull to minimize acoustic interactions with bone and to provide high-contrast photoacoustic images of the ICAs.^10,11 We demonstrated that patient-specific simulations have the potential to enable preoperative planning to determine appropriate placement of imaging system components.^10,11

Our vision for a presurgical workflow with the use of simulations is shown in Fig. 1, in direct comparison with a workflow without access to presurgical simulations. Without simulations, a surgeon may need to search and find optimal locations to place an ultrasound receiver to best receive the acoustic signals within the skull. This intraoperative process may be prone to the expenditure of time that could instead be spent operating on the patient and would increase the time that the patient is under anesthesia and on the operating table. In addition, the surgeon may have to abandon the use of intraoperative photoacoustic image guidance during the procedure if a viable receiver location is not determined in time, thereby sacrificing the benefits of real-time photoacoustic-based guidance during the surgery.¹⁰ In contrast, the simulations that we propose have the potential to preoperatively identify patient-specific ultrasound receiver locations.¹⁰ The surgeon may then use the simulation results and the patient’s preoperative CT scan to construct a time-efficient surgical plan for each operation. With this step completed, the neurosurgery may then be safely executed with the benefit of real-time, intraoperative photoacoustic image guidance.

Fig. 1

Surgical workflow with and without the use of patient-specific simulations. Without photoacoustic simulations, a surgeon will likely have to test multiple ultrasound receiver locations and their feasibility. Simulations have the potential to efficiently identify these locations before the surgical procedure.

Our group previously identified the ocular cavity, temporal region, and nasal cavity as three optimal ultrasound receiver locations, using compressional-wave-only simulations.^10,11 However, the shear waves known to propagate within dense media such as bone were excluded from our initial demonstrations. As an alternative, elastic wave simulations, which include both compressional and shear waves, are expected to more accurately represent the physical acoustic process. However, elastic wave simulations, such as the k-Wave elastic simulations based on the classical Kelvin–Voigt absorption model, are time consuming and memory intensive.¹⁷ Time and memory costs reduce the ease and likelihood of clinical translation.

This paper presents a comparison of compressional and elastic photoacoustic k-Wave simulations to investigate which simulation type (i.e., compressional or elastic) is required for preoperative surgical planning. The remainder of this paper is organized as follows. Section 2 details the simulation methods and quantitative metrics used for the comparison. Section 3 presents the resulting comparisons. Section 4 discusses the implications of the results with respect to the compressional wave simulations that were used in our experimental validation studies, with regard to the vision shown in Fig. 1 and with an eye toward reducing barriers to clinical translation. Finally, Sec. 5 concludes the paper.

2.

Methods

2.1.

Simulation Configuration

Three-dimensional (3D) photoacoustic k-Wave^17–19 simulations were performed after converting the CT volume of a human cadaver skull into heterogeneous, volumetric maps of corresponding density, compressional and shear wave sound speeds, and compressional and shear wave absorption prefactors. Figure 2 shows an example axial slice from the CT volume of a human cadaver skull. Homogeneous volumes were additionally modeled as the average density, sound speed, and absorption values of brain tissue to provide baseline photoacoustic images that do not contain heterogeneous tissue effects, such as aberration, attenuation, scattering, and reverberation. The simulated tissue properties for the heterogeneous and homogeneous cases are reported in Table 2. A two-dimensional (2D) cross-section of the propagating wave at various time points for the homogeneous and heterogeneous cases is shown in Fig. 3, with more details available in Video 1. The computational grid was defined with a symmetric voxel size of $0.3\mathrm{mm}\times 0.3\mathrm{mm}\times 0.3\mathrm{mm}$. Acoustic wave propagation was simulated with a sampling frequency of 82 MHz (i.e., time increment of ${1.12}^{-8}\mathrm{s}$) on a NVIDIA Quadro RTX 6000 GPU (Table 1).

Fig. 2

Axial slice from the CT volume of the human cadaver skull demonstrating a 2D cross-section of the 3D simulation configuration. Spherical photoacoustic sources were placed within the LCA and at distances of 6 to 13 mm from the LCA to represent the tip of a surgical instrument (a distance of 6 mm is shown). Green lines illustrate locations of independently placed ultrasound transducers.

Fig. 3

2D cross-section of acoustic wave propagation at various time points (i.e., $t_0$ through $t_3$) during (a) the homogeneous and (b) heterogeneous 3D simulations. There is no bone present in the homogeneous simulation, as represented by the grayscale background. In both cases, the acoustic wave propagates spherically outward from the initial pressure distribution (location of left ICA). In the heterogeneous simulations, acoustic interactions with cranial bone cause distortions (i.e., aberrations, attenuation, scattering, and reverberations) in the waveform, which ultimately degrade image quality. See Video 1 for the evolution of wave propagation for these two simulations (Video 1: MP4, 0.734 MB [URL: https://doi.org/10.1117/1.JBO.26.7.076006.1]).

Table 1

Transcranial photoacoustic simulation parameters.

Simulation parameter	Values
Transducer
Pitch	0.3 mm
Height	13.6 mm
Kerf	0.0 mm
Elements	64
Bandwidth	84.5%
Center frequency	3.6 MHz
Field of view	90 deg
Simulation grid
Voxel size	$0.3\times 0.3\times 0.3{\mathrm{mm}}^3$
Sampling frequency	82 MHz

Table 2

Simulated tissue properties of the homogeneous and heterogeneous volumetric maps.18,20–25 The two shear wave properties of sound speed in brain and absorption power law prefactors in bone and brain are not explicitly known. However, these properties can be estimated as approximately half the corresponding compressional speed of sound and approximately double the corresponding compressional absorption power law prefactor.24,26

		Speed of sound (m/s)	Absorption power law prefactor19,27 (dB MHz1.18/cm)	Density (kg/m3)
Bone	Compression	1300 to 349218,20	1.36 to 2.5022	812 to 277018,21
	Shear	650 to 174624	2.71 to 5.00a	812 to 277018,21
Brain	Compression	151925	0.2122,23	104525
	Shear	800a	0.75a	104525

^aEstimated shear wave values.

Phased array ultrasound transducers with 0.3 mm pitch, 13.5 mm height, 0 mm kerf, and 64 elements (as reported in Table 1) were positioned to receive transcranial photoacoustic signals from three acoustic windows: (1) the left ocular cavity, (2) the left temple region, and (3) the nasal cavity, as shown in Fig. 2. A spherical photoacoustic source with a diameter of 0.3 mm (i.e., a point source) or 4 mm (i.e., the diameter of a typical adult carotid artery) was placed in the location of the left carotid artery (LCA), as shown in Fig. 2. Background absorption was not modeled in the photoacoustic source distributions for the following two reasons. First, background absorption can be mitigated by carefully selecting the illumination wavelength to preferentially excite the oxygenated hemoglobin over surrounding tissues.²⁸ Second, our previous observations in transcranial photoacoustic imaging of the cadaveric ICAs did not display evidence of background absorption.¹⁰

To investigate photoacoustic image characteristics when the tip of a surgical instrument is also in the imaging plane, a 2-mm-diameter spherical photoacoustic source was positioned within the surgical site at distances of 6 to 13 mm from the 4-mm-diameter LCA target. These distances were measured from the centers of the LCA and instrument targets. Considering this measurement and the best possible system resolution described in Sec. 2.2, the minimum distance simulated was 6 mm. The distances of the LCA source relative to the center of the transducers located on the ocular, nasal, and temple regions were 7.40, 8.02, and 5.03 cm, respectively.

3.

Results

Table 3 compares the execution time and memory usage of the 3D compressional and elastic simulations, performed with the three independently placed transducers shown in Fig. 2. The last row of this table reports the mean time and memory usage across the three simulations to demonstrate that elastic simulations required $\sim 4.6\times$ more time and $2.7\times$ more memory than compressional simulations.

Table 3

Comparison of execution times and memory usage for the compressional and elastic k-Wave simulations.

Figures 4(a) and 4(b) show simulated photoacoustic channel data from the heterogeneous compressional and elastic simulations of the 0.3-mm-diameter LCA source obtained from the temple acoustic window. Although the initial wavefront shapes across elements look similar between the two simulation types, the compressional simulation has greater reverberations in the photoacoustic signal subsequent to the initial wavefront. To further appreciate the differences in the wavefronts between the simulation types, Fig. 4(c) shows the difference image of the elastic channel data subtracted from the compressional channel data. The wavefront in the compressional simulation has a greater amplitude than the wavefront in the elastic simulation, as indicated by the positive values along the wavefront in the difference image. Specifically, the maximum amplitudes of the wavefronts were $9{\mathrm{e}}^{-5}$ and $6{\mathrm{e}}^{-5}$ Pa for the compressional and elastic simulations, respectively.

Fig. 4

Simulated (a) compressional and (b) elastic photoacoustic channel data of the LCA obtained with the temporal acoustic window. (c) Difference image of the elastic channel data subtracted from the compressional channel data.

Figure 5(a) shows simulated photoacoustic images generated from the heterogeneous simulations of the 4-mm-diameter LCA source. The shape, location, and visibility of the targets generally agree between the compressional and elastic simulation pairs for each transducer location. Figure 5(b) quantifies the target visibility and detectability of the heterogeneous simulated images displayed in Fig. 5(a). There were minimal differences in contrast, SNR, and CNR measurements, and the gCNR measurements were equivalent. Therefore, the visibility and detectability of these photoacoustic image targets ranges from similar to equivalent for compressional and elastic wave simulations. As noted by Kempski et al.,³³ there is no benefit to improving particular methods when the gCNR is already near its maximum value of 1.0, and in this case, the method is the type of simulation used for presurgical planning.

Fig. 5

(a) Simulated photoacoustic images of the LCA obtained with the ocular, nasal, and temporal acoustic windows, from left to right, respectively. The top and bottom rows show compressional and elastic simulation results, respectively. The single target and multiple background ellipsoidal ROI are outlined in blue and orange, respectively. Images are displayed with the 10 dB dynamic range. (b) Corresponding mean ± one standard deviation of contrast, SNR, CNR, and gCNR^32,33 measurements.

Figures 6(a) and 6(b) show the $-6\mathrm{dB}$ contours of the point target images from the homogeneous and heterogeneous simulations, respectively. The best resolution of the system was 1.4 mm, obtained from the shortest cross-section of the $-6\mathrm{dB}$ contour of the homogeneous compressional simulation image obtained from the temporal region. Figure 6(c) shows the calculated areas of the contours. The contour areas of the baseline homogeneous simulated images were 9.33, 9.80, and $6.60{\mathrm{mm}}^2$ for the ocular, nasal, and temple acoustic windows, respectively, for compressional simulations and 10.27, 11.05, $7.08{\mathrm{mm}}^2$, respectively, for elastic simulations. The areas of the corresponding simulated images obtained with the heterogeneous skull model were 29.78, 38.26, and $24.22{\mathrm{mm}}^2$, respectively, for compressional simulations and 28.76, 35.29, and $24.21{\mathrm{mm}}^2$, respectively, for elastic simulations. Therefore, the resolution of these photoacoustic images is 0.45 to $1.24{\mathrm{mm}}^2$ better for the homogeneous compressional wave simulations than the homogeneous elastic wave simulations. However, in the presence of bone, the elastic wave simulations have 0.33 to $3.35{\mathrm{mm}}^2$ better resolution than compressional wave simulations for the heterogeneous case.

Fig. 6

Point spread function $-6\mathrm{dB}$ contours for the (a) homogeneous and (b) heterogeneous point target simulations. (c) Corresponding areas of the PSF $-6\mathrm{dB}$ contours. Homogeneous simulation measurements provide baseline resolution measurements without the negative effects of tissue heterogeneity.

Figure 7(a) shows an annotated photoacoustic image from a heterogeneous simulation of the LCA and a surgical instrument tip, obtained with the ocular acoustic window. The ground truth target-to-instrument distance is 11.27 mm, and the measured target-to-instrument distance is 11.18 mm in this image. Figure 7(b) shows image-based measurements of multiple target-to-instrument distances as a function of ground truth distances. The largest distance error between measurements and ground truth was obtained with the nasal cavity images, measuring 1.24 and 1.16 mm for the compressional and elastic wave simulations, respectively.

Fig. 7

(a) Heterogeneous simulated photoacoustic images of the two sources representing the LCA and instrument tip, obtained with the ocular acoustic window at a relative distance of 11.27 mm. The target-to-instrument distance is the length of the line between the source centers in the image. (b) Comparison of measured target-to-instrument distances as a function of the true target-to-instrument distances for compressional (○) and elastic ($\mathrm{\Delta }$) measurements in the ocular, nasal, and temporal acoustic windows. The ideal 1:1 relationship is shown as the dashed black line, with ± 1.4 mm boundaries, representing the best possible resolution of the system, indicated by the dashed gray lines.

The error ranges for the ocular and temporal regions were 0.09 to 0.91 mm and 0.14 to 0.88 mm for the compressional and elastic simulations, respectively. When comparing the target-to-instrument distance errors between the compressional and elastic simulations, the differences in errors were 0.03 to 0.09, 0.07 to 0.45, and 0.03 to 0.22 mm for the ocular, nasal, and temporal regions, respectively. Therefore, the relative target-to-instrument distances can be measured with accuracy within the 1.4 mm best possible resolution of the system from both compressional and elastic wave simulations.

4.

Discussion

This paper details our investigations of k-Wave elastic wave simulations as a method for presurgical planning of transducer placement during transcranial photoacoustic-guided neurosurgery. In particular, we assessed comparative differences in wave propagation, image quality (i.e., target visibility, target detectability, and target size estimation), and target-to-instrument tip distances (i.e., distances between a critical blood vessel and surgical instrument tip) using photoacoustic channel data and images generated from compressional and elastic transcranial simulations, resulting in three key observations. First, wavefronts in the channel data were similar in terms of general wavefront shape, but they differed in amplitude and reverberations between the compressional and elastic simulations, as shown in Fig. 4. Second, target visibility and detectability were qualitatively and quantitatively similar between the compressional and elastic simulations for the three acoustic windows tested, as shown in Fig. 5. Third, the accuracy of target-to-instrument tip distances measured using either compressional or elastic simulations was generally similar between elastic and compressional simulations and was within the best possible resolution of the system, as shown in Fig. 7. Although the latter two observations (based on results presented in Figs. 5 and 7) may seem counterintuitive with respect to shear waves being the primary acoustic component during acoustic propagation in bone, we hypothesize that the similarity of target visibility, target detectability, and target-to-instrument distances between the two simulation types is due to the strategic choice of transducer location ¹⁰ (i.e., the length of acoustic interaction with bone is minimized along the pathway from the LCA source to transducers) and the relatively minimal thickness of bone with respect to the distance between the source and each transducer.

Our three key observations indicate that it is likely sufficient to utilize the less time-consuming, less memory-intensive k-Wave compressional simulations for presurgical planning, confirming that appropriate complexity was included when comparing simulation to experimental results from a human cadaver.¹⁰ The arguably minimal 0.33 to $3.35{\mathrm{mm}}^2$ area differences in PSFs between k-Wave compressional and elastic simulations indicate that future investigations surrounding the effects of PSF differences could be useful prior to completely eliminating elastic simulations from the presurgical planning process. In particular, for the three acoustic windows tested, the areas of the $-6\mathrm{dB}$ contours of the PSFs achieved with the compressional wave simulations was up to $1.24{\mathrm{mm}}^2$ better than that of the elastic wave simulations for the homogeneous case and up to $3.35{\mathrm{mm}}^2$ worse than that of the elastic simulations for the heterogeneous case, as shown in Fig. 6. A neurosurgeon user study could possibly be performed to determine if these differences in image resolution would substantially alter a neurosurgeon’s presurgical plan.

Note that our observations are specifically catered to simulation results obtained with the k-Wave Toolbox software package. A comparison with multiple other simulation packages^34–36 was beyond the scope of this paper. Using the k-Wave simulation package, we demonstrated for the first time that it is likely sufficient to sacrifice the accuracy of elastic wave simulations by solely relying on compressional wave simulations for presurgical planning. Although this sacrifice affects the photoacoustic waveforms sensed by the simulated receivers (as shown in Fig. 4) and image resolution (as shown in Fig. 6), it may not impact the image quality needed to make surgical decisions based on target detectability and target-to-instrument distance. Thus, we now have a better understanding of this pathway to advance intraoperative transcranial photoacoustic imaging of the ICAs toward surgical use.

Regarding the benefit of shorter simulation times, though minimally invasive neurosurgeries are rarely emergency procedures, surgeons waiting on simulations to complete before operating on a patient might create bottlenecks in presurgical planning. This issue could be a barrier to entry for this technology in hospitals that might not have access to or financial resources for powerful GPUs. Without the GPUs that we used to accelerate our simulations, elastic wave simulations could take several hours to complete and ultimately disrupt the presurgical workflow. When utilizing compressional wave simulations, which require $\sim 4.6\times$ less time to simulate, surgeons can avoid these presurgical bottlenecks, regardless of their accessibility to GPUs. In addition, reduced computational memory places less strain on computational resource requirements for presurgical simulations and can enable surgeons to perform multiple patient-specific simulations at once. Thus, when implemented with compressional wave simulations, our proposed presurgical planning step has the potential to rapidly identify patient-specific optimal transducer locations for incorporation into the surgical plan, in support of our vision outlined in Fig. 1.

In summary, presurgical identification of optimal photoacoustic imaging system component locations with compressional wave simulations offers three primary benefits. First, this presurgical identification would eliminate wasting valuable operating room time to search for and find a suitable transducer location, thereby reducing the time that the patient is under anesthesia, the total procedure duration, and the medical cost. The second benefit is removal of the potential barrier that the surgeon may be unable to identify a viable transducer location and therefore unable to use photoacoustic image-guidance during the procedure. The third benefit is identification of transducer locations that would minimize image quality degradation from the presence of bone and thereby produce photoacoustic images of the ICAs with the best image quality possible for the patient. Possible future work includes a comparison of simulation types when modeling more complex source distributions, such as the possibility for optical absorption from the deoxygenated blood in the endothelial tubes of the cavernous sinus, and followup neurosurgeon user studies.

5.

Conclusions

The work presented in this paper is the first to reveal that the less time-consuming and less memory-intensive compressional k-Wave simulations are likely sufficient for identification of optimal transducer locations for transcranial photoacoustic-guided surgery. This assessment is based on target visibility and detectability (e.g., target contrast, SNR, CNR, and gCNR) and relative source-to-instrument distance producing similar or identical quantitative measurements for compressional and elastic simulations. These results have multiple implications for reducing barriers for clinical translation of simulation-based photoacoustic-guided surgery.