Photoacoustics can image spreading depolarization deep in gyrencephalic brain

Spreading depolarization (SD) is a self-propagating wave of near-complete neuronal depolarization that is abundant in a wide range of neurological conditions, including stroke. SD was only recently documented in humans and is now considered a therapeutic target for brain injury, but the mechanisms related to SD in complex brains are not well understood. While there are numerous approaches to interventional imaging of SD on the exposed brain surface, measuring SD deep in brain is so far only possible with low spatiotemporal resolution and poor contrast. Here, we show that photoacoustic imaging enables the study of SD and its hemodynamics deep in the gyrencephalic brain with high spatiotemporal resolution. As rapid neuronal depolarization causes tissue hypoxia, we achieve this by continuously estimating blood oxygenation with an intraoperative hybrid photoacoustic and ultrasonic imaging system. Due to its high resolution, promising imaging depth and high contrast, this novel approach to SD imaging can yield new insights into SD and thereby lead to advances in stroke, and brain injury research.

www.nature.com/scientificreports www.nature.com/scientificreports/ temporal resolution in a seconds range for a sub-millimeter imaging field 12 and has a sub-millimeter penetration depth. TPM is therefore mostly used in small animal models. LS imaging or LS flowmetry images changes in cerebral blood flow in single vessels 22 . It is complementary to the larger field of view IOS imaging 21 which images reflectance changes of light in one 24 or two 23 narrow bands. IOS has a potentially sub-second temporal resolution 28 and micron spatial resolution, while again being diffusion limited to a sub-millimeter penetration depth and no depth resolution. NIRS, in contrast to the other optical techniques, is no imaging technique but employs point measurement probes 25 or optrode strips 26 to monitor millimeter scale areas similar to electrodes.
Functional magnetic resonance imaging (fMRI) with blood oxygen level dependent (BOLD) or diffusion weighted contrasts is the only modality that has been used to image the hemodynamic response of SD deep in brain 7 . Substantial drawbacks besides the complex imaging setup are the poor spatiotemporal resolution 29,30 and low contrast 7 when compared to optical or electrical measurements.
Overall, it can be concluded that the imaging methods proposed to date either feature high spatiotemporal resolution (IOS, TPM, LS) or are capable to provide depth-resolved information on SD beyond the surface (fMRI, implanted electrodes), but cannot provide both. To address this bottleneck, we investigate photoacoustic (PA) imaging as a possible high-resolution imaging technique for measuring SD deep in the gyrencephalic brain. Near infrared (NIR) light can penetrate deep into tissue, is scattered and gets diffused, thereby losing spatial information after a fraction of a millimeter. Photoacoustics 31 is capable of imaging beyond this sub-millimeter optical diffusion limit through the PA effect 32 ; light is delivered as a nanosecond laser pulse and where it is absorbed, it causes sudden thermoelastic expansion which in turn gives rise to acoustic waves. These acoustic waves emitted by the PA effect are in the ultrasound spectrum and therefore scatter much less than NIR light in tissue and can be detected by ultrasound (US) probes. Reconstructing their origin yields PA images. A multispectral stack of such images can be processed to reconstruct images of estimated tissue oxygenation that feature the spatiotemporal resolution and imaging depth of US combined with the optical contrast of NIRS. Multispectral photoacoustic imaging has shown to image blood oxygenation and perfusion in a variety of applications [33][34][35][36][37] . In the context of brain imaging, however, the application of PA has been restricted to lissencephalic brains 30 , which have been scanned with a range of PA imaging systems, including linear array probes 38 . The potential of PA imaging for monitoring SD remains to be investigated.
Rapid neuronal depolarization and repolarization causes tissue hypoxia 12 . Therefore, our work is based on the assumption that the imaging of hemodynamic changes with photoacoustics enables the monitoring of SD deep inside the tissue. We hypothesize that multispectral PA imaging is able to image SD induced hemodynamic changes in the entirety of the cortical gray matter of a gyrencephalic brain. For the purposes of this study, we measure an estimation of blood oxygenation (sO 2 ) and total hemoglobin (HbT). Our imaging concept, which is illustrated in Fig. 1, relies on a hybrid photoacoustic ultrasonic (PAUS) imaging system which combines (1) an US research system featuring a linear US transducer with a center frequency of 7.5 MHz and broad acoustic response 39 with (2) a near infrared (NIR) fast tuning optical parametric oscillator (OPO) laser 40 (see Methods). The system operates in an interleaved PAUS imaging mode, acquiring multispectral PA sequences with corresponding US images for each PA image. The multispectral PA image stream is the source for the hemodynamic contrast information. The registered US images provide a reference for stabilizing the PA image with a motion compensation approach 41 , this fixes the anatomical positions of ROIs during prolonged recording (see Methods). The concurrent US imaging is also used as anatomical reference for the physician (e. g. in needle guidance during the stimulation). Each multispectral PA image stack is converted into an image of estimated sO 2 and HbT. For this, the abundances of Hb and HbO in each pixel are estimated with a spectral unmixing method based on a non-negative constrained linear least squares solver. Our estimates are = + HbT Hb HbO (total hemoglobin) and = sO HbO 2 /HbT (blood oxygenation). www.nature.com/scientificreports www.nature.com/scientificreports/

experiments and Results
Two experiments were performed with our PAUS system to investigate whether the monitoring of tissue oxygenation with PA enables the detection and monitoring of SDs in the entire depth of the cortical gray matter of a gyrencephalic brain. In both experiments, brain activity was monitored with ECoG using a standard subdural electrode strip (Fig. 1).
The aim of the initial wave experiment was to investigate if the hemodynamic response of the brain to an induced SD can be imaged with multispectral PA. We performed the experiment in an uninjured brain. To analyze tissue hemodynamics before, during and after the occurence of SD, we took continuous PAUS measurements starting 24 min before the the first potassium chloride (KCl) stimulation (see Fig. 1) and ending one hour after the stimulation. After the experiment we cut sagittal surgical slices from the extracted brain to relate the acquired PA and US images to the brain morphology as seen on the exposed tissue. As shown in Fig. 2 we were able to image PA signal colorred of vasculature up to a depth of approximately 1 cm, which allowed us to image the entire cortical gray matter in the field of view of the imaging plane.
By estimating sO 2 in each pixel of our reconstructed multispectral images we observed a single wave of hypoxia spreading from the point of KCl stimulation through the tissue at a speed of approximately 5 mm/min. The estimated sO 2 for two sub-surface regions of interest (ROI) is plotted in Fig. 3a to illustrate this wave. sO 2 in a wide field of view during the same time frame is shown in Supplemental Video 1 played at a factor 100 timelapse. The wave of hypoxia coincides with the ECoG measurements on two electrodes in the proximity whose signals are plotted in Fig. 3b; they clearly show a single SD wave moving through the cortex, while the other electrodes on both hemispheres showed no change in activity. Figure 3c-g shows the change in estimated sO 2 in the region around the stimulation as hypoxia propagating through the tissue followed by an increase in sO 2 over the baseline.
The purpose of the cluster experiment was to investigate the hemodynamic changes during SD clusters with PAUS. To this end, we repeatedly stimulated the brain with KCl, before we started the PAUS measurement. We stopped these stimulations when we observed the occurrence of clustered SDs both electrically with ECoG and optically with IOS as an additional state of the art reference for surface measurement of SD. Once these SDs had www.nature.com/scientificreports www.nature.com/scientificreports/ subsided we positioned the PAUS probe and started the measurement with a single new KCl stimulation after a baseline recording period of 15 min (details see Methods).
As shown in Fig. 4a,b as well as Supplemental Video 2, there again was no change in the sO 2 estimation in the imaged sagittal plane in the left hemisphere during the baseline period. After KCl stimulation, we observed repetitive waves of hypoxia propagating through the imaging plane, followed by an overcompensation in sO 2 propagating through the cortex to up to a depth of approximately 2 to 5 mm below the brain surface. Figure 4c illustrates one such wave propagating from left to right during one minute as a change in sO 2 estimation. The speed of the waves was measured as 3-9 mm/min between ROI 1 and 2. ECoG measurements on the left hemisphere shown  www.nature.com/scientificreports www.nature.com/scientificreports/ in Fig. 5b indicate a SD cluster with the same frequency and speed of the sO 2 changes (Fig. 5a). As was the case in the initial wave experiment no change in ECoG activity in the right hemisphere was observed.
In addition to the sO 2 estimation from spectral unmixing we estimated the total hemoglobin (HbT) for the cluster experiment; this is visualized in Supplemental Video 3. The changes of HbT in the ROI are shown in Fig. 5c, where ROIs 3 and 4 seem to exhibit low frequency vascular fluctuations (LF-VF) 42 which appear to be depressed after SD 11,43 .

Discussion
We investigated the imaging of SDs based on the concept of PA imaging. Our approach involves simultaneous US and multispectral PA imaging for time-resolved reconstruction of tissue oxygenation in sagittal image slices. Two in vivo porcine experiments with our PAUS system provide the following evidence suggesting that our concept allows for the detection and monitoring of SD.
(1) Hypoxia consistent with ECoG: By estimation of sO 2 , we observed pronounced drops in estimated sO 2 after KCl stimulation (cf . Figs 3a and 4b). This local hypoxia lasted for around 30 seconds and was followed by an overcompensation or return to baseline sO 2 . These changes were consistent with the occurrence of SD in ECoG. The indicators we used to identify SD in ECoG were based on consensus 15 : A characteristic abrupt DC shift followed by a longer lasting positivity, and a reduction in amplitudes of spontaneous AC activity. Both of which needed to spread with a speed of 1.5-9.5 mm/min between electrodes and not cross hemispheres.
(2) Transient increase in blood volume: By estimation of HbT, we also observed a so-called normal hemodynamic response -a pronounced transient increase in blood volume (hyperemia) which was followed by a mild long-lasting oligemia (Fig. 5c). (3) Speed of wave propagation: Both the changes in sO 2 and HbT propagated through the gray matter at speeds of 3-9 mm/min. This is consistent with speed of SD reported in the literature as 1.7-9.2 mm/min 11 or 1.5-9.5 mm/min 15 in gyrencephalic brain (3-9 mm/min in porcine brain 13 ). (4) Low-frequency vascular fluctuations: We also observed changes in low-frequency vascular fluctuations (LF-VF) 43 (Fig. 5c). The observed LF-VF "display[ed] a spreading suppression in a similar fashion to that of SDs" in ECoG (see 26 ). Note that LF-VF were only visible in the vicinity of larger vessels, which was the rationale for placing ROIs 3 and 4 in such regions. www.nature.com/scientificreports www.nature.com/scientificreports/ We conclude from these observations that our measurements clearly support our initial hypothesis and suggest that PAUS is able to image SD as a change in sO 2 .
While such changes in estimated sO 2 can be detected with our approach, it is worth mentioning that the absolute estimated values are only local correlates to the real physiological value. This is due to location dependent (mostly depth dependent) fluence effects. The measured PA signal is proportional to the product of optical absorption and light fluence. Fluence is highly dependent on the distribution of chromophores and scatterers throughout the tissue. As a consequence, quantification of chromophore concentrations from measured signal is an ill-posed inverse problem and subject to ongoing research. As our methodology relies on the local analysis of relative changes, we do not consider this challenge a problem for this work. In general, PA measurements cannot be quantitatively compared in different locations. The signals are, however, sufficient to estimate local changes in sO 2 44 . The different amplitude of both sO 2 and HbT changes for various ROIs are therefore likely artifacts and better quantification methods are necessary to investigate this aspect further 37,45 .
A related issue is that the results from HbT estimation are more susceptible to absolute changes in light fluence over all measured wavelengths, compared to the sO 2 based measurements. Such a change in illumination geometry can occur due to SD related brain swelling 46 which causes a slow shift in the absolute HbT signal as shown in Fig. 4c. Lower fluence, generally in higher depths, also cause the signal-to-noise ratio to deteriorate, which can be observed when comparing ROI 1 with ROI 3 in Fig. 5a.
In contrast to all other methods proposed to monitor SD to date, our approach has the unique advantage that it features both high resolution and high imaging depth. While it is not suitable for imaging the entire gyrencephalic brain, penetration depth is sufficient to image the entire thickness of the cortical grey matter. The thickness of the perfused gray matter in a porcine brain is usually less than 5 mm. With blood absorption as a contrast, the far less perfused white matter will only exhibit signal in larger vessels. Such a vessel can be seen approximately at 1 cm depth in Fig. 2 and more clearly in a depth of about 1 cm in the Supplemental Videos 2 and 3 as well as Fig. 4a. We were not able to reliably image deeper structures, due to low light fluence.
Furthermore, the simultaneous PA and US imaging proved to be useful for anatomical orientation during the intervention (i. e. in needle guidance for KCl stimulation).
With these promising results, we see a potential use of PA imaging for SD characterization i. e. during pharmacological trials on the gyrencephalic brain. As the thickness of the human cerebral cortex is comparable, usually averaging 2.5 mm and not exceeding 5 mm 47 , PA imaging would be ideally suited for the study of SD in patients, as well. While, PA imaging cannot currently penetrate through an intact human skull 30 , PA imaging could for example be used postoperative to study SD in stroke patients 48 with a hemicraniectomy 49 .
Our pilot study strongly suggests that photoacoustics could become a valuable tool for detection, imaging, and monitoring SD. Due to its high spatiotemporal resolution this approach can be used to more precisely study where (i. e. which neuron layer) SDs originate and how they propagate, thus adding to our understanding of the nature of SD and its contribution to brain injuries and disease progression. Methods PAUS imaging system. The custom built hybrid PAUS imaging system is based on a 128 channel ultrasound data acquisition system (DiPhAs, Fraunhofer IBMT, St. Ingbert, Germany) with a 128-element linear US transducer operating on a center frequency of 7.5 MHz and broad acoustic response (L7-Xtech, Vermon, Tours, France). Due to its low level application programming interface (API) access, the system allows for raw data access and an interleaved PAUS imaging mode. This interleaved mode acquires US data from several shots after each PA data acquisition. The data acquisition (DAQ) module is combined with a fast tuning OPO laser cart (Phocus Mobile, Opotek, Carlsbad, USA) which yields 690 nm-950 nm, 5 ns long laser pulses with a pulse repetition rate of 20 Hz and a per laser pulse power of up to 50 mJ. The wavelength of each laser pulse can be tuned in between shots, allowing for real time multispectral acquisition sequences. Laser fiber bundles ending in two line arrays are attached to the transducer by a 3D-printed frame including acrylic windows for the laser output. For each experiment the entire probe was wrapped in a sterile ultrasound probe cover and gold leaf of sub-micron thickness was placed between the active area of the US transducer and the probe cover to reduce artifacts created by light absorption in the US transducer. Gold leaf was chosen for its high near infrared reflectance. For live imaging and recording all APIs to the system were integrated in the Medical Imaging Interaction Toolkit (MITK) software framework and the MITK workbench application was used throughout the intervention to control the PAUS system, configure the image acquisition, and show live PA and US imaging streams. During our experiments we visualized both streams with 15-20 fps using delay and sum (DAS) beamforming for an imaging depth of 4 cm with 256 reconstructed lines. For the initial wave experiment we imaged the wavelength sequencence (735 nm, 756 nm, 850 nm, 900 nm) selected to distinguish Hb and HbO 50 . Because we added an estimation of HbT in the cluster experiment we instead imaged the isosbestic point of Hb and HbO at 798 nm for further reference, leading to the wavelength sequencence (760 nm, 798 nm, 858 nm). Image reconstruction. The raw radiofrequency (rf) PA data acquired during the experiments was matched with the laser pulse energies recorded by a pyroelectric sensor (Ophir PE25-C, Ophir Optronics, North Andover, USA) built in the laser system (Phocus Mobile, Opotek, Carlsbad, USA) and matched with the wavelengths of the laser pulses measured by a spectrometer (HR2000+, Ocean Optics, Dunedin, USA). The wavelengths of the pulses were measured independently of the imaging system to account for calibration errors. The rf PA slice was then corrected for the corresponding pulse energies. The recorded PAUS data was already beamformed live during the experiment to reduce the system load writing to disk. A single US slice was recorded after each PA slice. The US image was a compounded image averaged from US data acquired at five angles, equidistant from +10 deg to −10 deg and beamformed to 256 lines using a DAS algorithm with boxcar apodization. To convert the acquired rf PA slices into meaningful images suitable for multispectral analysis, the slices were beamformed with www.nature.com/scientificreports www.nature.com/scientificreports/ a reference DAS implementation 51 using Hanning apodization to 512 lines. B-Mode images with isotropic pixel spacing of 0.075 mm were formed with a Hilbert transform based envelope detection filter. US B-Mode images were formed in the same way, only adding a subsequent logarithmic compression.
Motion compensation. The PA images obtained after beamforming are corrected for inter-frame motion introduced by breathing, pulse or swelling. The later is especially relevant as spreading depolarization is closely linked to brain swelling 46 . The data is stabilized by motion correction (1) in order to enable a more stable spectral unmixing and (2) to assure that a given pixel location corresponds to the same anatomical ROI. To correct for the inter-frame motion an optical flow based method is employed. The optical flow of each US image relative to the first US image in the entire recording is estimated using an algorithm proposed by Farnebäck 52 . The flow estimated from the US B-mode image is then used to warp the corresponding PA B-mode image.
Experimental data analysis. Because of the slow propagation of SD wave fronts we averaged over ten motion corrected frames of the same wavelength and still have the 1 s temporal resolution of IOS. Spectral unmixing of those image sequences was then performed using a python implementation of a non-negative constrained linear least squares solver (scipy.optimize.nnls). In all figures and supplemental material plots and videos one PA data point is averaged over ten frames and then averaged over the ROIs. Speeds of SD wavefronts were obtained by measuring the time between the local minima of sO 2 ROI 1 and ROI 2. The positions of ROIs 1 and 2 in both experiments were chosen at 1 mm depth and 7.5 mm apart, in the center of the reconstructed image stream. ROIs 3 and 4 were chosen deeper and close to larger vessels to investigate the LF-VF effect which, as discussed, can only be observed there. ROIs are otherwise representative of the entire data set as can be seen in the Supplemental Videos. To take only perfused tissue into account, sO 2 estimates were masked to only include pixels where the unmixed HbT value was larger than the mean background noise plus two standard deviations.

Animals.
Protocols for all experiments were approved by the institutional animal care and use committee in Karlsruhe, Baden-Wuerttemberg, Germany (Protocol No. 35-9185.81/G-174/16). All experiments were performed in accordance with the relevant guidelines and regulations. Female German Landrace swines of 31 and 33 kg were premedicated with Midazolam (Dormicum 0.5-0.7 mg/kg) and Azaperone (Stresnil 4 mg/kg) intramuscularly. After premedication, two venous lines were placed in the ear veins, and propofol (Disoprivan 5-7 mg/kg) was administered intravenously to facilitate the intubation. The animals were then intubated and mechanically ventilated and the pressure controlled ventilation was adapted to a respiration rate of 12-20/min, a flow of 2.5 l O 2 / min, 2.0 l air/min, FiO2 35% and volume 7-10 ml per kg. The maintenance of anesthesia required inhalational anesthesia with isoflurane (Isosthesia 0.6-1.0%) and intravenous midazolam at a continuous dose of 0.5-0.7 mg/ kg/h via perfusion and maintained throughout the entire experiment. If a wakening reaction occurred, a bolus of propofol (Disoprivan 5-7 mg/kg) was administered. Temperature was monitored with a rectal probe. A 4-Fr catheter was placed in the right femoral artery for permanent monitoring of the mean arterial blood pressure (Raumedic AG, Helmbrechts, Germany). Capillary oxygen saturation was monitored from one ear. Arterial blood gases were obtained in the animal used for the initial wave experiment. Ringer's solution was given intravenously over 8-12 h, to compensate for intraoperative bleeding, urinary output and insensible losses. The two animals used in this study were used primarily for this project. After finishing the protocol the animal used for the cluster experiment was used for other unrelated studies.
Surgery. Animals were fixed in a stereotactic frame (Standard Stereotaxic Instruments, RWD Life Science, Shenzhen, China) and an extensive craniotomy with excision of the dura mater was performed, to view the subarachnoidal space bilaterally. Initially, the brain surface was immersed for 30 to 40 min in a standard lactated Ringer's solution with an elevated K+ concentration (7 mmol/l), as preconditioning for SD induction, as proposed by Bowyer et al. 53 for the KCl model of SD. EcoG was performed with two strips of 5 electrodes each (Ad-tech, Racine, Wisconsin, USA) that were placed at the lateral margins of the craniotomy below the dura mater and above the parietal cortex. A camera for IOS imaging and its corresponding light sources were mounted above the stereotactic frame. After preconditioning a 5-10 mm deep paraffin pool was filled over the exposed cortex, to reduce the diffusion of the KCl stimulation. When necessary, paraffin was withdrawn and new paraffin was added. The preparation time was 4-5 h before the KCl stimulations started. A gel pad (Aquaflex Ultrasound Gel Pad, Parker Laboratories, Fairfield, USA) was cut in shape of the exposed brain surface and placed in the paraffin pool. The custom designed PAUS probe (see Methods, PAUS Imaging system) was placed on top the gel pad and fixed relative to the frame. For the initial wave experiment the gel pad and PAUS imaging system was placed before the initial stimulation. For the cluster experiment the gel pad and system was placed after the initial KCl stimulations and the accompanying IOS imaging was performed. With the help of live US imaging, it was positioned to image a sagittal plane of the left hemisphere approximately 1 cm from the midline. For the cluster experiment we waited until any residual SD from prior stimulation subsided in the ECoG monitoring. Only then did we start recording PAUS data in a sagittal plane for 15 min as a baseline. After a sufficient baseline recording, spreading depolarization was triggered using 2-5 μl of 1 mol/l KCl solution with a Hamilton syringe. The stimulus needle was guided using the live PAUS image streams visualized in MITK.
Electrocorticography. ECoG recording with the subdural electrodes was perfomed in 10 active channels, using the Powerlab 16/SP analogue/digital converter coupled with the LabChart-7 software (ADInstruments,