空间与单细胞转录组学揭示大豆对豆薯层锈菌侵染的防御反应中两种独特细胞状态及其空间协同机制

本综述通过整合空间转录组学（Visium）与单核RNA测序（snRNA-seq）技术，系统解析了大豆（Glycine max）应答亚洲大豆锈病（ASR）病原体（Phakopsora pachyrhizi）侵染过程中的空间异质性免疫调控机制。研究发现感染区域与周边区域存在截然不同的转录组状态：感染区虽存在病原体但防御基因表达较低，而周边区虽病原体稀少却表现出强烈的非自主性防御反应（cell non-autonomous defense）。研究进一步通过加权基因共表达网络分析（hdWGCNA）鉴定出一个关键免疫应答模块（M3），其在应激细胞中特异性激活。该工作为理解植物－病原互作的空间协调机制提供了新视角，并为作物抗病育种提供了潜在靶点。

　　

1 引言

植物作为固着生物，必须适应包括生物和非生物胁迫在内的局部环境刺激。与动物利用移动免疫细胞清除病原体不同，植物依赖与病原体接触的局部细胞来检测和启动免疫反应。虽然批量RNA测序（RNA-seq）能够在全转录组水平评估植物对病原体感染的反应，但植物免疫反应的空间协调机制仍不清楚。本研究通过空间和单核转录组学实验，捕捉大豆对豆薯层锈菌（Phakopsora pachyrhizi）引起的亚洲大豆锈病（ASR）感染的反应空间模式。

2 材料与方法

2.1 植物生长与病原接种

10天龄Williams 82大豆植株接种豆薯层锈菌分离株FL-07（90-110K孢子/ml），于4天和7天后采集样本。单核RNA测序使用14天龄93Y21大豆幼苗，喷洒接种GA-05分离株（约100K孢子/ml），并于1、3、5天后采集受感染单叶。

2.2 单核RNA测序

2.2.1 组织收集

约40 mg新鲜ASR感染大豆叶片组织速冻于液氮，-80°C保存直至处理。

2.2.2 核分离

采用改良的已发表方案。冷冻组织在含有裂解缓冲液的预冷培养皿中切碎，使用Dounce组织研磨器匀浆，经20 μm和10 μm细胞筛过滤，离心后重悬核沉淀。

2.2.3 10x Genomics Chromium文库制备

分离的核用血细胞计数器计数，每个样本约一万个核加载到Chromium Next GEM芯片上，按照制造商说明制备文库。

2.2.4 测序

Chromium Next GEM文库在Illumina NovaSeq 6000上测序，读长按10x Genomics推荐，每个样本测序深度为3亿簇。

2.3 Visium空间转录组学

2.3.1 组织收集

叶片从植株上切下并进一步切成小矩形，置于填充Sakura OCT包埋剂的包埋模具中，定向为横切或平切面，于液氮冷却的异戊烷中速冻，-80°C保存。

2.3.2 冰冻切片和荧光成像

样本块移至Leica CM1950冰冻切片机， chamber温度设为-13°C。切10 μm厚切片，在Leica DM1000上观察真菌结构。为增加冰冻切片完整性，采用改良的 plastering技术，使用1:1稀释OCT的10 mM磷酸盐缓冲盐水（PBS）。切片用小麦胚芽凝集素Alexa Fluor 488缀合物（1:2000）在PBS中4°C染色过夜，80%甘油封片，在Leica Upright DM5000上用20x空气物镜荧光成像。

2.3.3 切片固定、染色和成像

Visium载玻片从-80°C取出，置于37°C载玻片加热器上1分钟。然后置于冰上塑料托盘中。预冷至-20°C的Farmer固定液（3份100%乙醇:1份冰醋酸）加至载玻片，静置2分钟。倾去固定液，加预冷至-20°C的100%乙醇，静置2分钟。倾去乙醇，干燥2分钟。加1:10稀释的1%甲苯胺蓝（溶于1%硼酸钠溶液），静置2分钟。倾去染液，加蒸馏水静置2分钟。倾去水，干燥5分钟后成像。载玻片在Leica Thunder Imager上用10x/0.30NA空气物镜和K5单色相机透射光成像。

2.3.4 10x Genomics Visium组织优化

组织优化根据制造商说明（CG000238 Rev F）进行。10 μm厚切片置于Visium组织优化载玻片上，按2.3.3处理成像，然后测试不同透化时间（5-15分钟）。使用Visium空间组织优化试剂盒试剂进行逆转录和组织去除。处理后在Leica Thunder Imager上荧光成像。根据最大荧光信号和最低信号扩散选择10分钟透化时间，用于后续Visium基因表达载玻片处理。

2.3.5 10x Genomics Visium基因表达文库制备

固定在Visium基因表达载玻片上的10 μm厚切片按制造商说明处理，使用Visium空间基因表达试剂盒和文库构建试剂盒分别进行cDNA扩增和文库构建。每个切片37°C透化10分钟，53°C逆转录45分钟。第二链合成和变性后，从单个切片回收cDNA片段，根据回收cDNA片段的定量PCR（qPCR）优化步骤确定的Cq值扩增17个循环。酶切、末端修复、接头连接和大小选择后，使用10x Genomics Dual Index TT Set A试剂盒中的单独条形码引物扩增所得cDNA插入片段，测序前合并。

2.3.6 测序

在Illumina NovaSeq 6000和NextSeq 550系统上按10x Genomics推荐的读长和深度（每个点50,000簇）进行测序。使用来自Read 1测序数据的单独空间条形码和唯一分子标识符（UMI）序列进行空间定位和独特分子事件表征，而Read 2测序数据用于确定原始cDNA片段和转录本的分子身份。

2.4 Visium空间和Chromium snRNA-seq数据的处理与分析

2.4.1 读段比对

使用10x Genomics Space Ranger软件（v2.0.1）和10x Genomics Cell Ranger软件（v8.0.1）将原始读段处理为每个条形码Visium点和每个条形码Chromium核的计数值。两个实验的读段同时比对到豆薯层锈菌基因组（Phapa1, GCA_025201825.1）和大豆Williams 82（Wm82v4）基因组。使用R 4.3.2中的Seurat 5.1.0包进行单核和空间数据分析。

2.4.2 样本过滤

手动对齐每个Visium载玻片上捕获区域周围的基准框，使用10x Genomics Loupe Browser软件（v8.1.2）选择组织覆盖点。此后，12个样本5,085个点上的33,994个基因（28,789个来自大豆，5,196个来自豆薯层锈菌）用于进一步分析。

选择高质量核（300 < UMI < 10,000; 250 < 基因 < 6000; 线粒体或叶绿体读段少于10%）用于进一步分析。使用DoubletFinder（v2.0.4）识别和去除双联体。过滤后，6个样本48,358个核上的44,478个基因（41,028个来自大豆，3,450个来自ASR）用于进一步分析。

2.4.3 标准化和聚类

空间和单核样本分别以相同方式处理，产生两个独立的Seurat对象。使用Seurat SCTransform标准化每个样本，然后对所有捕获区域或所有单核样本进行典型相关分析（cca）整合。使用整合数据对象，所有点或细胞使用30个主成分（PCs）和Louvain算法聚类，空间和单核对象的分辨率分别为0.5和0.25。还用30个PCs计算均匀流形近似和投影（UMAP）降维。

从完整单核对象中提取应激细胞类型簇，其中Louvain聚类以0.25分辨率运行。此后，识别所得亚簇的标记基因，并与完整数据集中识别的细胞簇比较。保留原始UMAP确认亚簇细胞类型注释。

2.4.4 基因水平分析

使用Seurat函数FindAllMarkers识别簇标记基因，显著性由Wilcoxon秩和检验和调整P < 0.05定义。使用Seurat函数FindMarkers识别健康簇和感染簇之间的差异表达基因（DEGs），显著性由Wilcoxon秩和检验和调整P < 0.05定义。使用fgsea（v1.28.0）进行基因本体（GO）富集，按上述函数计算的所有avg_log2FC值排序。使用Seurat函数AddModuleScore总结最初由Cabre等人（2021）和Morales等人（2013）识别的15个ASR诱导基因的表达。

2.4.5 标准化伪批量分析

使用Seurat函数AggregateExpression对空间和单核数据集每个基因的原始计数求和。求和的原始计数提供给DESeq2（v1.42.1）进行默认标准化。

2.4.6 可视化

使用多个R包创建数据可视化。EulerR（v7.0.2）和SuperExactTest（v1.1.0）用于创建维恩图和相关统计。UpsetR（v1.4.0）用于创建Upset图。Pheatmap（v1.0.12）用于创建所有热图。

2.5 高维加权基因共表达网络分析

使用R 4.3.3中的hdWGCNA 0.3.03包进行高维加权基因共表达网络分析（hdWGCNA）。网络构建前，过滤基因至在至少1%细胞中表达的大豆基因（24,705个基因），过滤细胞至标记为健康或应激叶肉细胞。使用ConstructNetwork函数构建带符号共表达网络，minModuleSize设为50，detectCutHeight设为0.995，mergeCutHeight设为0.2。为选择软功率阈值，按Morabito等人描述使用TestSoftPowers函数执行参数扫描。

为总结每个共表达模块中基因的表达，计算模块特征基因（MEs）。MEs定义为模块基因表达矩阵的第一主成分，模块中心基因通过与模块特征基因的相关性选择。对每个模块与每个模块前10个中心基因的拓扑重叠执行以基因为中心的UMAP降维。为比较模块特征基因值的统计显著性，采用未配对双侧Wilcoxon秩和检验。

3 结果

3.1 空间和单核转录组学同时捕捉大豆和ASR转录组

为表征大豆-ASR互作期间转录组的空间和时间变化，用ASR孢子感染大豆单叶。生成感染4天后大豆叶片的平切面（PD）切片和感染4天及7天后叶片的横切面（CS）切片。每种样本类型使用四个连续切片进行10x Genomics Visium实验。组织切片置于Visium载玻片上后，RNA透化到载玻片上，载玻片上的点含有带空间条形码、UMI和poly(dT)引物的引物，以捕获带poly(A)尾的转录本。空间条形码允许将测序读段追溯回Visium载玻片上的特定点，然后可叠加并与载玻片上切片的明场图像比较，获得基因表达的空间信息。poly(dT)引物能够捕获带poly(A)尾的RNA分子，这些分子在大豆和ASR核中都很普遍，允许在空间转录组学实验中同时捕获大豆和ASR转录组。

12个Visium切片捕获的总大豆基因数从一万九千到三万二千不等，ASR基因数从一千到六千不等，表明充分覆盖了大豆和ASR的整个转录组。所有12个Visium切片的主成分分析显示，4天和7天样本形成 distinct簇，表明ASR和大豆基因的表达在感染时间过程中发生变化。分析每个点的ASR读段百分比发现，7天横切面的ASR读段百分比显著高于4天平切面或4天横切面，表明随着疾病进展，ASR感染在大豆叶片中扩展。每个切片中高ASR存在点的空间映射显示ASR感染在大豆叶片中的空间异质性。为分析所有切片不同点之间的关系，合并12个切片的所有点并生成UMAP。合并不同切片数据后，UMAP显示来自不同切片的点良好整合。检查每个点的ASR读段百分比发现，高ASR读段百分比的点在UMAP上紧密定位，表明ASR感染区域具有相似的基因表达谱。这些结果表明空间转录组学可以同时捕获大豆和ASR转录组，并识别ASR活跃感染的特定空间区域。

虽然空间转录组学实验捕获了基因表达以及原位空间定位信息，但Visium载玻片上的点分辨率为55 μm，并非真正的单细胞分辨率。因此，使用10x Genomics Chromium平台用单核RNA测序补充空间转录组学研究。从感染1、3、5天后ASR感染大豆叶片分离核，每个时间点两个生物学重复，并进行单核RNA测序。六个单核RNA测序样本中，捕获的大豆细胞数从五千到一万一千不等，检测到的大豆总基因数从四万二千到四万六千不等，表明单核RNA测序数据捕获了大豆细胞的大部分全转录组。所有六个单核RNA测序样本的PCA分析显示所有生物学重复聚集在一起，表明单核RNA测序的可重复性。与Visium平台类似，Chromium平台也使用poly(dT)引物捕获转录本，允许同时捕获大豆和ASR转录本，六个单核RNA测序样本中捕获了三千到五千个ASR基因。单核RNA测序样本中每个核的ASR读段百分比从1天到3天以及从3天到5天显著增加，证实了ASR感染随时间进展。为探索六个单核RNA测序样本中所有核基因表达之间的关系，生成所有核的UMAP。合并数据后，来自不同样本的核均匀分布在UMAP上，表明良好整合，允许在特定核簇中检查不同样本间的基因表达变化。在单核RNA测序数据集中，一小簇ASR核显示每个核高ASR读段百分比，再次表明单核RNA测序实验成功捕获了大豆和ASR转录组。

这些结果表明空间转录组学和单核RNA测序都可以同时捕获宿主大豆和病原体ASR的转录组。空间转录组学数据允许可视化ASR感染中心，突出大豆叶片中ASR感染