Three-dimensional reconstruction method suitable for leguminous root nodule organ single cell space transcriptome

Abstract

The invention relates to a method suitable for three-dimensional reconstruction of a leguminous root nodule organ single cell space transcriptome, which comprises the steps of integrating single cell data and space group data and then integrating X-ray imaging data, optionally, the single cell data is obtained through an snRNA-seq technology, and optionally, the space group data is obtained through a Stereo-seq technology. The method is used for performing three-dimensional gene expression map construction on the plant root nodules, the advantages of the three technologies are combined together, and a user is allowed to explore root nodule cell population composition and a gene expression mode in a three-dimensional mode.

Description

TECHNICAL FIELD

[0001] The present application belongs to the field of biotechnology and imaging technology, and specifically relates to a method suitable for three-dimensional reconstruction of single-cell spatial transcriptome of plant root nodule, for example, legume plant, which is a novel method combining X-ray microscopy imaging technology with single-cell data and spatial transcriptomics, and is used for constructing three-dimensional gene expression map of plant root nodule under non-destructive conditions. BACKGROUND

[0002] X-ray CT scanning is used as a powerful tool for obtaining complex three-dimensional organ cell structure data. Since the attenuation of X-rays depends on density, it is particularly suitable for plant imaging due to the ubiquitous intercellular space in plant organs. X-ray CT scanning generates a three-dimensional image of a plant organ or tissue from a stack of X-ray images obtained when rotating the sample. Its advantages include: large penetration depth depending on photon energy, scalable field of view and resolution from 500 nm to 1 mm, minimal sample preparation and non-invasive real-time non-destructive imaging.

[0003] However, existing X-ray CT scanning imaging technology is mostly used for anatomical structure analysis, and lacks the ability to integrate with molecular expression information. Traditional spatial transcriptome technology relies on frozen sections or laser capture microdissection, which destroys the three-dimensional structure of the tissue, and the spatial information is limited to two-dimensional level.

[0004] Therefore, there is an urgent need for a non-invasive, high-resolution three-dimensional imaging method to integrate the fine anatomical structure of plant organs with spatially resolved gene expression information to study the spatiotemporal regulation mechanism of plant organ development process. SUMMARY

[0005] In view of the above problems, the present application integrates single-cell data and spatial group data, and further integrates the overall organ cell characteristics of X-rays to construct a 3D transcriptome map of root nodule single-cell spatial group, which combines the advantages of the three technologies together, allowing users to explore the composition of root nodule cell groups and gene expression patterns in three dimensions.

[0006] Specifically, the present application provides the following technical solutions.

[0007] In one aspect, the present application provides a method suitable for three-dimensional reconstruction of single-cell spatial transcriptome of legume plant, for example, soybean root nodule organ, which comprises the step of integrating single-cell data and spatial group data, and then integrating X-ray imaging data.

[0008] In some embodiments, the single-cell data is obtained by snRNA-seq technology.

[0009] In some embodiments, the spatial group data is obtained by Stereo-seq technology.

[0010] In some embodiments, the obtaining of the X-ray imaging data comprises a step of drying the sample.

[0011] In some embodiments, the snRNA-seq technique comprises the following steps:

[0012] a1. Single nuclei suspension, emulsion breaking, mRNA capture bead collection, reverse transcription, cDNA amplification and purification by titration microfluidic encapsulation;

[0013] a2. The resulting cDNA product is fragmented into 250-400 bp fragments;

[0014] a3. Indexing sequencing library construction for sequencing in paired-end mode, where each sequencing read comprises a 30 bp Read 1, a 100 bp Read 2 gene sequence and a 10 bp sample index barcode read;

[0015] a4. Preprocessing of snRNA-seq data.

[0016] In some embodiments, the 30 bp Read 1 comprises 10 bp of cell barcode 1, 10 bp of cell barcode 2 and 10 bp of unique molecular identifier.

[0017] In some embodiments, preprocessing of snRNA-seq data includes filtering raw sequencing reads and demultiplexing by barcode assignment; aligning the obtained reads to a soybean reference genome, e.g., Gmax_ZH13_V2.0, using Spliced Transcripts Alignment to a Reference (STAR), annotating to gene sets with PISA; automatically identifying valid cells using the “barcodeRanks” function of the DropletUtils package, removing background beads based on the number distribution of unique molecular identifiers per cell and reads with less than 500 UMIs; then calculating gene expression for cells using PISA and creating a gene x cell matrix for each library; merging expression data from different libraries of the same tissue using the “Merge” function of the R package Seurat, followed by normalization using the “NormalizeData” function with default parameters; after log-normalization by the “ScaleData” function, selecting 1500 highly variable genes using the “FindVariableFeatures” function; for subsequent clustering and visualization, using a principal component analysis (PCA)-based dimensionality reduction algorithm to extract multiple, e.g., 40, PCs, using the “FindNeighbors” function for Lovain clustering with a resolution of 1; projecting the clustering results to two-dimensional space using the “RunUMAP” function; identifying differentially expressed genes using the “FindAllMarkers” function.

[0018] In some embodiments, filtering raw sequencing reads includes removing reads with average base quality score less than 4, more than 2 bases with quality score less than 10, containing N bases, or inappropriate barcodes.

[0019] In some embodiments, the Stereo-seq technique comprises the following steps:

[0020] b1. Sample preparation and embedding;

[0021] b2. Cryosection preparation;

[0022] b3. RNA quality detection of cryosections;

[0023] b4. Chip processing;

[0024] b5. Tissue fixation;

[0025] b6. Tissue permeabilization treatment;

[0026] b7. In situ reverse transcription reaction;

[0027] b8. Tissue removal and product recovery;

[0028] b9. Reverse transcription product recovery and amplification;

[0029] b10. Reverse transcription product purification;

[0030] b11. Reverse transcription product fragmentation and amplification;

[0031] b12. Dual-fragment screening;

[0032] b13. Sequencing;

[0033] b14. Stereo-seq data preprocessing.

[0034] In some embodiments, in the RNA quality detection step of the frozen section, RIN≥7 is considered as a qualified sample.

[0035] In some embodiments, in the tissue permeation treatment step, the permeation time is 8 minutes.

[0036] In some embodiments, in the dual-fragment screening step, the quality control standard is that the yield is greater than 300 ng, and the main peak of the fragment distribution is between 400-600.

[0037] In some embodiments, in the sequencing step, the library concentration requires Qubit quantification≥2nM, and sequencing is performed in paired-end mode, wherein each sequencing read contains a 50 bp Read 1 and a 100 bp Read 2, and the spatial coordinate information is matched by barcode.

[0038] In some embodiments, Stereo-seq data preprocessing includes determining coordinate identification by mapping Stereo-seq raw data to the coordinates of the in situ capture chip, removing reads with UMI quality scores less than 10, aligning the remaining reads to a reference genome such as Gmax_ZH13_V2.0, and annotating to a gene set, and finally generating an expression profile matrix containing coordinate identification information for subsequent analysis.

[0039] In some embodiments, since the plant tissue only covers part of the surface of each Stereo-seq chip, a custom script is used to extract the exact area of the spatial transcriptome data, align the spatial transcriptome data with the ssDNA image, then determine the boundary according to the sample edge, and then convolve the spatial expression profile matrix into an N×N size spatial expression profile block and name it binN, N=50, and further adopt a ssDNA-based cell segmentation technique to convolve the spatial expression profile matrix into a cell-level spatial expression profile block.

[0040] In some embodiments, the ssDNA-based cell segmentation technique convolves the spatial transcriptome matrix into cell-level spatial transcriptome patches includes using the Patchify package in Python to segment the ssDNA image into sub-images of the same size; using the cyto2 algorithm in the Cellpose package to identify cells in each sub-image; using the pathchify package to merge the cell segmentation results of each sub-image, and output the merged cell segmentation results.

[0041] In some embodiments, when mapping the Stereo-seq raw data to the coordinates of the in situ capture chip, one mismatch is allowed to tolerate sequencing or PCR errors.

[0042] In some embodiments, integrating the single-cell data and spatial group data, i.e., integrating the pre-processed snRNA-seq data and Stereo-seq data.

[0043] In some embodiments, integrating the pre-processed snRNA-seq data and Stereo-seq data includes using the “FindTransferAnchors” and “TransferData” functions of the R package Seurat to calculate the proportion of cell clusters determined by snRNA-seq that are located to each cluster determined by spatial transcriptomics, completing the UMAP projection from snRNA-seq to cluster embedding of Stereo-seq by the “Knn” function, obtaining the spatial position of the clusters determined by snRNA-seq, and annotating them according to their anatomical structure, and using the R package Clustertree to identify other cell clusters based on the lineage relationship of the previously annotated cell clusters.

[0044] In some embodiments, the step of re-integrating the X-ray imaging data comprises extracting tissue images in the circular field of view of the serial scan pictures using the Segment Editor module in 3D Slicer; based on the extracted tissue images of the serial sections, constructing a 3D Mesh of the outer contour of the tissue using 3D Slicer and smoothing it, and storing the final 3D Mesh as an obj format; based on the first 40 principal components, clustering the single-cell transcriptome data using the Louvain algorithm (resolution = 1) to identify the cell type features of the tissue; mapping the Stereo-seq data onto the single-cell data using the FindTransferAnchors and TransferData functions in Seurat, and determining the cell type of each Squarebin of the Stereo-seq data according to its mapping score; calculating the Z-axis position of the corresponding 3D Mesh of the Stereo-seq according to the section position in the Stereo-seq serial section experiment record; using the TrakEM2 elastic registration algorithm in ImageJ to register the Stereo-seq serial sections with the X-ray images at the corresponding Z-axis positions, thereby determining the spatial position of the Stereo-seq serial sections in the 3D Mesh; and finally displaying the information such as gene expression and tissue type information on the Stereo-seq serial sections at the accurate spatial positions in the tissue 3D Mesh, thereby obtaining a three-dimensional reconstruction atlas of the nodule organ spatial transcriptome.

[0045] In some embodiments, the step of drying the sample comprises the following:

[0046] c1. 50% formaldehyde-acetic acid-ethanol (FAA) vacuum infiltration (-0.08 MPa) until no air bubbles come out of the sample, place the sample on ice, fix the sample overnight at room temperature, perform gradient ethanol dehydration in a 15 ml centrifuge tube, the gradient is 50%→70%→80%→90%→95%→100%→100%→100% ethanol, 15 min for each gradient, and place on ice;

[0047] c2. Pour out 100% ethanol and replace with 3:1 volume ratio of ethanol: acetone for 30 min, and place on ice;

[0048] c3. Replace with 3:2 volume ratio of ethanol: acetone for 30 min, and place on ice;

[0049] c4. Replace with acetone for 30 min;

[0050] c5. Repeat acetone replacement for 30 min, place the sample in a sample basket into a critical point drying instrument, replace with liquid CO2, slowly exhaust after warming up to 31.1°C, and complete critical point drying.

[0051] In another aspect, the present application provides a sample preparation method suitable for plant nodule X-ray imaging, comprising the following steps:

[0052] c1. 50% formaldehyde-acetic acid-ethanol (FAA) vacuum infiltration (-0.08 MPa) until no bubbles come out of the sample, place the sample on ice, fix the sample overnight at room temperature, perform gradient ethanol dehydration in a 15 ml centrifuge tube, the gradient is 50%→70%→80%→90%→95%→100%→100%→100% ethanol, 15 min for each gradient, and place on ice;

[0053] c2. Replace 100% ethanol with 3:1 volume ratio of ethanol:acetone for 30 min, and place on ice;

[0054] c3. Replace with 3:2 volume ratio of ethanol:acetone for 30 min, and place on ice;

[0055] c4. Replace with acetone for 30 min;

[0056] c5. Repeat acetone replacement for 30 min, place the sample in a sample basket and put it into a critical point drying instrument, replace with liquid CO2, slowly exhaust after warming up to 31.1℃, and complete critical point drying.

[0057] Definitions

[0058] X-ray imaging technology: an advanced characterization technology based on the penetration of X-rays for three-dimensional imaging, a non-destructive 3D imaging technology that can reveal the complex structure inside the sample.

[0059] RIN: full name for RNA integrity number, an important indicator for judging RNA integrity. The highest RIN value is 10, representing very good RNA integrity, and the lowest value is 0, representing very poor RNA integrity.

[0060] Neuroglancer: a high-performance WebGL viewer and visualization framework developed by Google's Connectomics team, used to view three-dimensional volume data. It supports multiple data sources and can display arbitrary (non-axis-aligned) volume data slice views as well as 3D networks and line segment models (skeletons).

[0061] HDF5: a file format that supports hierarchical data storage, allowing users to store multiple datasets and metadata in the same file.

[0062] 3D Slicer: an open-source software for image analysis, visualization, and processing. It supports multiple imaging modalities and provides powerful three-dimensional reconstruction, segmentation, registration, and quantitative analysis functions.

[0063] Registration: refers to the matching of geographical coordinates of different image patterns obtained by different imaging methods in the same area, including geometric correction, projection transformation and uniform scale processing.

[0064] TrakEM2 software: a program for morphological data mining, three-dimensional modeling and image stitching, registration, editing and annotation.

[0065] Coordinate mapping: refers to the process of converting a point or a set of points from one coordinate system to another in two-dimensional or three-dimensional space.

[0066] snRNA-seq: a technique for determining RNA molecules in the nucleus. It obtains the transcriptome information of the nucleus by separating the nucleus from the cell, then performing RNA extraction, transcript amplification and sequencing. During data acquisition, strict control and screening of data quality are required.

[0067] Beneficial effects

[0068] The method of the present application has the following advantages:

[0069] 1. Non-destructive: X-ray imaging avoids physical sectioning, preserving sample integrity for subsequent dynamic observation, that is, X-ray imaging technology can realize dynamic observation of cells in intact samples.

[0070] 2. Whole organ integration: supports data correlation from organ whole to single cell level.

[0071] 3. Improved sample preparation and X-ray imaging technology provide a powerful method for imaging of cell tissue in soybean nodule development.

[0072] 4. Although Stereo-seq enables sequencing-based spatial transcriptomics to measure gene expression levels in the whole transcriptome at subcellular resolution, it is still difficult to reconstruct organ 3D single cell expression profiles in current plant research. The combination of X-ray imaging technology and Stereo-seq in the present application can obtain 3D single cell expression profiles, which is beneficial to display the spatial group gene expression profile from a new perspective.

[0073] 5. Spatial group data is less but has spatial positioning, single cell data is larger, cell clustering is more clear, but cell loses spatial positioning information, and it is difficult to annotate cell groups. Integrating the advantages of the two, and adding the whole organ cytology characteristics of X-ray, constructing the 3D transcriptome atlas of nodule single cell spatial group, can combine the advantages of the three technologies together, allowing users to explore the composition of nodule cell groups and gene expression patterns in three dimensions. BRIEF DESCRIPTION OF DRAWINGS

[0074] Figure 1 The X-ray imaging results of soybean nodule. A: 3D view of the whole nodule. B: Cell morphology display of longitudinal section of nodule. C: Cell morphology display of transverse section of nodule.

[0075] Figure 2 The results of permeation time test of spatial group chip of nodule section are shown. The position of nucleic acid is in situ in the tissue after permeation for 8 minutes, and obvious nucleic acid diffusion is observed after permeation for 10 minutes, so that the permeation time of 8 minutes is the best permeation time of the nodule in this period, which is used for subsequent formal experiment conditions.

[0076] Figure 3 The structural information of spatial group section of selected nodule at each level is shown.

[0077] Figure 4 The results of spatial group data quality evaluation of nodule are shown.

[0078] Figure 5 The spatial bin division based on Stereo-seq data is shown, and 12 cell groups collectively show spatial positioning. The resolution is set to bin 50.

[0079] Figure 6 The spatial bin division based on Stereo-seq data is shown, and a single cell group shows spatial positioning. The resolution is set to bin 50.

[0080] Figure 7 The cell clustering and annotation based on snRNA-seq data are shown. One dot represents one nucleus, a total of 10113 nuclei, and a total of 51921 soybean coding genes are detected, of which there are an average of 1982 mRNAs per nucleus, and the number of genes is 1183.

[0081] Figure 8 The specific marker genes of each cell group of snRNA-seq data are shown.

[0082] Figure 9 The cell segmentation based on spatial group ssDNA image is shown, and the data of nodule single cell spatial group is displayed at the cell level resolution.

[0083] Figure 10 The spatial positioning of a single cell group of snRNA-seq data is shown.

[0084] Figure 11 The combined spatial positioning of the total cell group of snRNA-seq data is shown.

[0085] Figure 12Registration and integration based on soybean nodule image and spatial group data are shown, and spatial group data of root tip can be inquired from cell clustering, gene expression and cell localization. DETAILED DESCRIPTION

[0086] In order to make the objects, technical solutions and advantages of the present application clearer, the present application is further described in detail below with reference to specific embodiments and drawings.

[0087] Example 1: Construction of three-dimensional single-cell spatial transcriptome of soybean nodule

[0088] I. X-ray imaging

[0089] When X-ray imaging is performed, the sample needs to be fully dried and the structure integrity of the cells needs to be maintained, and the morphology of the cells does not change and does not collapse. After searching the prior art, the inventors found that there is no plant sample preparation method for X-ray imaging in the prior art, especially for young nodule with high water content (>80%), direct drying will cause cell collapse and structure damage.

[0090] To solve the technical problem, the inventors unexpectedly found that the nodule sample prepared by the method meets the sample requirements of X-ray imaging, and high-resolution X-ray three-dimensional images meeting the requirements can be obtained.

[0091] The application of X-ray imaging technology in plant tissue imaging can open up new ideas for researchers, and the current X-ray three-dimensional imaging technology requires a long scanning time, and the imaging cost of a sample is about 20-30 thousand yuan. The sample preparation method of the present application will provide a reference for other researchers and save a lot of time and cost.

[0092] In addition, although there are individual laboratory researchers who have observed legume root nodules using X-ray microscopy technology in the prior art, the results obtained by them are overall imaging of the root system, and the resolution does not reach the cell level (for example, see Teramoto et al. High-throughput three-dimensional visualization of root system architecture of rice using X-ray computed tomography, Plant Methods, 2020 16:66; Van Dusschoten, D., Metzner, R., Kochs, J., Postma, J. A., Pflugfelder, D., Bühler, J., Schurr, U., Jahnke, S. (2016). Quantitative 3D analysis of plant roots growing in soil using magnetic resonance imaging. Plant physiology, 170(3), 1176-1188). The preparation method of the present application can make the resolution of X-ray imaging reach the cell level.

[0093] Root nodule sample preparation: vacuum infiltration (-0.08 MPa) with 50% formaldehyde-acetic acid-ethanol (FAA) until no bubbles come out of the sample, place the sample on ice, fix the sample overnight at room temperature, perform gradient ethanol dehydration (50%→70%→80%→90%→95%→100%) in 15 ml centrifuge tubes, each for 15 min, on ice. Replace 100% ethanol with 3:1 volume ratio of ethanol: acetone for 30 min, on ice; replace with 3:2 volume ratio of ethanol: acetone for 30 min, on ice; replace with acetone for 30 min; repeat acetone replacement for 30 min, place the sample in a sample basket and put it into a critical point drying instrument, replace with liquid CO2(flow rate 10 L / min for 1 h), slowly exhaust after warming up to 31.1°C (supercritical state), complete critical point drying.

[0094] After ensuring that the sample is completely dry and the morphology is intact, it is placed in the X-ray microscope for imaging. The X-ray microscope used is Zeiss Xradia 620 (voltage 30-160 kV, 25 W; spatial resolution 0.5 microns; 0.4x, 4x, 20x, 40x objective) of Carl Zeiss Management Co., Ltd.

[0095] II. Spatial group sequencing

[0096] 1. Sample preparation and embedding

[0097] Fresh / frozen soybean root nodules were placed in OCT embedding medium, placed on ice, and vacuumed for 10 minutes. Directional embedding was performed in embedding boxes using OCT embedding medium, and the samples were frozen on dry ice. RNase-free operation was performed throughout the process.

[0098] 2. Frozen section preparation

[0099] The constant temperature freezing microtome was set to -20°C, and the section thickness was 10 μm.

[0100] 3. RNA quality detection

[0101] 10-20 10 μm thick tissue sections were stored in 1.5 ml EP tubes pre-cooled at -20°C, and total RNA extraction (Beijing Huayueyang Biological, 0416) and quality detection (RIN≥7 was considered as a qualified sample) were performed. The remaining tissue was continuously embedded with OCT, and after the quality inspection was confirmed, the sections were subjected to formal spatiotemporal transcriptome experiment.

[0102] 4. Chip processing

[0103] The reagents and consumables used were provided by Huada (Stereo-seq Transcriptome Kit T, 201KT13114; STOmics Accessory Kit, 1000033700; Stereo-seq PCR Adaptor, 301AUX001).

[0104] Carefully pick up the spatially resolved sequencing chip (Stereo-seq Chip T Carrier (1 cm * 1 cm), Cat. No. 200CT13114, see Chen, A., Liao, S., Cheng, M., Ma, K., Wu, L., Lai, Y., Qiu, X., Yang, J., Li, W., Xu, J., et al. (2022). S patiotemporal transcriptomic atlas of mouse organogenesis using DNA nanoball patterned arrays. Cell 185: 1777-1792. https: / / doi.org / 10.1101 / 2021.01.17.427004) with tweezers and place it in a new 24-well plate. Wash twice with 400 μΐ of 0.1x SSC; cover the chip surface with 100 μΐ of NF-H20 (with RNase inhibitor) for 1 min to remove RNase, and rinse once with 400 μΐ of NF-H20 (without RNase inhibitor); dry the chip surface, and the patch is ready. Carefully cover the frozen section on the chip and quickly bake the section at 37 °C for 3 min.

[0105] 5. Tissue fixation

[0106] Place the chip in a pre-cooled methanol solution at -20 °C for 30 min; remove the fixed tissue chip from the methanol and place the chip in a culture dish with parafilm; after the methanol evaporates and the chip surface dries, wash the chip once with 100 μΐ per chip of 0.1x SSC and keep the chip surface moist; after the 0.1x SSC is absorbed, take a 1 / 20 RI fluorescence photo, select the FITC channel on the microscope, 10x lens, and take a full picture. Figure 3 After the liquid on the chip surface is absorbed, add 100 μΐ per chip of 0.1x SSC, absorb the liquid on the chip surface, and immediately proceed to the next step of permeabilization.

[0107] 6. Tissue permeabilization

[0108] Freshly prepare the permeabilization reagent 0.01N HC1 (100 μl 0.1N HC1 + 900 μl NF-H2O), prepare the dilution of permeabilization enzyme according to the number of chips, each chip uses 1.1 μl permeabilization enzyme added to 110 μl 0.01N HC1, mix well with pipette, remember not to shake. Put the prepared permeabilization reagent in a 37°C incubator for 3 minutes; transfer the chip culture dish to a 37°C incubator, add 100 μl / chip of preheated and equilibrated permeabilization reagent for 8 minutes (preliminary experiments found that 8 minutes of permeabilization in situ information is good, 12 minutes of permeabilization has obvious nucleic acid diffusion, such as Figure 2 shown); aspirate the permeabilization reagent on the surface of the chip, add 100 μl / chip of 0.1x SSC; aspirate the liquid on the surface of the chip and transfer it to a new culture dish, immediately add the next reaction solution.

[0109] 7. In situ reverse transcription reaction

[0110] Gently add 100 μl of RT Mix equilibrated to room temperature to the chip, and react in a 42°C incubator for 1 hour.

[0111] 8. Tissue removal and product recovery

[0112] Aspirate the RT Mix on the surface of the chip, wash the chip once with 0.1x SSC; add the tissue removal buffer, and react in a 37°C incubator for 30 minutes. Rinse the chip with 0.1x SSC and repeat the washing twice for the next reaction.

[0113] 9. RT product recovery and amplification

[0114] cDNA release: add 400 μl / well of cDNA Release Mix to the chip, and seal it in a 55°C incubator for 3 hours, shaking every hour.

[0115] cDNA recovery: recover the liquid in the reaction well to a new 1.5 ml centrifuge tube, and add 350 μl / well of NF-H2O to wash the reaction well chip and recover it.

[0116] 0.8x beads purification: mix the recovered liquid with VAHTSTM DNA Clean Beads equilibrated to room temperature, shake well, incubate at room temperature for 10 minutes; after centrifugation, place the EP tube on a magnetic stand for 3 minutes, remove the supernatant after the liquid is clear; add 1 ml of freshly prepared 80% ethanol (presently prepared), stand for 30 seconds; discard the supernatant and repeat the above steps; discard the supernatant and air dry at room temperature for 5 minutes; add NF-H2O with a volume of 42 μl to dissolve, mix well and stand at room temperature for 5 minutes, centrifuge for a moment, and stand on the magnetic stand for 3-5 minutes, then recover 42 μl of supernatant after the liquid is clear.

[0117] Prepare PCR Mix: 42 μΐ of the above recovered sample, supplemented with NF-H20 to make up the volume to 42 μΐ, add 58 μΐ of PCR Mix, total volume 100 μΐ, and perform PCR reaction; determine the concentration of PCR product.

[0118] 10. cDNA purification

[0119] Mix the PCR reaction solution with VAHTSTM DNA Clean Beads (VAZYME), shake to mix, incubate at room temperature for 10 min; after centrifugation, place on a magnetic stand for 3 min, remove the supernatant; add 1 ml of 80% ethanol, stand for 30 s; discard the supernatant and repeat the steps; discard the supernatant and dry at room temperature for 3 min; add 42 μΐ of NF-H20 to restore the volume, mix well and stand for 5 min, centrifuge briefly, place on a magnetic stand for 3 min, and recover the supernatant; take 1 μΐ of cDNA sample to determine the concentration and cDNA fragment distribution.

[0120] 11. cDNA fragmentation and amplification

[0121] Prepare the fragmentation and amplification system, mix 100 μΐ of PCR amplification product and VAHTSTM DNA Clean Beads (VAZYME), shake to mix, incubate at room temperature for 5 min; after brief centrifugation, place on a magnetic stand for 3 min, transfer the supernatant; mix 20 μΐ of VAHTSTM DNA Clean Beads with the supernatant, shake to mix, incubate at room temperature for 5 min; after brief centrifugation, place on a magnetic stand for 3-5 min, aspirate and discard the supernatant; add 200 μΐ of 80% ethanol, stand for 30 s, aspirate and discard the supernatant carefully; centrifuge briefly, separate on a magnetic stand, and aspirate the liquid at the bottom of the tube; dry the magnetic beads at room temperature for 5 min; add 20 μΐ of NF-H20 to restore the volume, stand at room temperature for 5 min, centrifuge briefly, place on a magnetic stand for 3 min, and transfer the supernatant to a new PCR tube; take 1 μΐ of the screened product to determine the concentration.

[0122] 12. Two-stage fragment screening

[0123] The screened product is detected for fragment distribution using an Agilent 2100 High Sensitivity Chip. Quality control standards: yield greater than 300 ng, and main peak of fragment distribution between 400-600.

[0124] 13. Sequencing preparation

[0125] Library concentration: Qubit quantification > 2 nM; Sequencing on DNBSEQ T1 sequencing platform, MGISEQ-2000 of BGI-research or DNBSEQ-T7 sequencing platform (Beijing) in paired-end mode. Each sequencing read contains a 50 bp Read 1 and a 100 bp Read 2. Spatial coordinate information is matched by barcode.

[0126] 14. Stereo-seq data preprocessing

[0127] Stereo-seq raw data were determined coordinate identities (CIDs) by mapping to the coordinates of the in situ capture chip (allowing one mismatch to tolerate sequencing or PCR errors). After removing reads with UMI quality scores below 10, the remaining reads were aligned to the reference genome (Gmax_ZH13_V2.0) and annotated to the gene set (version 2.1). Finally, a transcriptome matrix containing CID information was generated for subsequent analysis. Since the plant tissue only covers part of the surface of each Stereo-seq chip, a custom script was used to extract the exact area of the spatial transcriptome data. For this purpose, the spatial transcriptome data were aligned with the ssDNA image, and then the boundary was determined according to the sample edge (an irregularly arranged unequal interval straight line on the surface of the Stereo-seq chip, which is not planted with DNA oligonucleotides for mRNA capture). Due to the characteristics of mRNA captured on a single DNA nanoball (DNB) at a subcellular resolution (500 nm), the spatial transcriptome matrix was convolved into an N x N size spatial transcriptome block (named binN, N = 50). Further, the spatial transcriptome matrix was convolved into a cell-level spatial transcriptome block using ssDNA-based cell segmentation technology. Using the Patchify package in Python (v3.12.10), the ssDNA image was segmented into sub-images of the same size (default length and width are both 512 pixels), saving computing resources and facilitating calculation; the cyto2 algorithm in the Cellpose package was used to identify cells in each sub-image; the patchify package was used to merge the cell segmentation results of each sub-image, and the merged cell segmentation results were output.

[0128] III. Single nucleus sequencing

[0129] 1. Soybean nucleus extraction

[0130] Nuclei were extracted from root nodules of soybean for snRNA-seq. Fresh tissues harvested from plants were immediately chopped in a petri dish with a scalpel blade, and 2 ml of pre-chilled homogenization buffer (Sigma, CELLYTPN1-1KT) was added. Then, 1 ml of the mixture was transferred to a 1 ml Dounce homogenizer (Kimble), frozen in a freezer, and homogenized for 10-20 times. The homogenized mixture was then filtered through a 40-um cell strainer into a 1.5-ml centrifuge tube, centrifuged at 1260 g for 5 min at 4 °C to pellet the nuclei, and then resuspended in blocking buffer containing 1% bovine serum albumin (BSA) and 0.2 U / µL of RNase inhibitor. Finally, the nuclei were re-pelleted by centrifugation at 600 x g for 5 min at 4 °C, and then resuspended in 1x phosphate-buffered saline (PBS) containing 1% BSA for subsequent library preparation.

[0131] 2. snRNA-seq library construction and sequencing

[0132] snRNA-seq libraries were constructed using the DNBelab C series high-throughput single-cell system (BGI-research). The snRNA-seq library preparation steps: single nuclei suspension by titration microfluidic encapsulation, emulsion breaking, mRNA capture bead collection, reverse transcription, cDNA amplification and purification. The resulting cDNA product was fragmented into 250-400 bp fragments. The indexed sequencing library was constructed using the Stereo-seq library construction kit (BGI-Research, 101KL114) according to the standard protocol, and quantified using the Agilent Bioanalyzer 2100 and Qubit single-stranded DNA (ssDNA) assay kit (Thermo Fisher Scientific). Sequencing was performed on the MGISEQ-2000 or DNBSEQ-T7 sequencing platform (Beijing) of BGI-research in paired-end mode. Each sequencing read includes a 30 bp Read 1 (including 10 bp cell barcode 1, 10 bp cell barcode 2, and 10 bp unique molecular identifier, UMI), a 100 bp Read 2 gene sequence, and a 10 bp sample index barcode read.

[0133] 3. snRNA-seq data preprocessing

[0134] snRNA-seq data were processed using open source software (https: / / github.com / MGI-tech-bioinformatics / DNBelab_C_Series_HT_scRNA-analysis-software) to obtain read count matrix for each gene and each cell. First, raw sequencing reads were filtered (reads with average base quality score lower than 4, more than 2 bases with quality score lower than 10, containing N base or inappropriate barcodes were removed) and demultiplexed by barcodes assignment. Then, obtained reads were aligned to soybean reference genome (Gmax_ZH13_V2.0) using Spliced Transcripts Alignment to a Reference (STAR; v2.7.9) and annotated to gene sets using PISA (version 2.1). Valid cells were automatically identified using the “barcodeRanks” function of DropletUtils tool based on the UMI number distribution of each cell to remove background beads and cells with UMI counts less than 500. Gene expression of cells was then calculated using PISA and gene x cell matrix was created for each library. Expression data of different libraries from the same tissue were merged using the “Merge” function of R package Seurat (v4.1.0) and subsequently normalized using the “NormalizeData” function with default parameters. After log-normalization by the “ScaleData” function, 1500 highly variable genes were selected using the “FindVariableFeatures” function. For subsequent clustering and visualization, 40 PCs were extracted using principal component analysis (PCA)-based dimensionality reduction algorithm, and Lovain clustering was performed using the “FindNeighbors” function with resolution set to 1. The clustering result was projected to two-dimensional space using the “RunUMAP” function. Differentially expressed genes were identified using the “FindAllMarkers” function.

[0135] IV. Integration and annotation of single-nuclei sequencing data and spatially-resolved sequencing data

[0136] Cell types were annotated by integrating snRNA-seq data and Stereo-seq data. First, we used the “FindTransferAnchors” and “TransferData” functions of the R package Seurat (v4.1.0) to calculate the proportion of cell clusters determined by snRNA-seq that were located to each cluster determined by spatial transcriptomics. Projection from UMAP of snRNA-seq to cluster embedding of Stereo-seq was done by the “Knn” function. In this way, we obtained the spatial location of clusters determined by snRNA-seq and annotated them according to their anatomical structure. Some other cell clusters were identified based on the lineage relationship of previously annotated cell clusters using the R package Clustertree (v1.8.7).

[0137] V. 3D single-cell spatial ensemble expression profiling

[0138] We used X-ray microscopy for multiscale high-resolution imaging (1.5 microns) as a reference for reconstructing the 3D structure of the sample. The tissue images in the circular field of view of the continuous scanning pictures were extracted using the Segment Editor module in 3D Slicer (v5.8.1). Based on the extracted tissue images of the continuous sections, we used 3D Slicer to construct the 3D Mesh of the outer contour of the tissue and smooth it, storing the final 3D Mesh as an obj format. Based on the first 40 principal components, we used the Louvain algorithm to cluster the single-cell transcriptomic data (resolution = 1) to identify the cell type features of the tissue. We used the FindTransferAnchors and TransferData functions in Seurat (v4.1.0) to map the Stereo-seq data to the single-cell data, and according to the mapping score of each Squarebin of Stereo-seq data, we determined its cell type. According to the section position in the Stereo-seq continuous section experiment record, we calculated the Z-axis position of the corresponding 3D Mesh of Stereo-seq. Using the TrakEM2 elastic registration algorithm in ImageJ (v1.53c), we registered the Stereo-seq continuous sections with the X-ray images of the corresponding Z-axis position, thereby determining the spatial position of the Stereo-seq continuous sections in the 3D Mesh. Finally, we can display the gene expression and tissue type information on the Stereo-seq continuous sections at the accurate spatial position within the 3D Mesh of the tissue.

[0139] Results and analysis of results

[0140] Nodules were fixed with formalin-acetic acid-ethanol (FAA, 50%) and dehydrated in a series of ethanol concentration gradient, and then dried using a CO2 dryer. We obtained completely dried and intact nodule samples. Further, we used a Zeiss Xradia 610 Versa X-ray microscope for multiscale high-resolution imaging (1.5 microns) as a reference for reconstructing the 3D structure of the sample Figure 1 ). The pre-experiment found that the nodule tissue section was well-informed after 8 minutes of permeabilization, and there was obvious nucleic acid diffusion after 12 minutes of permeabilization Figure 2 ). We performed serial sectioning of the entire nodule and selected 21 frozen sections for spatial group sequencing Figure 3 ). After deconvolution, the nodule spatial expression profile was processed with bin50, and a total of 15454 bins were detected, with 4368 mRNAs detected in each bin, 2183 genes, and a total of 37086 genes detected Figure 4 ). According to the spatial positioning and gene expression pattern of the nodule, the Stereo-seq spatial group data showed 12 cell groups with specific spatial positioning Figure 5 ). After the extraction of single nuclei from the same batch of samples, nodule single-cell sequencing was performed, and a total of 10113 nuclei were obtained. Through the snRNA-seq data analysis process, 10,113 high-quality nuclei were divided into 26 cell types Figure 7 ). Most cell groups had cell-specific marker genes, indicating that the grouping was clear and reasonable Figure 8 ). Among the 26 cell types, 16 were mapped to the specific cell positioning of the spatial cluster, so we annotated the single-cell groups based on spatial positioning Figure 9 and 7). The merging of each cell group with clear positioning showed clear nodule cell group positioning Figure 10). Further, tissue images in the circular field of view of the serial scan pictures were extracted using the Segment Editor module in 3D Slicer (v5.8.1). Based on the extracted tissue images of serial sections, a 3D Mesh of the tissue outer contour was constructed using 3D Slicer, and was smoothed, and the final 3D Mesh was stored as an obj format. Based on the first 40 principal components, the single-cell transcriptome data was clustered using the Louvain algorithm (resolution = 1) to identify the cell type signature of the tissue. The Stereo-seq data was mapped to the single-cell data using the FindTransferAnchors and TransferData functions in Seurat (v4.1.0), and the cell type of each Squarebin in the Stereo-seq data was determined according to the mapping score. According to the section position in the Stereo-seq serial section experiment record, the Z-axis position of the Stereo-seq corresponding 3D Mesh was calculated. Using the TrakEM2 elastic registration algorithm in ImageJ (v1.53c), the Stereo-seq serial sections were registered with the X-ray images at the corresponding Z-axis position, thereby determining the spatial position of the Stereo-seq serial sections in the 3D Mesh. Finally, we can display the gene expression and tissue type information on the Stereo-seq serial sections at the accurate spatial position in the tissue 3D Mesh.

[0141] We constructed a single-cell data 3D model of nodule based on the tissue cell model of X-ray microscope, the section spatial group data of Stereo-seq, and the single-nucleus transcriptome data information, and these data were stored in the Soybean Organ Transcriptome Atlas (SOTA) (https: / / db.cngb.org / stomics / soybean / nodule.3d / ). The model allows users to explore the composition of nodule cell groups and gene expression patterns in three dimensions.

[0142] REFERENCES

[0143] Fan J., Shen Y., Chen C., Chen X., Yang X., Liu H., Chen R., Liu S., Zhang B., Zhang M., Zhou G., Wang Y., Sun H., Jiang Y., Wei X., Yang T., Liu Y., Tian D., Deng Z., Xu X., Liu X., and Tian Z. (2025). A large-scale integrated transcriptomic atlas for soybean organ development. Mol. Plant. 18, 669-689.

[0144] The above-described embodiments of the present application are merely intended to further illustrate the purposes, technical solutions and beneficial effects of the present application, and should not be used to limit the present application. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present application should be included in the protection scope of the present application.

Claims

1. A method for three-dimensional reconstruction of the spatial transcriptome of single cells of leguminous plants such as soybean nodules, characterized in that: The method includes the steps of integrating single cell data and spatial group data, and then integrating X-ray imaging data. Optionally, the single cell data is obtained by snRNA-seq technology. Optionally, the spatial group data is obtained by Stereo-seq technology. Optionally, the X-ray imaging data is obtained including the step of drying the sample.

2. The method according to claim 1, characterized in that The snRNA-seq technique comprises the following steps: a1. Single nucleus suspension, emulsion rupture, mRNA capture bead collection, reverse transcription, cDNA amplification and purification by titrating microfluidics packages; a2. The resulting cDNA product is fragmented into 250-400 bp fragments; a3. Construct an indexed sequencing library and perform sequencing in paired-end mode, wherein each sequencing read includes a 30 bp Read 1, a 100 bp Read 2 gene sequence, and a 10 bp sample index barcode read. Optionally, the 30 bp Read 1 includes a 10 bp cell barcode 1, a 10 bp cell barcode 2, and a 10 bp unique molecular identifier. a4. Preprocess the snRNA-seq data.

3. The method according to claim 2, characterized in that The snRNA-seq data were preprocessed by filtering the raw sequencing reads and demultiplexing them by barcode assignment; the obtained reads were aligned to the soybean reference genome, such as Gmax_ZH13_V2.0, using Spliced ​​Transcripts Alignment to a Reference (STAR), annotated to the gene set using PISA; the "barcodeRanks" function of the DropletUtils tool was used to automatically identify valid cells, and background beads and reads with less than 500 UMI counts were removed based on the distribution of the number of unique molecular identifiers per cell; PISA was then used to calculate the gene expression of the cells and create a gene × cell matrix for each library; the expression data of different libraries of the same tissue were merged using the "Merge" function of the R package Seurat, and then normalized using the "NormalizeData" function with default parameters; after logarithmic normalization using the "ScaleData" function, 1500 highly variable genes were selected using the "FindVariableFeatures" function; for subsequent clustering and visualization, a dimensionality reduction algorithm based on principal component analysis (PCA) was used to extract multiple, for example 40, PCs, and the "FindNeighbors" function was used for Lovain clustering with a resolution of 1; the clustering results were projected into two-dimensional space using the "RunUMAP" function; and the "FindAllMarkers" function was used to identify differentially expressed genes.

4. The method according to claim 3, characterized in that Filtering of raw sequencing reads included removing reads with an average base quality score lower than 4, more than 2 bases with a quality score lower than 10, containing N bases, or inappropriate barcodes.

5. The method according to any one of claims 1 to 4, characterized in that The Stereo-seq technique includes the following steps: b1. Sample preparation and embedding; b2. Frozen section preparation; b3. RNA quality testing of frozen sections. Optional: RIN ≥ 7 is considered acceptable. b4. Chip processing; b5. Tissue fixation; b6 tissue permeabilization, optionally, the permeabilization time is 8 minutes; b7. In situ reverse transcription reaction; b8. Reverse transcription product recovery and amplification; b9. Purification of reverse transcription products; b10. Reverse transcription product interruption and amplification; b11 double fragment screening, optionally, the quality control standard is a yield greater than 300ng, the main peak of the fragment distribution is between 400-600; b12. Sequencing, optionally with a library concentration of ≥ 2 nM using Qubit quantification, and paired-end sequencing, where each sequencing read contains a 50 bp Read 1 and a 100 bp Read 2, with spatial coordinate information matched by barcodes; b13. Stereo-seq data preprocessing.

6. The method according to claim 5, characterized in that Stereo-seq data preprocessing includes mapping the raw Stereo-seq data to the coordinates of the in situ capture chip to determine the coordinate identifiers, removing reads with UMI quality scores lower than 10, aligning the remaining reads to the reference genome, such as Gmax_ZH13_V2.0, and annotating them to gene sets, and finally generating an expression matrix containing coordinate identifier information for subsequent analysis. Optionally, one mismatch is allowed when mapping the raw Stereo-seq data to the coordinates of the in situ capture chip to tolerate sequencing or PCR errors. Optionally, since plant tissue only covers a portion of the surface area of ​​each Stereo-seq chip, the expression matrix is ​​not fully mapped. domain, extracting precise regions of spatial transcriptome data using a custom script, aligning the spatial transcriptome data with the ssDNA image, and then determining boundaries based on sample edges, convolving the spatial expression profile matrix into N×N spatial expression profile blocks and naming them binN, where N=50, and further using ssDNA-based cell segmentation technology to convolve the spatial expression profile matrix into cellular-level spatial expression profile blocks; optionally, the ssDNA-based cell segmentation technology to convolve the spatial expression profile matrix into cellular-level spatial expression profile blocks includes using the Patchify package in Python to segment the ssDNA image into sub-images of equal size; The cyto2 algorithm in the Cellpose package is used to identify cells in each sub-image; the pathchify package is used to merge the cell segmentation results identified in each sub-image and output the merged cell segmentation results.

7. The method according to any one of claims 1 to 6, characterized in that Integration of single-cell data and spatial group data involves integrating pre-processed snRNA-seq data with Stereo-seq data. Optionally, integration of pre-processed snRNA-seq data with Stereo-seq data involves calculating the proportion of cell clusters identified by snRNA-seq localized to each cluster identified by spatial transcriptomics using the "FindTransferAnchors" and "TransferData" functions of the R package Seurat. Cluster embedding from UMAP projection of snRNA-seq to Stereo-seq is performed using the "Knn" function to obtain the spatial positions of clusters identified by snRNA-seq and annotate them according to their anatomical structures. Additional cell clusters are identified based on the lineage relationships of previously annotated cell clusters using the R package Clustertree.

8. The method according to any one of claims 1 to 7, characterized in that The steps of reintegrating X-ray imaging data include using the Segment Editor module in 3D Slicer to extract tissue images in the circular field of view of continuous scanning images; based on the extracted tissue images of continuous slices, use 3D Slicer to construct a 3D mesh of the outer contour of the tissue, smooth it, and store the final 3D mesh in obj format; based on the first 40 principal components, use the Louvain algorithm to cluster the single-cell transcriptome data (resolution = 1) to identify the cell type characteristics of the tissue; use the FindTransferAnchors and TransferData functions in Seurat to map the Stereo-seq data to the single-cell data, and determine its cell type based on the mapping score of each Squarebin of the Stereo-seq data; calculate the Z-axis position of the Stereo-seq corresponding 3D Mesh based on the slice position in the Stereo-seq continuous slice experiment record; use the TrakEM2 elastic registration algorithm in ImageJ to align the Stereo-seq continuous slices with the X-ray images of the corresponding Z-axis position to determine the Stereo-seq continuous slices in 3D The spatial position of the nodule in the Mesh; finally, the information on the Stereo-seq continuous slices, such as gene expression and tissue type information, is displayed at the accurate spatial position within the tissue 3D Mesh, thereby obtaining a three-dimensional reconstruction map of the spatial transcriptome of the nodule organ.

9. The method according to any one of claims 1 to 8, characterized in that The steps for drying the sample include the following: c1. Vacuum infiltrate the sample with 50% formaldehyde-acetic acid-ethanol (FAA) at -0.08 MPa until no bubbles appear. Place the sample on ice and fix it overnight at room temperature. Dehydrate the sample in a 15 ml centrifuge tube using a gradient of ethanol: 50% → 70% → 80% → 90% → 95% → 100% → 100% → 100% ethanol. Each gradient is 15 minutes long and the sample is kept on ice. c2. Pour off the 100% ethanol and replace with 3:1 ethanol:acetone by volume for 30 minutes. Place on ice. c3. Replace with ethanol:acetone in a 3:2 volume ratio for 30 minutes and place on ice. c4. Replace with acetone for 30 min; c5. Repeat the acetone replacement for 30 minutes. Place the sample in a basket and place it in the critical point drying instrument. Replace the sample with liquid CO2. Raise the temperature to 31.1°C and slowly exhaust the air to complete critical point drying.

10. A sample preparation method suitable for X-ray imaging of plant root nodules, comprising the following steps: c1. Vacuum infiltrate the sample with 50% formaldehyde-acetic acid-ethanol (FAA) at -0.08 MPa until no bubbles appear. Place the sample on ice and fix it overnight at room temperature. Dehydrate the sample in a 15 ml centrifuge tube using a gradient of ethanol: 50% → 70% → 80% → 90% → 95% → 100% → 100% → 100% ethanol. Each gradient is 15 minutes long and the sample is kept on ice. c2. Pour off the 100% ethanol and replace with 3:1 ethanol:acetone by volume for 30 minutes. Place on ice. c3. Replace with ethanol:acetone in a 3:2 volume ratio for 30 minutes and place on ice. c4. Replace with acetone for 30 min; c5. Repeat the acetone replacement for 30 minutes. Place the sample in a basket and place it in the critical point drying instrument. Replace the sample with liquid CO2. Raise the temperature to 31.1°C and slowly exhaust the air to complete critical point drying.

Images