Adult C57BL/6 mice (8 weeks old, 25–35 g), provided by Animal Center of Southern Medical University (Guangzhou, China), were used in the present study. The animals were housed for at least 1 week before surgery under controlled environmental conditions with ambient temperature of 25°C, relative humidity of 65%, and 12/12 h light-dark cycle. The animals were free access to gain food and water. The protocol was reviewed and approved by the Institution Animal Care and Use Committee of The Second Affiliated Hospital of Guangzhou Medical University (NO. B2019-037).

The mice were anesthetized using isoflurane (1–2%/oxygen). A cortical stroke model was made by occlusion of the distal middle cerebral artery using a previously reported method (Wen et al., 2019). In brief, a 1 cm skin incision was made between the ear and eye (usually right side). The temporal muscle was removed from the skull. A hole was made by the drill right above the MCA and the artery was coagulated by the electrocoagulation forceps proximal and distal to the bifurcation. The temporal muscle was relocated to cover the burr hole. Suture the wound and place the animal in a nursing box at 32°C to recover from anesthesia and return it to the cage.

Magnetic resonance imaging (MRI) scanning was performed 6 h after stroke onset using a Bruker Biospec 7.0 T system (PharmaScan, Bruker Biospin, Rheinstetten, Germany) with a mouse brain array coil and a transmit only volume coil. The anesthetized animals were secured within the cradle by tooth and ear bars, and a mouse head four-channel phased array surface receiver coil was placed on the head. Body temperature was maintained at 37 ± 0.5°C during the MRI scanning procedure by a closed circuit thermal jacket.

T2 weighted scans using a fast-spin echo sequence: echo time (TE) 33 ms, repetition time (TR) 8,000 ms, field of view (FOV) 30 mm × 30 mm, acquisition matrix 512 × 512, acquiring 0.4 mm thick slices. A four-shot spin-echo planar imaging DWI scan (TE 30.5 ms, TR 8,000 ms, FOV 25 mm × 25 mm, acquisition matrix 96 × 96, three directions x, y, z, B-values = 0 and 1,000 s/mm², 50 contiguous thick, 0.4 mm thick). The PWI can be performed non-invasively by tagging protons in the arterial blood supply with an inversion pulse by using a T2*-weighted echo-planar sequence: TE 76 ms, TR 11,482 ms, FOV 24 mm × 24 mm, acquisition matrix 128 × 128, 1 excitation, repetition 50 times. All MR images were processed on a commercial workstation (ParaVision Acquisition 6.0.1).

As shown in Figure 1A, the spatial maps were constructed based on the spot barcodes derived from the ST data. The gene expression quantities were then connected to the corresponding spots. Global Moran’s I was calculated to determine whether it was clustered. If not, the gene was aborted for further analysis. Next, the ROI was circled in the ST maps and the spot barcodes of the ROI were obtained. Then, local Moran’s I was used to determine the CO type, which was group into five categories: non-significant, H-H, L-L, high-to-low (H-L), and low-to-high (L-H) spots (Figure 1B). The number of H-H and/or L-L spots were used for the following statistics. Finally, visualization of DEGs was performed using Hot Spot analysis with Arcmap 10.8.

Global Moran’s I measures the spatial autocorrelation based on both the locations and the quantities of features. A value of 0 for Moran’s I indicated no autocorrelation. In the brain section, no gene showed no autocorrelation. Z-score and P-value were used to evaluate the significance of global Moran’s I (Figure 2). Based on the global Moran’s I, z-score and P-value, the spatial pattern were group into: clustered, random and dispersed (Figure 2). If global Moran’s I was close to + 1 and z-score was more than 1.65, the pattern was clustered, which means elevated gene expression had similar elevation values close to each other like Nrgn (Figure 2). If global Moran’s I was close to −1 and z-score was less than −1.65, the pattern was dispersed, which means dissimilar values were next to each other like a checkerboard (Figure 2).

Spatial transcriptomic data of normal brain section was obtained from the 10X Genomics dataset. The brain hemisphere was divided into sensorimotor cortex, basal ganglia, cingulate gyrus, entorhinal cortex, hippocampus, thalamus, and hypothalamus (Figure 3A). The ratio of H-H spots in the sensorimotor cortex was calculated and is shown in Figure 3E. Then we determined the proportion of H-H spots (more than 50% in Figure 3E) in the other brain regions (Figure 3F). For instance, Arpp19 mostly enriched in the sensorimotor cortex while it also was found in the cingulate gyrus and entorhinal cortex, which was verified in ST map (Figures 3B–D,F). DEGs in other brain regions were also identified and typical genes are shown in Figure 3G. All these results indicated that saSpatial was able to detect DEGs of ROI within a section.

Spatial analysis of spatial transcriptomics was then applied to determine the DEGs in the brain section from ischemic stroke. Ischemic stroke was induced in a mouse model by occlusion of the distal middle cerebral artery. MRI scanning was performed 6 h after stroke onset to determine the ischemic penumbra (mismatch of perfusion- and diffusion-weighted imaging) and core. Penumbra has been defined as brain tissue at a risk of infarction (Albers et al., 2018). The brains were then harvested for ST sequencing. Ischemic core and penumbra and the corresponding area in the normal brain section were marked with reference to HE and MRI images (Figure 4B).

Considering numerous gene expressions were similar in the ischemic penumbra and core, which was also proved by the following data (Figure 4A), there were three comparisons of the number of H-H spots: between the ischemic penumbra or core and the corresponding area in the normal brain sections; between the ischemic penumbra or core and the normal region in the ischemic stroke brain section; and between ischemic core and penumbra (Figure 4A). Bubble plots of the top 20 DEGs for each comparison identified by saSpatial are shown in Figure 4A. For instance, Gm42418 was highly expressed in the ischemic hemisphere while there was no significant difference between ischemic ore and penumbra. In contrast, Cdk5r2 was significantly increased in the ischemic core but not in the penumbra (Figure 4C). Hence, Cdk5r2 but not Gm42418 was a marker gene for ischemic core. As for ischemic penumbra, saSpatial identified the up-regulated DEGs like Hspa1a, which was like a thin loop around the core (Figure 4C).

Substantial attention was paid to the up-regulated expression, while expression of numerous genes were halted by ischemia. In contrast to the up-regulated DEGs as described above, the number of L-L spots was used for detection of down-regulated DEGs. Bubble plots of the top 20 DEGs for each comparison are shown in Figure 5A. There were a lot of genes down-regulated in ischemic area while most of them showed no significant differences between ischemic core and penumbra such as Vps8 and Sparc (Figure 4A). For the ischemic core, Lrrc77 was a marker gene (Figure 5B). As for the ischemic penumbra, down-regulated DEGs like Prkcb were identified (Figure 5C).

Differentially expressed genes detected by saSpatial were compared with that detected by FindMarkers (Figure 6). For the ischemic core, FindMarkers found only 210 DEGs while there were 1,021 DEGs identified by saSpatial (Figure 6A and Supplementary Table 1). The majority of DEGs detected by FindMarkers could be also found by saSpatial (Figure 6A) while numerous genes like Itm2a, Foxp1, Cbln2, Mtus2, Nab2, and Cnp were not detected by FindMarkers (Figure 6B). For the ischemic penumbra, FindMarkers can only find 44 DEGs and the majority of DEGs were missed (Figure 6C). saSpatial found that there were 497 potentially valuable DEGs in the ischemic penumbra (Figure 6C) such as Gng13, Cyr61, Fosb, Gadd45g, S100b, and Gadd45b, which were confirmed by ST maps (Figure 6D).

The sensorimotor cortex is composed of six layers. DEGs in each layer were revealed by saSpatial and the top 3 DEGs in each layer are shown in Figure 7A. Layer V is the internal pyramidal layer and contains large pyramidal neurons, the axons of which form the corticospinal tract. saSpatial showed that Ighm was a marker gene in Layer V (Figure 7B).

Different cortical layers responded differently to ischemic insult. Six layers of ischemic core and that in the normal brain section were marked with reference to the HE images (Figure 8A). The fraction of DEGs in each layer is shown in Figure 8B. The ST and Hotspot maps are shown in Figure 8C. For example, Cbln4, which was not found in the normal brain, was up-regulated in the ischemic core and enriched mainly in the Layer III. Scn4b was up-regulated in the Layer V while it was down-regulated in other layers (Figures 8B,C).

Spatial analysis of spatial transcriptomics was based on the ST data. Hence, the resolution of saSpatial relied on the raw ST data. A spot in ST data represented genes of an area including about 10 cells. Thus, saSpatial cannot reveal the cellular types. Secondly, CO types of spots mainly consider the effect of adjacent spots while some genes were affected by remote ones. Finally, the DEGs identified by saSpatial needs to be verified by the additional biological experiment.