2025, Vol. 55
再生是指动物恢复失去或受伤的组织和器官的过程[1], 尽管强弱不同, 但大多数动物都具有再生能力。一些动物在受伤后, 其身体的组织碎片可以再生出完整的生物体, 这种能力被称为全身再生(Whole-body regeneration), 具备全身再生能力的动物能够重建起完整的身体构造(Body plan), 补充所有缺失的细胞类型, 并整合新、旧组织以重新启动所有的生物体功能[2]。相比之下, 其他具备再生能力的动物仅能再生出部分的器官、组织或细胞类型, 这种能力被称为结构再生(Structural regeneration)。
再生的分子和细胞机制目前已经在部分动物模型中得到研究, 揭示出再生过程的多样性[1]。部分脊椎动物依靠胚基具备了结构再生的能力, 从心脏、肢体和晶状体等器官到复杂的多组织结构, 如鳍、四肢和趾尖等在某些脊椎动物中均可再生[3]。虽然脊椎动物只表现出结构再生的能力, 但全身再生却在无脊椎动物中广泛存在[4]。属于扁形动物门(Platyhelminthes)的涡虫具备全身再生能力, 其身体被切成多段后, 每一段依然能再生出完整的体轴和器官[5]。同时, 全身再生的能力也在非双侧对称动物中广泛存在, 水螅(Hydra)属于刺胞动物(Cnidaria), 具备极强的再生能力[6]。对具有全身再生能力动物的研究表明, 虽然再生过程整体上呈现出多样性, 但是也存在部分相似之处, 例如, Wnt信号通路在上述物种组织正确再生过程中均起到了重要作用[4, 7-8]。
无肠目(Acoela)在过去的研究中常被归类为扁形动物门的一个目。与扁形动物门的涡虫类似, 无肠目动物也具有极强的全身再生能力。近年来对无肠目物种大分子数据的分析以及对形态特征的评估都表明, 无肠目不属于扁形动物门, 而是在进化上位于双侧对称动物的基部[9]。澳洲异尾涡虫(Heterochaerus australis)属无肠目(Acoela)盘旋科(Convolutidae)异尾属(Heterochaerus)。常与虫黄藻(Zooxanthella)和虫绿藻(Zoochlorellae)共生, 其解剖结构简单, 身体前端具有平衡囊、额腺和眼点等感觉器官, 没有消化系统与排泄系统, 雌雄同体, 无明显的生殖腺[10]。与其他无肠目动物类似, 澳洲异尾涡虫同样具有全身再生能力, 将其横切为头部和尾部两段组织后, 各自又能再生出完整的身体结构[10]。目前对无肠目动物再生过程的研究表明, Wnt与Bmp信号通路对其再生过程中体轴的正确建立起到了关键作用。Wnt信号控制其前后轴的再生, Bmp信号控制其背腹轴的再生, 这些通路的扰动会导致其再生过程中产生异常的表型[4]。同时, 在对尾部再生过程中Wnt-3基因及其他Wnt配体的研究中, 还建立了无肠目尾部组织再生启动的基因调控网络[11]。
无肠目动物基于其重要的进化地位以及全体再生的能力, 可能具备了双侧对称动物原始的再生调控网络。目前, 对于无肠目动物再生过程的研究主要集中于再生关键基因的功能及调控关系, 而再生作为一个复杂的生物过程往往基于细胞凋亡、增殖和分化等多种生物过程以及复杂的基因调控关系, 对再生过程中基因调控网络的解析有助于进一步理解复杂的再生过程[1]。基于此, 本研究对无肠目动物澳洲异尾涡虫头部组织与尾部组织多个再生时期进行了转录组测序, 使用加权基因共表达网络分析(Weighted gene co-expression network analysis, WGCNA)鉴定出头部组织与尾部组织再生过程中共有的关键基因共表达模块, 并进一步揭示了澳洲异尾涡虫再生过程中头部组织与尾部组织共有的基因调控网络的级联表达模式, 为研究双侧对称动物原始的再生调控网络提供参考。
1 材料与方法 1.1 澳洲异尾涡虫的培养及再生样品采集将本研究所用的澳洲异尾涡虫置于光照恒温培养箱中培养, 设置光暗周期为光照16 h∶黑暗8 h, 其中光照环境温度21 ℃, 黑暗环境温度19 ℃, 培养在盛有过滤灭菌海水的玻璃培养皿中, 每3 d投喂一次新鲜孵化的丰年虾幼虫, 最后一次喂食丰年虾幼虫3 d后, 开展切割再生实验。
挑选体积相似的成年澳洲异尾涡虫, 于身体中央部位横切, 以获得头部(前端)和尾部(后端)两部分组织, 分别放置在新培养皿中, 于光照恒温培养箱中培养。于再生0 h、3 h、6 h、12 h、1 d、2 d、4 d和7 d时, 分别取头部与尾部的再生组织各一块, 置于EP管中, 并放置于液氮中速冻, 将速冻组织保存于-80 ℃冰箱中等待送检。
1.2 RNA测序和质量控制将澳洲异尾涡虫各再生时期的速冻组织样品放于干冰中, 寄送至诺禾致源公司(天津)进行转录组测序。待样品RNA质检合格后, 制备测序文库。由于通过高通量测序得到的原始数据含有测序接头序列以及低质量序列, 为保证后续生物信息学分析所用数据的质量, 对原始数据进行过滤, 以去除被接头污染的碱基数大于5 bp的、低质量的(质量值Q≤19的碱基占总碱基的比例超过50%)及含N比例大于5%的Reads。最终得到质量较高的干净序列(Clean reads), 用于后续的数据分析。
1.3 转录组分析基于转录组测序质控后的Clean reads数据, 使用FastX_toolkit在线软件(https://github.com/agordon/fastx_toolkit)中fastq_quality_filter命令, 对预处理后的测序数据进行质控, 去掉含碱基N的基因序列, 并保留超过80%的碱基质量值大于20的高质量Reads。使用STAR在线软件(https://github.com/alexdobin/STAR)将Reads比对到经本课题组测序组装的澳洲异尾涡虫参考基因组上, 根据该基因组的GFF注释文件, 使用featureCounts命令对比对到澳洲异尾涡虫参考基因组上各基因的Reads数进行统计。由于转录本长度的差异会影响其表达量的计算, 因此需要计算基因长度并根据基因中所有Reads数目之和计算出每个基因的表达量。为减少不同样本的基因表达水平差异, 采用每千个碱基的转录每百万映射读取的转录本(Transcripts perkilobase million, TPM)作为归一化方法, 对表达量进行归一化处理, 得到头部组织与尾部组织各基因在不同再生时间点中的表达谱。TPM的计算方式如下: 将每个基因的表达量除以该样本中所有基因的表达量之和, 再乘以标准化因子100 000。该方法可以消除不同样本之间基因个数差异的影响, 使得各个基因的表达量可以直接进行比较[12]。
1.4 基因差异表达分析本研究以再生0 h的头部和尾部组织的表达谱数据作为对照组数据, 以头部组织与尾部组织再生3 h、6 h、12 h、1 d、2 d、4 d和7 d的表达谱数据作为实验组数据, 使用edgeR软件[13]分别计算头部组织和尾部组织再生各时期的差异表达基因, 设置|log2(Fold change)|≥2, 并且将阈值q＜0.05的基因作为筛选显著差异表达基因。
1.5 WGCNA基因共表达网络分析本研究使用WGCNA软件的一系列函数, 分别针对头部组织与尾部组织在8个再生时期中21 653个基因的表达数据, 构建各自再生过程的基因共表达网络。首先, 使用goodSamplesGenes函数进行异常值检测, 并确定适当的软阈值β(头部组织β=6, 尾部组织β=9), 以构建一个无尺度的网络。为清除背景噪音和伪关联带来的误差, 将邻接矩阵转为拓扑重叠矩阵(Topological overlap matrix, TOM), 同时利用函数dissTom=1－OM对拓扑重叠矩阵取逆, 得到相异性矩阵。其次, 利用函数hclust对相异矩阵进行层次聚类, 采用动态切割(Dynamic tree cut)法切割基因聚类树, 该过程可以把表达模式相近的基因合并在同一个分支, 每个分支代表一个共表达模块[14]。然后, 为识别头部和尾部组织中差异表达基因富集的模块, 使用超几何分布检验对每个模块进行差异表达基因的过表达分析, 并通过Benjamini-Hochberg法校正p值[15]。最后, 使用Cytoscape软件进行基因共表达网络的可视化[16]。
1.6 实时荧光定量PCR分析分别随机选取头部与尾部组织再生0 h、3 h、6 h、12 h、1 d、2 d、4 d和7 d的澳洲异尾涡虫各3只, 使用Trizol试剂(供应商: Invitrogen, CA, USA)提取各组织样品总RNA, 使用RNase-freeDNaseⅠ(供应商: TaKaRa, Dalian, China)处理样品以除去基因组DNA。采用微量分光光度计和1.5%琼脂糖凝胶电泳检测RNA的浓度和完整性。使用反转录M-MLV试剂盒(供应商: TaKaRa, Dalian, China)合成cDNA第一链, 反转录在PCR仪中进行, 具体程序如下: 30 ℃, 10 min; 72 ℃, 15 min; 冰上冷却。合成的cDNA置于-20 ℃保存备用。
根据基因差异表达分析结果, 在头部与尾部组织再生各时期的共有的差异表达基因中, 各选出一个显著差异表达基因, 以验证基因差异表达分析结果的准确性。同时根据基因共表达网络分析结果选出部分关键枢纽基因, 以验证基因的共表达关系。针对鉴定出的各基因序列, 利用Primer Premier 5软件设计荧光定量PCR引物(见表 1—3), 并选择Actin为内参基因。建立实时荧光定量PCR反应体系(总体系为20 μL): 1 μL的cDNA溶液(浓度20 ng/μL); 10 μL SYBY Color qPCR Master Mix; 正向和反向引物(浓度均为2 μmol/L)各4 μL。各样品分别设3个生物学重复以及3个技术重复。使用Roche Light Cycler 96实时荧光定量PCR仪进行扩增。设置反应条件: 94 ℃, 10 min; 94 ℃, 15 s; 60 ℃, 1 min, 40个循环; 94 ℃, 15 s, 1个循环。通过分析产物的溶解曲线检验引物的特异性, 采用2-ΔΔCt方法计算目的基因的相对表达量, 将数据结果用GraphPad Prism 8软件进行单因素方差分析及绘图, 所有数据均以平均值±均值标准误(Mean±SEM)来表示。
|
表 1 代表性差异表达基因荧光定量PCR引物 Table 1 Representative differentially expressed genes primers for qRT-PCR |
|
|
表 2 头部组织枢纽基因荧光定量PCR引物 Table 2 Hub genes primers for qRT-PCR of head tissue |
|
|
表 3 尾部组织枢纽基因荧光定量PCR引物 Table 3 Hub genes primers for qRT-PCR of tail tissue |
澳洲异尾涡虫完整个体的切割位置以及再生过程中形态特征变化如图 1所示: 再生0~12 h, 虫体形态特征变化不大; 再生2 d, 头部与尾部组织均再生出可辨认的新组织结构; 再生7 d, 头部与尾部组织基本完成再生, 从而具备了完整的虫体结构, 但新生的组织结构还需要进一步生长, 以达到与原组织类似的生长状态。
|
((a)完整虫体形态特征及切割位置; (b)头部组织各再生时期虫体形态特征; (c)尾部组织各再生时期虫体形态特征。(a) Morphological characteristics and cutting positions of the whole body; (b) Morphological characteristics of the body at each regeneration period of head tissue; (c) Morphological characteristics of the body at each regeneration period of tail tissue.) 图 1 澳洲异尾涡虫再生过程形态特征 Fig. 1 Morphological characteristics of H. australis during different stages of regeneration |
将共16个澳洲异尾涡虫头部组织与尾部组织再生过程样品的Clean Reads比对到本课题组测序组装的澳洲异尾涡虫参考基因组上, 经TPM归一化处理后, 计算出澳洲异尾涡虫头部组织与尾部组织再生过程的基因表达谱数据, 转录组数据质控与比对统计结果如表 4和5所示。
|
|
表 4 头部组织各样品转录组质控结果统计 Table 4 Statistical result of quality control for the samples of head tissue |
|
|
表 5 尾部组织各样品转录组质控结果统计 Table 5 Statistical result of quality control for the samples of tail tissue |
基于澳洲异尾涡虫头部组织与尾部组织各再生时期的表达谱数据计算各再生时期的差异表达基因。从头部组织再生3 h、6 h、12 h、1 d、2 d、4 d和7 d样品中分别筛选出2 044、1 683、2 452、674、2 013、2 638、4 121个差异表达基因, 总计再生过程中有8 721个差异表达基因, 其中, 再生1 d的差异表达基因数最少(共674个), 再生7 d时差异表达基因数最多(共4 121个)。尾部再生3 h、6 h、12 h、1 d、2 d、4 d和7 d组分别筛选出3 294、4 646、2 814、821、3 450、3 440、4 040和2 216个差异表达基因, 总计再生过程中共有10 686个差异表达基因, 其中, 再生1 d的差异表达基因数最少(共821个), 再生6 h时差异表达基因数最多(共4 646个)。
由于头部组织与尾部组织再生过程中需要重新建立的身体结构不同, 为找到在头部与尾部组织再生过程中均发挥作用的关键再生基因, 本研究建立了头部与尾部组织再生过程中共有的差异表达基因集, 该基因集中共6 696个基因, 其中可能包含澳洲异尾涡虫再生过程中关键基因。
为验证头部与尾部组织再生过程中共有的差异表达基因集的可信度, 从头部和尾部组织再生过程共有的差异表达基因中, 于再生3 h、6 h、12 h、1 d、2 d、4 d和7 d每个时期各选出一个显著差异表达基因(见表 1)进行实时荧光定量PCR(见图 2)。结果显示, 头部和尾部组织各再生时期选出的差异表达基因代表的表达量均显著高于该组织再生0 h时该基因的表达量, 且同一基因在头部组织中的相对表达量与在尾部组织中的相似, 从而验证了差异表达基因分析结果的准确性。
|
(***表示显著性P＜0.001；****表示显著性P＜0.000 1。所有数据均以Mean±SEM表示。***indicates significance P < 0.001; ****indicates significance P < 0.000 1 compared to the control group. The data are presented as mean±standard error of the mean.) 图 2 头部(a)和尾部(b)组织再生各时期共有差异表达基因相对表达量 Fig. 2 Relative mRNA expression level of representative differentially expressed genes both in the regeneration process of head(a) and tail(b) tissue |
本研究分别对澳洲异尾涡虫头部与尾部组织在再生0 h、3 h、6 h、12 h、1 d、2 d、4 d和7 d的样品基因表达谱数据，使用WGCNA软件建立各自再生过程的基因共表达网络(见图 3)，其中不同的颜色代表不同的基因共表达模块。头部组织的样品中鉴定出8个基因共表达模块，分别命名为H1—H8，各模块所含基因的数量从688到5 550不等；尾部组织的样品中同样鉴定出8个基因共表达模块，分别命名为T1—T8，各模块所含基因的数量从1 213到6 380不等。
|
图 3 头部(a)和尾部(b)组织再生过程基因共表达网络聚类树 Fig. 3 WGCNA cluster tree for the regeneration process of head(a) and tail(b) tissue |
头部与尾部组织再生过程中共有的差异表达基因在一定程度上代表了澳洲异尾涡虫再生过程中关键基因。为分别鉴定头部与尾部组织再生过程中的关键基因共表达模块，使用超几何检验，对头部与尾部组织的基因共表达模块进行了共有再生差异表达基因的富集分析，并通过Benjamini-Hochberg法调整p值，以进行多次检验校正[15](见图 3)。分析结果显示，H7模块与T8模块为头、尾共同差异表达基因所富集的模块。
2.5 再生关键模块鉴定对头部组织与尾部组织再生过程的基因共表达网络进行模块相似性比较，以鉴定澳洲异尾涡虫再生过程中关键的基因共表达模块(见图 4)。比较结果显示，富集头部与尾部组织再生共有差异表达基因的H7和T8模块同样也为头与尾部组织再生过程基因共表达网络中相似性最强的模块，有超过一半的高连通度基因在这两个模块中共同存在。作为头部与尾部组织再生过程中相似性最强的共表达模块，H7和T8模块中的大部分基因在头部和尾部组织再生过程中均存在共表达的调控关系，这暗示了这些基因所组成的基因调控网络在澳洲异尾涡虫头部与尾部组织再生过程中均起到了重要作用。
|
(表中数值代表两基因共表达模块间共有的基因数。The values in the table represents the number of genes shared between the two WGCNA modules.) 图 4 头部和尾部组织再生过程基因共表达模块相似性分析 Fig. 4 Similarity analysis of WGCNA modules for the regeneration process of head tissue and tail tissue |
结合差异表达基因富集分析(见图 3)与基因共表达模块相似性分析(见图 4)的结果，本研究以H7与T8模块为澳洲异尾涡虫再生关键模块进行后续分析。分别对H7与T8模块中的基因进行KEGG功能富集分析(见图 5)，结果表明，H7模块基因显著富集的信号通路(p＜0.05)与T8模块相似，均富集到了同细胞能量平衡维持功能有关的Ampk信号通路[17]、同细胞异常增殖抑制功能有关的p53信号通路[18]以及同DNA复制和细胞周期相关的信号通路。
|
图 5 T7(a)和H8(b)模块KEGG富集分析 Fig. 5 KEGG enrichment analysis of module H7(a) and T8(b) |
为进一步解析H7和T8模块的基因调控网络，分别基于头部组织的H7模块与尾部组织的T8模块，构建了各自的基因共表达网络(见图 6)。从H7模块的5 550个基因中进一步鉴定出7个子基因共表达模块，各模块的基因数从279到1 531不等。从T8模块的6 380个基因中进一步鉴定出8个子基因共表达模块，各模块的基因数从169到1 591不等。
|
图 6 头部组织H7(a)和尾部组织T8(b)模块基因共表达网络聚类树 Fig. 6 WGCNA cluster tree for the module H7 in head tissue(a) and T8 in tail tissue(b) |
依据基因注释以及表达模式，从每个子模块中选取连通度前25%的基因作为该模块的枢纽基因(Hub gene)进行分析。为直观展示这些基因的表达模式，绘制了基因表达量热图(见图 7和8)。鉴于头部组织的H7b、H7c、H7d模块和尾部组织的T8b、T8c模块中基因数量少、连通度普遍较小且没有明确的表达模式，故未对上述模块的基因进行后续分析。
|
图 7 H7子模块枢纽基因表达模式 Fig. 7 Expression pattern of picked hub genes in the submodule of H7 |
|
图 8 T8子模块枢纽基因表达模式 Fig. 8 Expression pattern of picked hub genes in the submodule of T8 |
头部组织H7a模块的枢纽基因大多在再生3 h时高表达，主要包括衔接因子相关蛋白复合体基因Adaptin、Map激酶相互作用的丝氨酸/苏氨酸激酶基因Mknk。H7f模块的枢纽基因主要在再生3和6 h时高表达，其中：再生3 h时高表达的基因有胞外信号调节激酶Erk、DNA复制的解旋酶/核酸酶Dna2、丝裂原活化蛋白激酶基因Mapk3、纤母细胞生长因子受体Fgfr3和生长激素调节蛋白基因Grtp1；再生6 h时高表达的基因有卷曲蛋白基因Frizzled7、跨膜受体蛋白基因Notch和衔接蛋白基因Sh2b1等。与H7f模块同样在再生3和6 h高表达的H7e模块主要的枢纽基因包括在再生3 h时高表达的β-连环蛋白基因β-catenin、骨形态发生蛋白基因Bmp2和细胞周期蛋白依赖性激酶基因Cdk1，以及在再生6 h时高表达的细胞周期素基因Cyclin-N、细胞周期素基因Cyclin-E、丝氨酸/苏氨酸蛋白激酶基因Pak4、F-盒蛋白基因F-box等。在上述模块之后，在再生12 h和1 d时高表达的H7d模块中的基因主要包括5-氧脯氨酸酶基因Oplah、pre-mRNA剪接因子基因Xab2、长链脂酰辅酶A合成酶基因Acsl4、α-酮戊二酸脱氢酶α亚基基因Bckdha和超氧化物歧化酶基因Sod等。
尾部组织再生过程中T8d、T8f和T8g模块最先在再生3 h时高表达，随后在再生6 h时T8e和T8h模块高表达，最后T8a模块中的基因在再生1 d时高表达。最先表达的T8d模块枢纽基因主要包括潜在转化生长因子结合蛋白基因Ltbp4、丝裂原活化蛋白激酶基因Map4k5和骨形态发生蛋白受体基因Bmpr2等。同样在再生3 h时高表达的T8g模块主要包括腺苷酸环化酶基因Ac9、介导TGF-β家族成员的信号转导基因Smad4。T8f模块主要包括早老素Presenilin2、跨膜受体蛋白Notch基因。在尾部组织再生6 h时高表达的T8e和T8h模块主要包括细胞周期蛋白基因Cyclin-A和Cyclin-H、细胞周期蛋白依赖性激酶基因Cdk1、转录因子基因E2f4、细胞分裂周期蛋白基因Cdc16以及胰蛋白酶基因Trypsin-1等。在再生1 d时高表达的T8a模块主要包括40S核糖体蛋白基因Rps2和Rps8、酰基转移酶基因Agpat2以及胃三酰甘油脂肪酶基因Lipf等。
将上述头部与尾部组织的子模块枢纽基因的共表达关系可视化(见图 9)，从而发现在澳洲异尾涡虫头部与尾部组织再生过程中：再生3 h时，伤口刺激再生相关的信号通路基因高表达，从而启动再生；再生3~6 h，同细胞周期和DNA复制相关的基因表达上调，从而促进细胞的分裂增殖；再生12 h~1 d，与代谢相关的基因表达上调，从而为新组织的再生提供能量。
|
(图中各节点代表H7和T8各子模块中所选枢纽基因。不同颜色代表不同高表达再生时期的基因：蓝色代表再生3 h高表达的基因；绿色代表再生3~6 h高表达的基因；橙色代表再生12 h~1 d高表达的基因。(a)头部组织H7子模块可视化基因共表达网络；(b) 尾部组织T8子模块可视化基因共表达网络。Each node in the figure represents the picked hub gene in the submodule of H7 and T8, and different colors represent genes at different regeneration periods of high expression level: Blue represents genes with high expression level at 3h regeneration; Green represents genes with high expression level at 3~6 h regeneration; Orange represents genes with high expression level at 12 h~1 d regeneration. (a)Network interactions in the submodules of H7; (b) Network interactions in the submodules of T8.) 图 9 H7、T8子模块枢纽基因可视化基因共表达网络 Fig. 9 Network visualization of hub genes interactions in submodules of H7 and T8 |
为验证上述枢纽基因在头部与尾部组织再生过程中是否存在共表达关系，分别从头部与尾部组织中选出了部分枢纽基因(见表 2和3)进行实时荧光定量PCR分析(见图 10和11)。结果显示，各时期选出的枢纽基因在再生过程中的表达模式符合本研究分析结果(见图 7和8)，即位于同一功能模块的枢纽基因具备相似的高表达时期，这进一步印证了基因共表达网络分析结果的准确性。
|
((a) 头部组织再生关键信号通路基因模块两个枢纽基因相对表达量；(b) 头部组织细胞周期相关基因模块两个枢纽基因相对表达量；(c) 头部组织代谢相关基因模块两个枢纽基因相对表达量。ns表示与对照组相比无显著性; *表示显著性P＜0.05;**表示显著性P＜0.01;*表示显著性P＜0.001;****表示显著性P＜0.000 1。所有数据均以标准差(Mean±SEM)表示。(a) Relative mRNA expression level of two hub genes from key signaling pathway genes module in head tissue regeneration; (b) Relative mRNA expression level of two hub genes from cell cycle related genes module in head tissue regeneration; (c) Relative mRNA expression level of two hub genes from metabolic related genes module in head tissue regeneration. ns indicates no significant difference; * indicates significance P < 0.05; ** indicates significance P < 0.01; *** indicates significance P < 0.001; **** indicates significance P < 0.000 1 compared to the control group. The data are presented as mean±standard error of the mean.) 图 10 头部组织H7模块枢纽基因mRNA相对表达量 Fig. 10 Relative mRNA expression level of hub genes of H7 module |
|
((a)尾部组织再生关键信号通路基因模块两个枢纽基因相对表达量；(b)尾部组织细胞周期相关基因模块两个枢纽基因相对表达量；(c)尾部组织代谢相关基因模块两个枢纽基因相对表达量。ns表示与对照组相比无显著性，*表示显著性P＜0.05，**表示显著性P＜0.01，***表示显著性P＜0.001，****表示显著性P＜0.000 1。所有数据均以标准差(Mean±SEM)表示。(a) Relative mRNA expression level of two hub genes from key signaling pathway genes module in tail tissue regeneration; (b) Relative mRNA expression level of two hub genes from Cell cycle related genes module in tail tissue regeneration; (c) Relative mRNA expression level of two hub genes from metabolic related genes module in tail tissue regeneration.ns indicates no significant difference; * indicates significance P < 0.05; ** indicates significance P < 0.01; *** indicates significance P < 0.001; **** indicates significance P < 0.000 1 compared to the control group. The data are presented as mean±standard error of the mean.) 图 11 尾部组织T8模块枢纽基因mRNA相对表达量 Fig. 11 Relative mRNA expression level of hub genes of T8 module |
澳洲异尾涡虫作为双侧对称动物，具有前后轴与背腹轴。当对其进行横切后，所产生的头部与尾部组织均能通过再生形成一个完整的个体。然而，由于再生过程中两者需要重建的身体部位不同，因此，两者在再生过程中形成了特异性的基因调控网络，以支持不同身体部位的再生。聚焦于上述特异性基因调控网络，目前已在无肠目动物以及其他具有全身再生能力的动物中展开了较为深入的研究。在对属于扁形动物门的地中海涡虫(Schmidtea mediterranea)[4]以及与澳洲异尾涡虫同属于无肠目动物的三带黑豹蠕虫(Hofstenia miamia)[11]再生过程的研究中发现，Wnt信号通路在尾部组织再生过程中起到了关键作用。在这两个物种再生过程中，抑制Wnt信号的表达将导致这两种动物原本应再生出尾部组织的部位再生出头部组织，从而形成双头的个体。相反，再生过程中Wnt信号的过表达则会导致两个物种再生出尾部组织，从而形成双尾的个体[4, 11]。
除特异性的基因调控网络外，头部与尾部组织再生过程中也必然存在共性的基因调控网络，这些基因在两者再生过程中均发挥作用，在一定程度上代表了澳洲异尾涡虫再生过程中的关键基因调控网络。本研究聚焦于澳洲异尾涡虫再生过程中的关键基因调控网络，对其再生过程中的头部与尾部组织进行了转录组测序，并通过WGCNA分析流程解析了澳洲异尾涡虫头部与尾部组织再生过程中共有的基因调控网络，确定了再生过程中的关键基因模块，这有助于挖掘澳洲异尾涡虫再生过程的潜在的分子机制。
本研究重点关注了头部与尾部组织再生过程中共有的差异表达基因，以减少头部与尾部组织再生过程中特异性的基因调控网络的影响。通过WGCNA分析流程，分别构建了头部与尾部组织再生过程中的基因共表达网络，鉴定出头部与尾部组织再生过程共有的差异表达基因集所富集的共表达模块H7和T8(见图 3)。同时，基因共表达模块相似性分析的结果表明，这两模块的相似性最强，即两模块的基因在头部与尾部组织再生过程中均存在共表达的调控关系(见图 4)，推测其在澳洲异尾涡虫再生过程中起到了关键作用。
在鉴定出澳洲异尾涡虫再生关键基因共表达模块H7和T8后，对两模块中的基因开展了KEGG功能富集分析，结果显示，H7和T8模块中的基因显著富集于Ampk信号通路、p53信号通路以及同DNA复制和细胞周期相关的信号通路中，这些信号通路的功能分别同细胞能量平衡维持、细胞异常增殖抑制以及细胞增殖有关[17-18]，从而共同参与维持澳洲异尾涡虫头部组织与尾部组织的正确再生。为能进一步解析H7和T8模块的基因共表达网络，分别对两者的基因共表达网络子模块进行鉴定，并分析这些子模块枢纽基因的表达模式及功能。
在头部组织H7的子模块中，枢纽基因主要在再生3 h~1 d高表达。头部组织再生3 h时H7a、H7e和H7f模块的枢纽基因高表达，这些基因主要参与了Mapk/Erk信号通路的调控。目前有研究表明Mapk/Erk信号通路参与了许多动物的再生过程，该信号通路在蝾螈(Cynops orientalis)视网膜[19]、斑马鱼(Danio rerio)心脏[20]和大鼠(Rattus norvegicus)肝脏[21]的再生过程中发挥着重要作用。同时，Mapk/Erk信号通路也参与了涡虫再生的一系列关键环节，在涡虫再生过程中，它参与了干细胞的创伤反应[22]、再生启动[23]、增殖[24]和分化[25]等。在对地中海涡虫的研究[26]中发现，该涡虫的伤口响应依赖于Erk信号通路的激活，并且这一通路在伤口修复中具有高度保守性。在该涡虫中Erk基因的表达在受伤部位的组织中，受伤约0.5 h后达到峰值，然后以约1 h的延迟在相邻的组织间传播。使用抑制剂进行的功能分析结果表明，Erk基因的表达与干细胞有丝分裂直接相关。同时，Erk基因还参与了伤口信号在伤口附近与伤口远端的组织中的传播，伤口远端组织的Erk信号的传播对于伤口反应至关重要，Erk信号激活的延迟或缺失会影响该涡虫的再生过程[26]。头部组织再生3 h时高表达的枢纽基因功能表明，澳洲异尾涡虫再生的关键通路在此时因伤口的诱导已经激活，从而开启与再生相关的一系列反应。
头部组织再生6 h时，H7e和H7f模块的枢纽基因高表达。其中Notch参与的Notch信号通路影响细胞增殖、分化和凋亡[27]，对涡虫的研究显示，该通路在涡虫神经再生过程中发挥重要作用[28]。细胞周期蛋白(Cyclin)是一类具有细胞周期特异性或时相性的蛋白质，它与周期素依赖性激酶(Cyclin-dependent kinase, Cdk)共同影响细胞周期的运行，不同的周期蛋白在细胞周期中表达的时期不同。以往研究证明，Cyclin-E在G1期过渡到S期中起重要作用，其通过与Cdk-2相互作用，在多种细胞谱系的发育过程中调控细胞的增殖、静止状态和细胞命运[29]。Cyclin-D同样在G1期到过渡S期中起重要作用，其通过与Cdk-4相互作用，在胚胎发生过程中调节体腔细胞谱系和肠细胞的增殖，同时对线虫的研究证明，其参与了性腺发生过程中的性别决定[30-31]。E3泛素蛋白连接酶Trim71能够与miRNA相互作用，从而促进胚胎干细胞增殖并维持其状态，以往对哺乳动物的研究表明，其与miRNA相互作用后能够结合Ago2蛋白形成复合物，抑制周期蛋白依赖激酶抑制因子1A(Cyclin dependent kinase inhibitor 1A, Cdkn1a)的转录，从而促进细胞G1期到S期过渡，加速胚胎干细胞的更新[32]。头部组织再生6 h时高表达的枢纽基因表明，此时有关细胞周期和细胞增殖的基因已经开始表达，这为新组织的补充创造条件。
头部组织再生12 h~1 d时，H7d模块的枢纽基因高表达，其中超氧化物歧化酶Sod能够催化超氧阴离子自由基歧化为过氧化氢和氧气，保护细胞免受活性氧中间体的损伤，从而维持细胞的正常功能和生物体代谢平衡。Oplah催化5-氧代-L-脯氨酸的裂解生成L-谷氨酸；Acsl4催化长链脂肪酸转化为其活性形式酰基辅酶A，从而使其可以进一步参与β-氧化[33-34]。该时期头部组织高表达的基因功能主要同抗氧化、蛋白质和脂质代谢有关，推测该时期再生组织的代谢强度增加，从而为新组织的再生过程提供能量。
在尾部组织T8的子模块中，枢纽基因主要在再生3 h~1 d高表达。T8d、T8f和T8g三个模块的枢纽基因在尾部组织再生3 h时高表达，这些基因主要参与了Wnt与Bmp信号通路的调控。目前对涡虫再生机制的研究表明，Wnt信号通路是涡虫尾部组织正确再生的关键信号通路，在再生过程中对建立正确的前后轴发挥着重要作用[4]。在涡虫再生过程中：抑制Wnt信号的表达会导致头部组织再生出新的头部组织，而不是尾部组织，从而产生双头的个体；相反，过表达Wnt信号则会导致尾部组织的错误再生，从而形成双尾的个体[11]。同时，在涡虫再生过程中，Bmp信号通路在建立正确的背腹轴方面发挥重要作用。在Bmp信号通路中，Bmp与Bmpr2结合，Bmpr2与Bmpr1结合并磷酸化Bmpr1，从而形成了1型和2型受体二聚体。磷酸化后的1型受体在C末端的两个丝氨酸残基上招募并磷酸化下游受体Smads1/5/8，磷酸化的Smad1/5/8(pSmad1/5/8)与Smad4结合形成三聚体复合物，并移位到细胞核。在细胞核内，该三聚体复合物直接调控下游靶基因的转录。在涡虫的再生过程中，Bmp在背侧表达，而干扰Bmp与Smad4的表达均会导致该涡虫再生出异常的体型[35]。尾部组织再生3 h，高表达的枢纽基因大多为澳洲异尾涡虫再生过程中关键信号通路的基因，但是具体的信号通路与头部组织该时期的具体通路有所不同。推测导致差异的原因如下：在涡虫再生过程中建立正确的体轴，头部与尾部组织在再生的起始阶段，它们的再生关键通路基因的表达模式便已产生了差异。
尾部组织再生6 h, T8e、T8f和T8h模块的枢纽基因高表达。有研究发现，Cyclin-H的表达和其表达蛋白的活性在整个细胞周期中保持稳定，它与Cdk7相互作用，参与细胞周期控制和RNA聚合酶2的转录[36]。Cyclin-A控制细胞周期G1期过渡到S期以及G2期过渡到M期，通过与Cdk1或Cdk2形成特异性丝氨酸/苏氨酸蛋白激酶全酶复合物来发挥作用[37]。转录因子E2f4通过E2识别位点与DP蛋白协同结合而形成Drtf1/E2f复合物，并在调控细胞周期从G1期向S期的过渡中发挥作用[38]。与头部组织的H7e和H7f模块类似，该时期高表达的枢纽基因的主要编码细胞周期蛋白以及同细胞周期和细胞增殖有关的蛋白，从而为新组织的再生创造条件。
尾部组织再生12 h~1 d，T8a模块的枢纽基因高表达。在甘油-3-磷酸生物作为前体参与甘油磷脂和甘油三酯的生物合成过程中，Agpat2通过在甘油骨架的sn-2位置加入酰基，从而将溶血磷脂酸转化为磷脂酸[39]。Lipf催化三酰甘油的水解，产生游离脂肪酸、二酰甘油、单酰甘油和甘油[40]。该时期尾部组织高表达的枢纽基因功能主要同蛋白质和脂质代谢有关，推测在该时期再生组织的代谢强度增加，从而为后续的再生过程提供能量。
总结H7和T8子模块枢纽基因表达模式和功能如下：澳洲异尾涡虫再生0 h~1 d，头部与尾部组织均呈现出再生关键信号通路、细胞周期、代谢与抗氧的枢纽基因级联表达，共同构成基因调控网络表达模式(见图 9)。
对于再生过程调控机制的解析一直以来都是科学领域的重要课题，目前在同样具有极强再生能力的涡虫研究中已取得较为领先的进展。在大多数动物中，损伤可以在远离伤口的未损伤组织中引起反应。在地中海涡虫再生过程中的研究发现，其纵向的体壁肌肉细胞通过Erk基因传播伤口信号以激活伤口远端组织的反应，并证明了远端组织反应的激活对于再生过程是必要的。在本研究中，Erk基因同样也是澳洲异尾涡虫再生3 h关键基因共表达网络中的一员(见图 9)，但其在澳洲异尾涡虫伤口远端组织反应激活中是否具有同样重要作用还需要进一步的研究。同时，大量的涡虫再生的研究成果为未来解析澳洲异尾涡虫再生调控网络提供了良好的基因参考。基于无肠目动物的进化地位，比较两者再生基因调控网络也能够为再生过程进化研究提供参考。
4 结语本研究基于澳洲异尾涡虫头部与尾部组织各再生时期的表达谱数据，通过基因差异表达分析，构建了头部与尾部组织再生过程共有的差异表达基因集，并使用WGCNA分析流程分别构建了头部与尾部组织再生过程中的基因共表达网络，还鉴定出澳洲异尾涡虫的再生关键基因共表达模块H7和T8。
通过功能富集分析发现，H7、T8模块中的基因显著富集到Ampk信号通路、p53信号通路以及同DNA复制和细胞周期相关的信号通路中，这些信号通路的功能分别与细胞能量平衡维持、细胞异常增殖抑制以及细胞增殖有关，从而共同参与维持澳洲异尾涡虫头部与尾部组织的正确再生。
同时，本研究基于H7和T8模块分别构建了两者的基因共表达网络，并针对两者子模块枢纽基因的功能及其在再生过程中表达模式进行了分析。结果显示，所有的子模块枢纽基因普遍在再生3 h—1 d期间高表达，同时在该时期内，头部与尾部组织均呈现出相似的基因调控网络表达模式，即再生关键信号通路、细胞周期、代谢相关的枢纽基因级联表达的基因调控网络表达模式。
综上所述，本研究基于WGCNA分析流程鉴定出在澳洲异尾涡虫的头部与尾部组织再生过程中同样起关键作用的基因共表达模块, 并进一步揭示了头部与尾部组织再生过程中共同存在的基因调控网络的级联表达模式。
| [1] |
Srivastava M. Beyond casual resemblance: Rigorous frameworks for comparing regeneration across species[J]. Annual Review of Cell and Developmental Biology, 2021, 37415-37440. (  0) |
| [2] |
Bely A E, Nyberg K G. Evolution of animal regeneration: Re-emergence of a field[J]. Trends in Ecology & Evolution, 2010, 25(3): 161-170. ( 0) |
| [3] |
Tsonis P A. Regeneration in vertebrates[J]. Developmental Biology, 2000, 221(2): 273-284. DOI:10.1006/dbio.2000.9667 ( 0) |
| [4] |
Srivastava M, Mazza-Curll K L, van Wolfswinkel J C, et al. Whole-body acoel regeneration is controlled by Wnt and Bmp-Admp signaling[J]. Current Biology, 2014, 24(10): 1107-1113. DOI:10.1016/j.cub.2014.03.042 ( 0) |
| [5] |
Ge X Y, Han X, Zhao Y L, et al. An insight into planarian regeneration[J]. Cell Proliferation, 2022, 55(9): e13276. DOI:10.1111/cpr.13276 ( 0) |
| [6] |
Hobmayer B, Rentzsch F, Kuhn K, et al. WNT signalling molecules act in axis formation in the diploblastic metazoan Hydra[J]. Nature, 2000, 407(6801): 186-189. DOI:10.1038/35025063 ( 0) |
| [7] |
Gurley K A, Rink J C, Sanchez Alvarado A. Beta-catenin defines head versus tail identity during planarian regeneration and homeostasis[J]. Science, 2008, 319(5861): 323-327. DOI:10.1126/science.1150029 ( 0) |
| [8] |
Lengfeld T, Watanabe H, Simakov O, et al. Multiple Wnts are involved in Hydra organizer formation and regeneration[J]. Developmental Biology, 2009, 330(1): 186-199. DOI:10.1016/j.ydbio.2009.02.004 ( 0) |
| [9] |
Philippe H, Brinkmann H, Copley R R, et al. Acoelomorph flatworms are deuterostomes related to Xenoturbella[J]. Nature, 2011, 470(7333): 255-258. DOI:10.1038/nature09676 ( 0) |
| [10] |
孙婷, 汪安泰. 澳洲异尾涡虫(无肠目, 盘旋科)新纪录及描述[J]. 动物学杂志, 2014, 49(3): 415-422. Sun T, Wang A T. Heterochaerus australis, a new record of turbellarian (Acoela, Convolutidae) from China[J]. Chinese Journal of Zoology, 2014, 49(3): 415-422. ( 0) |
| [11] |
Ramirez A N, Loubet-Senear K, Srivastava M. A regulatory program for initiation of Wnt signaling during posterior regeneration[J]. Cell Reports, 2020, 32(9): 108098. DOI:10.1016/j.celrep.2020.108098 ( 0) |
| [12] |
Li B, Dewey C N. RSEM: Accurate transcript quantification from RNA-Seq data with or without a reference genome[J]. BMC Bioinformatics, 2011, 12323. ( 0) |
| [13] |
Robinson M D, McCarthy D J, Smyth G K. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data[J]. Bioinformatics, 2010, 26(1): 139-140. DOI:10.1093/bioinformatics/btp616 ( 0) |
| [14] |
Langfelder Peter, Horvath Steve. WGCNA: An R package for weighted correlation network analysis[J]. BMC Bioinformatics, 2008, 9(1): 559. DOI:10.1186/1471-2105-9-559 ( 0) |
| [15] |
Benjamini Y, Hochberg Y. Controlling the false discovery rate: A practical and powerful approach to multiple testing[J]. Journal of the Royal Statistical Society, 1995, 57(1): 289-300. DOI:10.1111/j.2517-6161.1995.tb02031.x ( 0) |
| [16] |
Shannon P, Markiel A, Ozier O, et al. Cytoscape: A software environment for integrated models of biomolecular interaction networks[J]. Genome Research, 2003, 13(11): 2498-2504. DOI:10.1101/gr.1239303 ( 0) |
| [17] |
Garcia D, Shaw R J. AMPK: Mechanisms of cellular energy sensing and restoration of metabolic balance[J]. Molecular Cell, 2017, 66(6): 789-800. DOI:10.1016/j.molcel.2017.05.032 ( 0) |
| [18] |
Levine A J. p53: 800 million years of evolution and 40 years of discovery[J]. Nature Reviews Cancer, 2020, 20(8): 471-480. DOI:10.1038/s41568-020-0262-1 ( 0) |
| [19] |
Mizuno A, Yasumuro H, Yoshikawa T, et al. MEK-ERK signaling in adult newt retinal pigment epithelium cells is strengthened immediately after surgical induction of retinal regeneration[J]. Neuroscience Letters, 2012, 523(1): 39-44. DOI:10.1016/j.neulet.2012.06.037 ( 0) |
| [20] |
Liu P, Zhong T P. MAPK/ERK signalling is required for zebrafish cardiac regeneration[J]. Biotechnology Letters, 2017, 39(7): 1069-1077. DOI:10.1007/s10529-017-2327-0 ( 0) |
| [21] |
Araujo T G, de Oliveira A G, Tobar N, et al. Liver regeneration following partial hepatectomy is improved by enhancing the HGF/Met axis and Akt and Erk pathways after low-power laser irradiation in rats[J]. Lasers in Medical Science, 2013, 28(6): 1511-1517. DOI:10.1007/s10103-013-1264-y ( 0) |
| [22] |
Owlarn S, Klenner F, Schmidt D, et al. Generic wound signals initiate regeneration in missing-tissue contexts[J]. Nature Communications, 2017, 8(1): 2282. DOI:10.1038/s41467-017-02338-x ( 0) |
| [23] |
Yuan J, Wang Z, Zou D, et al. Expression profiling of planarians shed light on a dual role of programmed cell death during the regeneration[J]. Journal of Cellular Biochemistry, 2018, 119(7): 5875-5884. DOI:10.1002/jcb.26779 ( 0) |
| [24] |
Ermakov A M, Ermakova O N, Ermolaeva S A. Study of possible involvement of MEK mitogen-activated protein kinase and TGF-beta receptor in planarian regeneration processes using pharmacological inhibition analysis[J]. Ontogenez, 2014, 45(5): 355-360. DOI:10.7868/S0475145014050048 ( 0) |
| [25] |
Tasaki J, Shibata N, Nishimura O, et al. ERK signaling controls blastema cell differentiation during planarian regeneration[J]. Development, 2011, 138(12): 2417-2427. DOI:10.1242/dev.060764 ( 0) |
| [26] |
Fan Y, Chai C, Li P, et al. Ultrafast distant wound response is essential for whole-body regeneration[J]. Cell, 2023, 186(17): 3606-3618. DOI:10.1016/j.cell.2023.06.019 ( 0) |
| [27] |
Bolos V, Grego-Bessa J, de la Pompa J L. Notch signaling in development and cancer[J]. Endocrine Reviews, 2007, 28(3): 339-363. DOI:10.1210/er.2006-0046 ( 0) |
| [28] |
Dong Z, Huo J, Liang A, et al. Gamma-secretase inhibitor (DAPT), a potential therapeutic target drug, caused neurotoxicity in planarian regeneration by inhibiting Notch signaling pathway[J]. Science of the Total Environment, 2021, 781146735. ( 0) |
| [29] |
Richardson H E, O'Keefe L V, Reed S I, et al. A Drosophila G1-specific cyclin E homolog exhibits different modes of expression during embryogenesis[J]. Development, 1993, 119(3): 673-690. DOI:10.1242/dev.119.3.673 ( 0) |
| [30] |
Boxem M, van den Heuvel S. lin-35 Rb and cki-1 Cip/Kip cooperate in developmental regulation of G1 progression in C. elegans[J]. Development, 2001, 128(21): 4349-4359. DOI:10.1242/dev.128.21.4349 ( 0) |
| [31] |
Yanowitz J, Fire A. Cyclin D involvement demarcates a late transition in C. elegans embryogenesis[J]. Developmental Biology, 2005, 279(1): 244-251. DOI:10.1016/j.ydbio.2004.12.022 ( 0) |
| [32] |
Loedige I, Gaidatzis D, Sack R, et al. The mammalian TRIM-NHL protein TRIM71/LIN-41 is a repressor of mRNA function[J]. Nucleic Acids Research, 2013, 41(1): 518-532. DOI:10.1093/nar/gks1032 ( 0) |
| [33] |
Watanabe T, Abe K, Ishikawa H, et al. Bovine 5-oxo-L-prolinase: Simple assay method, purification, cDNA cloning, and detection of mRNA in the coronary artery[J]. Biol Pharm Bull, 2004, 27(3): 288-294. DOI:10.1248/bpb.27.288 ( 0) |
| [34] |
Schaffer J E, Lodish H F. Expression cloning and characterization of a novel adipocyte long chain fatty acid transport protein[J]. Cell, 1994, 79(3): 427-436. DOI:10.1016/0092-8674(94)90252-6 ( 0) |
| [35] |
Gavino M A, Reddien P W. A Bmp/Admp regulatory circuit controls maintenance and regeneration of dorsal-ventral polarity in planarians[J]. Current Biology, 2011, 21(4): 294-299. DOI:10.1016/j.cub.2011.01.017 ( 0) |
| [36] |
Tirode F, Busso D, Coin F, et al. Reconstitution of the transcription factor TFIIH: Assignment of functions for the three enzymatic subunits, XPB, XPD, and cdk7[J]. Molecular Cell, 1999, 3(1): 87-95. DOI:10.1016/S1097-2765(00)80177-X ( 0) |
| [37] |
Pagano M, Pepperkok R, Verde F, et al. Cyclin A is required at two points in the human cell cycle[J]. EMBO Journal, 1992, 11(3): 961-971. DOI:10.1002/j.1460-2075.1992.tb05135.x ( 0) |
| [38] |
Ma L, Quigley I, Omran H, et al. Multicilin drives centriole biogenesis via E2f proteins[J]. Genes & Development, 2014, 28(13): 1461-1471. ( 0) |
| [39] |
Haque W, Garg A, Agarwal A K. Enzymatic activity of naturally occurring 1-acylglycerol-3-phosphate-O-acyltransferase 2 mutants associated with congenital generalized lipodystrophy[J]. Biochemical And Biophysical Research Communications, 2005, 327(2): 446-453. DOI:10.1016/j.bbrc.2004.12.024 ( 0) |
| [40] |
Rogalska E, Ransac S, Verger R. Stereoselectivity of lipases. ⅡStereoselective hydrolysis of triglycerides by gastric and pancreatic lipases[J]. Journal of Biological Chemistry, 1990, 265(33): 20271-20276. DOI:10.1016/S0021-9258(17)30500-8 ( 0) |