Coregistration and Spatial Compounding of Optoacoustic Cardiac Images via Fourier Analysis of Four-Dimensional Data

Abstract

Volumetric optoacoustic tomography has been shown to provide unprecedented capabilities for ultrafast imaging of cardiovascular dynamics in mice. Three-dimensional imaging rates in the order of 100 Hz have been achieved, which enabled the visualization of transient cardiac events such as arrhythmias or contrast agent perfusion without the need for retrospective gating. The fast murine heart rates (400–600 beats per minute) yet impose limitations when it comes to compounding of multiple frames or accurate registration of multi-spectral data. Herein, we investigate on the capabilities of Fourier analysis of four-dimensional data for coregistration of independent volumetric optoacoustic image sequences of the heart. The fundamental frequencies and higher harmonics of respiratory and cardiac cycles could clearly be distinguished, which facilitated efficient retrospective gating without additional readings. The performance of the suggested methodology was successfully demonstrated by compounding cardiac images acquired by raster-scanning of a spherical transducer array as well as by unmixing of oxygenated and deoxygenated hemoglobin from multi-spectral optoacoustic data.

1. Introduction

Cardiovascular imaging is essential for phenotyping the structure and function of the mammalian heart and vascular system. Preclinical studies generally rely on mice to model arrhythmias, myocardial infarction and other conditions [1,2,3]. However, imaging the mouse heart is challenged by the small size of internal structures and the very fast beating rate (400–600 beats per minute) [4]. Gating is used in high-resolution magnetic resonance imaging (MRI) and X-ray computed tomography (CT) to overcome the big difference between the beating frequency and the temporal resolution of these modalities [5,6]. This is achieved by either triggering the acquisition of data after detection of a physiological event (prospective gating) or by continuous acquisition and subsequent data correlation with the cardiac cycle (retrospective gating). Self-gated cardiac imaging methods are particularly valuable for small rodents as they do not require additional electrocardiogram (ECG) readings, but they generally lead to long acquisition times [7]. On the other hand, real-time imaging modalities have been used for visualizing the heart dynamics. Pulse-echo ultrasound (US) is widely used for in vivo cardiovascular imaging in rodents [4,8] and intravital microscopy has facilitated single cell tracking in cardiac tissues [9].

Recently, volumetric optoacoustic (OA) tomography (OAT) has been suggested as a particularly suitable modality for imaging the murine heart in vivo [10]. OAT provides high-contrast images of optical absorption within biological tissues with a single laser pulse [11,12,13], thus avoiding the need of image compounding. High frame rates in the order of 100 Hz have then been achieved for four-dimensional (real-time three-dimensional) imaging of the heart [14]. This has facilitated the visualization of transient events, such as arrhythmias [15] or the perfusion of a contrast agent [16]. By measuring the time of arrival of a bolus of indocyanine green (ICG) in the right and left ventricles, respectively, it was possible to estimate the pulmonary transit time (PTT), defined as the time it takes for blood to travel though the pulmonary circulation. It was shown that this physiological parameter changes in mice subjected to myocardial infarction [16,17] or chronic hypoxic stress [15], thus defining new capabilities in cardiovascular research. Fully exploiting the cardiac imaging capacities of volumetric OAT may however require acquisition of independent image sequences. For example, the effective field of view (FOV) scales with the achievable resolution and generally needs to be enhanced to cover the entire heart and surrounding regions, which can be done e.g., by stitching the images acquired by scanning the spherical array used for collecting OA signals [18]. Also, estimation of the bio-distribution of spectrally-distinctive substances such as oxygenated and deoxygenated hemoglobin requires acquisition of images corresponding to several excitation wavelengths. Synchronization of sequential acquisitions of cardiac images is challenged by the fast murine heart rate, and alternative methods must be developed for cardiac gating.

Herein, we investigate on the capabilities of Fourier analysis of four dimensional data as a self-gated cardiac imaging approach for co-registration of independent volumetric OA image sequences of the heart. We investigate on the feasibility of stitching heart images acquired at different scanning positions and multi-spectral un-mixing of cardiac data based on this methodology. We discuss on the potential implications in cardiovascular research and other fields.

2. Materials and Methods

2.1. Optoacoustic Imaging System

The OAT system used for the acquisition of three-dimensional images of the mouse heart is depicted in Figure 1a. A custom-made matrix ultrasound transducer array (Imasonic, SaS, Voray, France) consisting of 256 piezocomposite elements (4 MHz central frequency, ~100% detection bandwidth) densely distributed along its spherical active surface (40 mm radius, 90° angular coverage) was used to collect OA signals [19]. An optical parametric oscillator (OPO) laser (tunable wavelength between 670 and 1250 nm) guided with a fiber bundle inserted through a central aperture (8 mm diameter) of the transducer array was used to provide optical excitation. In the experiments, the pulse repetition frequency of the laser was set to 100 Hz and the wavelength was tuned to 700, 800 and 900 nm. OA signals were digitalized by a custom-made data acquisition system (Falkenstein Mikrosysteme GmbH, Taufkirchen, Germany), triggered with the Q-switch output of the laser, at 40 megasamples per second and transferred through a 1 Gbps Ethernet connection to a PC for storage and further processing. Raster scanning of the transducer array along the x-y axes (being y the vertical direction and x a normal direction parallel to the mouse surface) was performed inside a water tank by two translation motorized stages (Model: RCP2-RGD6C-I-56P-4, IAI Inc., Shizuoka Prefecture, Japan) in order to achieve a larger FOV (Figure 1a). For this, a custom-made holder for the spherical array along with a counter weight were attached to the stages [20,21]. The holder enables translating the array in a direction normal to the mouse surface. In the experiments, the array was fixed at a position for which the heart was located at approximately the center of its spherical surface. For each position of the transducer array, a sequence of 1000 consecutive frames was acquired to fully cover multiple cardiac cycles. Scanning parameters were adjusted with a custom MATLAB interface (MathWorks Inc., Natick, MA, USA).

2.2. Image Reconstruction, Stitching and Unmixing

Reconstruction of three-dimensional OA images was performed as follows. The acquired signals were first deconvolved with the electric impulse response of the individual piezoelectric elements and band-pass filtered between 0.1 and 7 MHz. Subsequently, a graphics processing unit (GPU)-implemented back-projection reconstruction algorithm was used to separate reconstruct the images corresponding to each excitation wavelength [22]. After reconstruction, the image sequences were synchronized based on phase analysis (Section 2.3) and further processed. Specifically, the images corresponding to different positions of the spherical array were stitched (compounded) by considering the voxels with maximum intensity in each frame [23]. Spectral un-mixing was performed to estimate the bio-distribution of two specific chromophores, namely the oxygenated (HbO₂) and deoxygenated (HbR) hemoglobin, from OA images at multiple wavelengths. Specifically, least square fitting of OA images at different wavelengths on a pixel-by-pixel basis to the molar extinction coefficients of HbO₂ and HbR was performed for this purpose [24].

2.3. Heartbeat Frequency and Phase Estimation

The heart and breathing rates were first estimated by calculating the mean Fourier transform of the time profiles of each voxel of the reconstructed three-dimensional images. The fundamental frequency of the heart could clearly be differentiated from the breathing frequency and its harmonics (Figure 1b). Three-dimensional phase maps were then calculated by considering the phase of the Fourier transform for the estimated calculated heart frequency. This enables monitoring the state of the heart during the cardiac cycle assuming that the heart beats regularly and in a cyclical manner.

2.4. Spatial Compounding and Coregistration

Coregistration and spatial compounding was performed based on three-dimensional phase mapping for the fundamental frequency of the heart. The phase map for a given position of the array and/or excitation wavelength was first calculated and taken as reference. Specifically, the image for the central position of the raster scan was taken as a reference for spatial compounding multiple positions and the image for 800 nm wavelength was taken as a reference for coregistration of multi-spectral data. The phase difference with respect to other images sequences was estimated via least square minimization with respect to the reference frame. For this, an area of 5 × 5 × 5 mm³ was considered in the case of multi-spectral data, while an area of 2 × 2 × 2 mm³ corresponding to the overlapping region was taken for different scanning positions. Spatial compounding of multiple scanning positions was then performed with a rigid body transformation consisting of a translation matrix. All data processing procedures were implemented in MATLAB (MathWorks Inc., Natick, MA, USA).

2.5. Animal Handling

In vivo mouse imaging experiments were performed in accordance with the Swiss Federal Act on Animal Protection and were approved by the Cantonal Veterinary Office Zürich (ZH 161/18). A female athymic nude mouse (8 weeks old, Janvier Lab, Le Genest-Saint-Isle, France) was imaged. A specifically designed holder was used to keep the mouse in a stationary position inside a water tank with temperature stabilized to 34 °C using a feedback-controlled heating stick. The mouse was anesthetized by isoflurane anesthesia (4% v/v for induction and 1.5% during the experiments, Abbott, Cham, Switzerland) in an oxygen/air mixture (100/400 mL/min). Gas anesthesia was provided through a custom-made breathing mask with the head positioned above the water surface at all times. The eyes were covered with vet ointment (Bepanthen, Bayer AG, Leverkusen, Germany) to ensure protection from the laser light as well as to prevent eye dehydration during the scanning. The mouse fully recovered after the experiments and no effect on its well-being and behavior was observed.

3. Results

OA data were acquired sequentially at several array positions and wavelengths; hence spatio-temporal discrepancies exist between the corresponding image sequences. Determination of the heart and breathing rates as well as phase information of the mouse heart is then essential for image compounding. Two prominent peaks in the low-frequency range of the temporal Fourier transform for a sequence of 1000 frames can clearly be distinguished (Figure 1b, top). These correspond to cardiac (7.1 Hz) and respiratory (2 Hz) rates, which are within the normal physiological range. The high-frame-rate of the OAT system (100 Hz) facilitates distinguishing all harmonic components corresponding to cardiac and respiratory motion, which would appear as alias frequencies if acquired at lower rates. The fundamental heart frequency can also be easily identified in the temporal Fourier transforms of a shorter sequence consisting of 100 frames (Figure 1b, bottom), which facilitates reducing the data required for effective compounding. Figure 1c shows three-dimensional views of three reconstructed OA images for 800 nm excitation at consecutive scanning positions in the vertical direction. It is shown that only a limited region of the heart can be effectively covered at each position. The corresponding phase maps for a given slice (indicated in Figure 1c) are displayed in Figure 1d. The difference in color distributions represents the phase shift between different sequences. Stitching was enabled by correcting the phase differences between consecutive scanning positions. This is illustrated in Figure 1d, where the central scanning position was taken as reference. Different heart structures are labelled in the combined phase map.

Figure 2 displays maximum intensity projections (MIPs) along the depth direction of OA images of the mouse heart for 800 nm wavelength excitation. The enhanced FOV achieved by stitching multiple images is clearly illustrated by comparing the image corresponding to a single position of the array (Figure 2a) and the compounded image (Figure 2b). Both images correspond to a FOV of 16 × 16 × 10 mm³. Only a part of the heart is visible in Figure 2a, and the image appears empty for regions outside the effective FOV of the array. Figure 2b shows that the whole heart can be visualized after the spatial compounding procedure. This image effectively covers a FOV of 16 × 16 × 10 mm³, which makes possible to clearly visualize the whole heart and surrounding vascular networks. The outline of the heart, the left and right ventricles (LV and RV), the left and right atria (LA and RA) as well as the coronary artery (CA) are labeled in Figure 2b. Figure 2c showcases a series of consecutive full-view images exhibiting heart motion during the cardiac cycle, from which it was possible to distinguish representative behavior of the heart at systole and diastole. The heart motion is more clearly visualized in a movie available in the on-line version of the journal (see the supplementary) corresponding to a rotating view of the compounded image.

The results obtained for the co-registration of multi-spectra OA images are shown in Figure 3. Figure 3a shows the MIPs along the depth direction of the OA images for a given position of the array as a function of excitation wavelength. Two representative time points corresponding to systole and diastole are shown. Physiological structures are labeled on the images. Generally, it is shown that the OA signal intensity increases with wavelength for structures containing highly oxygenated blood (e.g., the LA). The un-mixing results are displayed in Figure 3b. The images were normalized to the maximum signal intensity of the two components after un-mixing. As expected, it is shown that HbO₂ is the dominant component in all heart chambers, particularly in the LA. HbR was mainly detected in the RV as well as in vascular networks surrounding the heart.

4. Discussion and Conclusions

The feasibility to efficiently combine independent sequences of cardiac OA images can significantly advance the capabilities of this modality for cardiovascular research studies. In this work, we have shown two representative examples. By combining the images acquired by raster-scanning the transducer array, it was possible to identify cardiac structures such as the heart chambers or the coronary artery. Being able to visualize the entire heart is essential for quantifying key functional parameters. For example, the blood volume in the ventricles needs to be estimated to assess indicators of pumping effectiveness such as the end systolic volume (ESV), end diastolic volume (EDV), stroke volume (SV) and ejection fraction (EF). Oxygen saturation measurements in the heart can also be very valuable. For example, the oxygen saturation in the right atrium is an optimal estimator of mixed venous oxygen saturation (SvO2), which can be used to assess insufficiencies in cardiac output and oxygen delivery [25,26]. Myoglobin and hemoglobin oxygen saturation values in the cardiac muscle can also potentially be used as bio-markers of heart diseases [27].

Image coregistration was facilitated by the periodicity of breathing and heart motion. Recently, we have shown that a sequence of volumetric OA images of the heart could be efficiently compressed via principal component analysis (PCA), which enabled reconstructing subsequent frames with only a few OA signals [28]. In this work, we have shown that the fundamental frequencies and higher harmonics of respiratory and cardiac cycles could clearly be distinguished even when using relatively short acquisition sequences of 100 frames (1 s). Combination of independent images can then potentially be done in a short time, which may enable studying dynamic events such as hemodynamic responses to a breathing gas challenge [29,30]. It appears challenging to use the suggested approach to study faster events such as cardiac arrhythmias, which may represent a limitation of the methodology. However, frame rates in the order of kHz have recently been achieved in OAT [31,32]. This can potentially facilitate the development of more advanced image compounding approaches. Fourier analysis has further been used to provide a useful means to stitch multiple images acquired during a hand-held scan with no available reference on the position of the transducer array [23]. The results obtained in this work suggest that it may be further possible to stich cardiac images acquired in a hand-held manner. This is of particular relevance considering the ongoing clinical translation of the OA technology [33,34,35]. Imaging the human heart is challenged by the strong acoustic distortion induced in the rib cage. It is also important to take into account that the relatively large depth at which the human heart is located represents a limitation due to the strong attenuation of light. However, efficient image stitching may enable finding “acoustic and optical windows” facilitating estimating cardiac parameters that cannot be measured with existing technologies.

The enhanced information obtained by combining multiple sequences of images is of particular relevance in cardiovascular biology but can also be useful in other biomedical research fields where OA is currently being used. Generally, real-time imaging is affected by breathing motion artefacts. Arterial pulsation and structural motion induced by the heart beat can also affect the acquired images. These effects are particularly prominent for high-resolution OAT systems [18]. A gating imaging approach enabling compensating for breathing and heart motion can enhance the resolution of whole-body imaging systems [21,36,37] as well as avoid interference effects e.g., in neuroimaging studies [38,39]. In conclusion, Fourier analysis demonstrated that three-dimensional OA images of the heart contain sufficient information for the coregistration of independent frame sequences where synchronization is not possible. A self-gated imaging approach can then easily be defined, which is expected to provide new insights into cardiovascular biology and other fields.