大鼠脊髓转录组发育特征及时间模式研究

工程7（2021）1592研究组织工程-文章大鼠脊髓转录组发育的时间模式和分子网络特征Jian Yanga，b，#，Lili Zhaob，#，Sheng Yib，#，Fei Dingb，Yumin Yangb，Yan Liub，Yongjun Wangb，MeiLiub，Chengbin Xueb，Lian Xub，Leilei Gongb，Xinghui Wangb，Yu Zhangb，Bin Yub，Guo-li Mingc，Xiang，Gua，b，a苏州大学医学院人体解剖学与组织胚胎学系，苏州215123b南通大学附属医院江苏省组织工程与神经损伤修复临床医学中心，江苏省教育部神经再生重点实验室，南通大学神经再生联合创新中心，江苏南通226001c宾夕法尼亚大学佩雷尔曼医学院神经科学系，费城，PA 19104，美国阿提奇莱因福奥文章历史记录：收到2019年2020年3月15日修订2020年6月24日接受在线预订2021年保留字：发展基因表达模式先天免疫脊髓转录组A B S T R A C T脊髓发育的分子网络特征是组织工程不可或缺的，在胎盘哺乳动物中仍然知之甚少，特别是在它们与重要生物过程如再生的关系方面。在这里，使用大规模的时间转录组分析大鼠脊髓从胚胎阶段到成年，我们表明，波动的RNA表达水平反映了高度活跃的转录调控，这可能会启动脊髓图案。我们还表明，microRNA（miRNA）和转录因子表现出马赛克的基础上，他们的表达模式，而差异选择性剪接事件揭示，选择性剪接可能是一个驱动力的节点的发展兰维尔。我们的研究还支持先天免疫和内在生长能力之间存在负相关。表观遗传修饰在不同发育阶段表现出各自的调节功能，而鸟嘌呤核苷酸结合蛋白（G蛋白）偶联受体（包括嗅觉受体（ORs））在轴突生长中可能发挥多效性作用。这项研究为研究脊髓发育提供了宝贵的资源，并补充了越来越多的单细胞数据集。这些发现也为开发新的组织工程策略提供了遗传学基础©2021 THE COUNTORS.Elsevier LTD代表中国工程院出版，高等教育出版社有限公司。这是一篇CC BY-NC-ND许可下的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）中找到。1. 介绍研究胎盘哺乳动物脊髓发育的机制仍然很重要在组织工程、神经科学和再生医学领域[1]。在脊椎动物中，促进脊髓再生的内在生长能力在胚胎发育期间最大[2，3]，在幼年期变得较弱，在成年期几乎消失[4，5]。然而，从胚胎期到成年期脊髓发育过程中发生的转录变化，以及它们与脊髓内在生长能力的关系，仍然很差*通讯作者。电子邮件地址：gming@pennmedicine.upenn.edu（G.- L. Ming），nervegu@ntu.edu.cn（X. Gu）。#这些作者对这项工作做出了同样的明白此外，在发育中的大鼠脊髓的长期基因表达模式的系统研究尚未进行。为了在组织水平上追踪脊髓的长期发育过程，我们收集了从胚胎第9天（E9 d）到出生后第12周（P12 w）的大鼠的整个脊髓用于RNA测序（RNA-Seq）分析。从总共12个时间点收集36个样本以研究相关基因表达动态。2. 材料和方法2.1. 大鼠脊髓样品制备成年、妊娠和未成熟的无特定病原体（SPF）Sprague Dawley（SD）大鼠由中国南通大学实验动物中心提供大鼠麻醉，https://doi.org/10.1016/j.eng.2021.10.0012095-8099/©2021 THE COMEORS.由爱思唯尔有限公司代表中国工程院和高等教育出版社有限公司出版。这是一篇基于CC BY-NC-ND许可证的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）。可在ScienceDirect上获得目录列表工程杂志首页：www.elsevier.com/locate/engJ. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921593-≥2..乙醚吸入和断头器（Pushin Zr-11;上海普信仪器技术有限公司，有限公司、中国）在采样日06：30用冰块冷却制备。氧饱和、冰冷的人工脑脊液（ACSF）的制备如下组成成分：130 mmol/L氯化钠，5mmol·mL-1KCl、2 mmol·L-1KH2PO4、1.5 mmol·L-1CaCl2、6 mmol·L-1MgSO4、10 mmol·L-1 葡 萄 糖 ， 10 mmol·L-1HEPES （ pH 7.2 ， 渗 透 压 305mOsm）。断头后，快速分离大鼠脊柱，并用新鲜冰冷的ACSF洗涤两次 。 在 解 剖 显 微 镜 （ Olympus SZ 51; Olympus Corporation ，Japan）下，在新鲜、冰冷的ACSF中解剖整个脊髓。剥离硬脊膜、蛛网膜和软脊膜。在总共12个时间点（E9天、E11天、E14天、E18天、E19天）从E9天至P12周年龄的大鼠收集脊髓样品。P1d、P3d、P1w、P2w、P3w、P4w、P8w和P12w），各具有三个生物学重复。将所有组织在液氮中快速冷冻并储存在180 °C下。总共分离36个脊髓组织用于产生RNA-Seq数据。2.2. RNA提取和文库制备及测序使用TRIzol试剂（Ambion，Inc.，USA）根据制造商的方案。使用Epicentre Ribo-Zero rRNA Removal Kit（Epicentre，USA）去除总RNA的核糖体RNA（rRNA），并使用Illumina HiSeq 2000系统（Illumina，Inc.，美国）[6]。使用Lagos-Quintana等人的工作中提出的变性聚丙烯酰胺凝胶电泳（PAGE）方法分离小的非编码RNA（ 18-[7] 的 文 件 。 然 后 使 用 TruSeq 小 RNA 文 库 制 备 试 剂 盒（Illumina，Inc.）并在Illumina HiSeq 2000平台上测序。2.3. RNA-Seq分析在下游分析之前去除具有衔接子序列、高含量（10%）的未知碱基和低质量读数（读数中Q10的接下来，通过使用短寡核苷酸比对程序2（SOAP2）去除映射到rRNA序列的读段来过滤rRNA污染物[8]。大鼠参考基因组和注释文件（版本Rn6.0）从加州大学圣克鲁斯分校（UCSC）基因组浏览器数据库（genome.ucsc.edu）下载。为了定量基因表达，使用Bowtie 2将RNA-Seq读数映射到基因[9]，并通过期望最大化（RSEM，一种用于从RNA-Seq数据估计基因和同种型表达水平的软件包）使用RNA-Seq进行定量[10]。使用每百万映射读段的外显子的每个内切酶的片段（FPKM）来归一化基因长度和文库大小。NOIseq[11]用于鉴定两个连续发育期之间的差异表达基因（DEG）。将两个样本之间对数倍数变化≥1且分歧概率≥0.8的基因定义为DEG。2.4. 阶段特异性表达基因为了检测在给定阶段特异性表达的基因，使用Yanai等人[12]描述的方法。根据以下公式计算分级组织特异性指数s方程式：PN1-^xnxn¼;阶段n和N表示本研究中的发育时间点的数量。阶段特异性表达基因（SSEG）定义为s≥0.9的基因。2.5. 加权基因共表达网络分析为了鉴定脊髓发育过程中的遗传模式，我们使用WGCNA R软件包[13]进行了加权基因共表达网络分析（WGCNA），其中使用了E9d至P12 w样本中表达的所有基因的归一化表达（FPKM）选择软阈值来确定功率参数。将模块与来自用于检索相互作用基因/蛋白质（STRING）数据库的搜索工具的用R.2.6. 功能富集分析使用注释、可视化和综合发现数据库（DAVID）BioinformaticsResources 6.7[16，17]（http://david.abcc.ncifcrf.gov/）评估功能富集2.7. 可变剪接分析使用TopHat[18]鉴定剪接点，并使用转录物剪接的重复多变量分析（rMATS）[19]检测两个相邻样品之间的差异可变剪接事件（DASE）（例如，E9 d vs E11 d ，E11 d vs E14 d，E14 d vs E18 d）。rMATS计算包含和跳跃亚型的相对丰度，并计算差异剪接的P值和假发现率（FDR）。显著差异剪接定义为FDR≤ 0.05。2.8. 长链非编码RNA分析为了鉴定已知的和新的长非编码RNA（lncRNA），首先使用TopHat将每个样品的RNA-Seq读段映射到基因组，然后用Cufflinks组装[20]。在Cufflinks中实现的Cuffmerge和Cuffcompare用于合并来自所有样品的转录物，并且将组装的转录物与已知的信使RNA（mRNA）和lncRNA进行比较。保留了标记为“i“、"j“、" u“、"x“和" o“的转录本编码潜力计算器（CPC）[21]和iSeeRNA[22]用于预测新转录本编码潜力。 lncRNA的定量和差异表达分析与RNA-Seq分析相同（第2.3）。2.9. 小RNA分析首先去除衔接子序列和低质量读段。针对数据库映射干净序列以去除rRNA-、小条件RNA（scRNA）-、小核仁RNA（snoRNA）-、小核RNA（snRNA）-和转移RNA（tRNA）-相关读段。将剩余序列映射到miRbase （ v21 ） 和 Rfam （ v11 ） 数 据 库 ， 分 别 使 用 Bowtie2 和cmsearch[23]注释和分类为调节小RNA（microRNA（miRNAs）和Piwi相互作用RNA（piRNAs））。使用每百万转录本（TPM）计算miRNA表达水平[24]。差异使用DEGseq检测表达的小RNA[25]。miRNAs编号：1N-1其中^xn¼最大值xnð1Þ名为. 对数倍数变化。1和Q值0.001被定义为sig-1≤n≤N.二、 ≥≤基于所有脊髓的标准化表达数据发展时间点。xn表示基因在显著地差异表达。匠心途径分析(IPA)用于鉴定潜在的miRNA调节剂。J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）159215942.10. 实时定量聚合酶链反应（qRT-PCR）使用TRIzol试剂分离总RNA，然后使用PrimeScript RT ReagentKit（TaKaRa Bio，China）逆转录。为了定量基因表达，使用SYBRPremix Ex Taq （ Tli RNaseH Plus; TaKaRa Bio ， 中 国 ） 在StepOnePlus RT-PCR系统（Applied Biosystems，美国）。用2-DDCT法定量相对mRNA表达。GAPDH用作内部对照。所有实验重复三次。2.11. 免疫荧光将大鼠脊髓在4 °C下在4%多聚甲醛中固定过夜。将样品横切成12μm的厚度。用磷酸盐缓冲盐水（PBS，0.01mol·L-1，25°C; 21-040-CVR; Corning Inc.，USA）各5分钟，使用Immunol染色封闭缓冲液（Beyotime，中国）在25 °C下封闭切片1小时。将切片与一抗在4 °C下孵育12 -16小时。用PBS冲洗三次后，加入第二抗体，并将样品在25 °C下孵育1小时。然后将切片用Hoechst 33342染色15分钟。每个实验重复三次。使用以下抗体：抗Dnmt 3a（1：1000; ab 188470; Abcam，UK）、抗NEUR 0 D4（1：100; 14610-1-AP; Proteintech，USA）、抗巢蛋白（1：200;MAB 353; Merck Millipore ） 、 抗 NeuN （ 1 ： 200; MAB 377;Merck Millipore，USA）、Cy 3绵羊抗兔免疫球蛋白G（IgG）（1：500; C2306; Sigma-Aldrich，USA）、Alexa Fluor 488驴抗小鼠IgG（1：500; A21202; Invitrogen，USA）。3. 结果和讨论3.1. 大鼠脊髓转录组的景观我们从转录组中鉴定了15646个mRNA和667个miRNA。随后，我们使用qRT-PCR和免疫荧光（IF）随机筛选了两种基因和两种miRNA（Dnmt 3a、Neurod 4、miR-29 c-3 p和miR-18 a-3 p）用于验证（图1A和1B）。1（a）和（b））。我们还比较了IF数据与Allen Brain Atlas （ABA ）数据库中保存的小鼠脊髓原位杂交（ISH）数据[26]。总的来说，我们的IF结果与ISH数据一致（ABA数据库中只有两个点，出生后第4天和第56天），至少在出生后几天（P1w和P8w，在图1（a）中用星号标记）。这些结果证实了我们数据的质量，适用于下游分析。当共表达基因参与相似的生物学途径或受相似的调节途径调节时，它们倾向于聚集在一起[27]。WGCNA基于无标度网络和软阈值来计算基因的关系并检测共表达模块[13，28]。我们应用WGCNA来研究生物功能模块，并鉴定了总共21个模块（M1功能富集分析的模块确定了四个不同的阶段，在脊髓发育的基础上，生物过程的条款（图。 1（d）和（e））。第1阶段，即从E9d到E14d，包括从神经管形成原始脊髓。我们发现模块M6中的基因参与胚胎发育和转录调控，而M14中的基因在细胞命运决定中发挥作用（附录A中的表S1和S2）。第2阶段（从E14 d至P1 w）包括脊神经系统的初步发育，反映在与神经分化、增殖、营养和轴突发育相关的基因的富集（图1（e），附录A中的表S3第三阶段（从P1w至P8w）的特征在于参与髓鞘形成、先天免疫系统和氧化还原过程的基因富集，而第4阶段（P8w至P12w）的这些阶段性变化与相应时期的发展特征相一致[29]。此外，我们采用基因表达聚类分析来研究整体基因表达模式，并确定了两种表达模式（图1和图2中的模式1和2）。 2（a））。模式1中的基因在胚胎和新生儿发育过程中高度表达（7481个基因），并在胚胎发育的一些关键途径中富集（附录A中的表S8）。模式2中的基因从婴儿期到成年期高度表达（6926个基因），并参与能量代谢、免疫发育和细胞增殖（附录A）。SSEG的分析有助于更好地理解与特定发育阶段相关的生物学[12]。我们鉴定了55个在E9 d特异性表达的基因，这些基因与肠上皮细胞分化、前/后模式特化和神经板形态发生相关（图2（a）和附录A中的表S10）。在E11d发现了43个与转录调节、神经元生成和组织稳态有关的SSEG，而在E14d发现了9个与蛋白质自磷酸化和轴突细胞成分有关的SSEG有趣的是，在P8w和P12w中鉴定的SSEG都富集了与免疫系统相关的术语，如先天免疫应答抗原加工和呈递（表S10）。对于基因调控因子，我们发现miRNA的表达模式在出生前后呈现两种趋势（图2（b）中的miRNA谱）。我们还观察到两种转录因子（TF）相关的表达模式，以新生儿期（P3d）为边界线（TF谱见图11）。 2（b））。在P3d前，TF在与命运决定、干细胞功能调节和神经发育相关的方面富集（附录A中的表S11）。P3d后，与细胞分化和先天免疫激活相关的TF表达更高（表S11）。我们还比较了miRNA和TF的表达模式，以了解它们之间存在的关系。我们注意到，大多数TF（约70%）在P3 d前表达相对较高，而大多数miRNA（约60%）的表达水平相对较低;然而，这种模式在P3 d后逆转，TF的值约为30%，miRNA的值约为40%（我们称之为这一发现表明，宏观调控系统在脊髓中具有协调作用。我们还注意到，在E9d和E11d，干细胞群体维持可能是活跃的（附录A中的图S1）。在P3d，参与多巴胺能突触传递的基因似乎比参与胆碱能突触传递的基因更早表达（附录A中的图S2），这表明多巴胺能神经元可能在脊髓中更早出现。研究发育过程中的能量和物质代谢在这篇文章中，我们关注了参与能量储备代谢过程、自噬和氧化磷酸化的生物过程的基因的表达模式（图1A和1B）。2（c）和（d），附录A表S12和S13）。对于能源储备Meta-曲线过程（图中左侧） 2（c）），来自分支a和b的基因是富含糖原生物合成过程的术语正调节（P0.01;表S12）。来自分支b和c的基因富集了周期性调节的术语正调控。腺苷3 ′，5 ′-一磷酸（cAMP）生物合成过程（P <0.01;表S12）。对于自噬的积极调节J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921595图1.一、基因表达验证和共表达网络。（a）通过qRT-PCR（左）和免疫荧光染色（右）验证Dnmt 3a和Neurod 4的表达白色加号（+）表示阳性信号。P1w和P8w中的星号（*）表示我们的数据与Allen Brain Atlas（ABA）数据库中的原位杂交（ISH）数据一致。 比例尺：75μm。（b）通过qRT-PCR验证miR-29 c-3 p和miR-18 a-3 p的表达。rno：Rattusnorvegicus（挪威大鼠）。的x轴(a)（左）和（b）表示顺序发育阶段。红色和蓝色的y轴分别显示RNA-Seq的FPKM和qRT-PCR中的倍数变化(c)21个共表达模块的树状图。不同的颜色表示不同的模块。(d)几个重要的共同表达模块。Hub基因标记为红色。颜色图例是表达式的行Z分数（e）四个发展阶段及其特点。代表的共表达模块显示在括号中（底部）。J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921596图二、大鼠脊髓发育过程中的典型表达模式（a）大鼠脊髓中的mRNA表达模式上图显示SSEG的分布（s≥0.9）;下图显示两种mRNA表达模式。（b）miRNA和转录因子（TF）的嵌合表达模式左侧和右侧热图分别显示TF和miRNA的表达模式（c）参与能量代谢的基因的表达模式左侧热图显示了参与能量储备代谢过程的基因的表达模式;右侧热图显示了参与自噬正调控的基因的表达模式参与自噬的四个主要标志物基因的表达趋势（多项式拟合）显示在右上方。右侧y轴根据标记有星号（*）的基因的表达水平进行(d)氧化磷酸化相关基因的表达。a、b、c、d、e、f、g和h是分支。颜色图例是表达式的行Z分数。过程（图2（c）右侧）中，分支e的基因在对葡萄糖饥饿过程的细胞应答中（P0.01;表S13）和线粒体的正调控中显著富集。P3 d前的吞噬过程（P0.01;表S13）。P3d后，分支d显著富集与嗜异性正调控相关的基因（P0.01;表S13），这是一种复合物，天然免疫的本质和维持组织内环境稳定的功能。有趣的是，与氧化磷酸化电子传递链相关的基因的谱显示了成年期高表达的趋势（区块f;图2）。 2（d）和附录A中的表S14）。这种模式可能有助于为成年人提供更多的能量来源3.2. RNA表达通过比较相邻时间点（例如，E9 d vs E11 d，E11 d vs E14 d）。有趣的是，我们在E11d与E14d组中鉴定了2438个DEG（图14）。 3（a），Fig. 表S3和表S15附录A）。从数字上看，这些DEG约占三分之二的DEG被上调，并参与突触组装、轴突发生和神经系统发育（图2中的红色条）。所有下调的DEG均参与细胞分裂和细胞周期（图3（b）中的青色条）。此外，这两个相邻时间点之间差异表达的lncRNA和mRNA的分布显示出与DEG相似的模式（图3（a），红星）。因此，我们把这种现象称为为了进一步研究导致RNA表达波动的潜在调控因素，我们对参与基因表达转录后调控的miRNA进行了额外的分析。下调的miRNA数量显示出大幅增加，并且在E14天前占所有差异表达miRNA的76%（附录A中的表S16）。使用miRWalk2.0数据库的分析[30]显示，在J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921597图3.第三章。RNA波动和潜在的调控因子。（a）RNA波动。差异表达的mRNA和lncRNA的数量分别由红线和蓝线表示显著下调的piRNA和miRNA的百分比分别由橙色和绿色线表示（b）在E11d和E14d阶段之间差异表达的基因的基因本体红色和青色条分别表示上调和下调基因的富集。(c)上图显示了参与增强子结合的57个TF的表达模式的热图（颜色图例是表达的行Z分数）。下图显示了指示从热图中描绘的TF推断的潜在基因调控模型的示意图（d）有助于RNA波动和神经系统发育的潜在调节因素SWI/SNF：转换/蔗糖不可发酵; Wnt：无翅/int-1类。E11 d组与E14 d组相比，下调的miRNA显著富集于RNA聚合酶II启动子转录延伸的负调控（P0.01）、DNA模板转录的负调控（P0.01）和转录调控区DNA结合（P0.01）（附录A中表S17）。这一数字的大幅上升，修饰的miRNA可以降低mRNA降解的水平（图3（a））。已知piRNA通过Piwi-piRNA相互作用介导转录沉默我们发现下调的piRNA的比率从E9 d到E14 d增加（图3（a），橙色线;附录A中的表S18）。除了miRNAs和piRNAs，我们还鉴定了271个作用于DEG上游的TF。J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921598（附录A中的表S19），其中57个参与增强子结合的分子功能（附录A中的表S20）。表达聚类分析表明，在RNA表达的上述波动之前，57个TF中的约80%以高水平表达（图3（c）中的热图）。蛋白质编码基因通过相邻lncRNA的顺式调节是一种常见的过程[32]。我们进一步鉴定了在13个TF附近表达的29个lncRNA，其中27个显示与其相邻TF共表达（图S4和表S21附录A）。有趣的是，靶向这13个TF的miRNA主要以低水平表达。此外，这13个TF可以基于它们的分子功能分为两部分：A部分（例如，Ezh 2、Epc 1、Ruvbl 2、Smarca 4和Smarcc 1），其由参与表观遗传重编程和染色质重塑的TF组成;和B部分（例如，Tcf 3、Isl1、Tead 1、Tead2和Nkx 2 -5），其包含在细胞分化和神经发育的起始中起作用的TF（图3（d））。这些潜在的上游调节因子（miRNA、lncRNA和TF）可能协同促进所观察到的mRNA波动。这种波动可能参与脊髓的图案化和脊髓组织的区域化3.3. 选择性剪接事件从E9d到P12w，我们鉴定出5种主要的可变剪接形式（A3SS、A5SS观察到的DASE模式表明，大多数为SE事件（占所有事件的51.4%;附录A中的表S22）[19]。DASE的两个峰值与波动期和新生儿期一致（图4（a））。经历SE剪接事件的基因在生物学过程中富集，如神经元投射发育（P0.01，E11 d vs E14 d）和基于微管的运动（P0.05，E18 d vs P1d），而经历RI剪接事件的基因在生物学过程中富集，如动作电位期间的膜去极化（P0.01，E11 d vs E14 d）和核质的细胞成分（P0.01，E18 d vs P1d）。其他三个剪接事件（MEX、A5SS和A3SS）随时间推移表现出稳定的曲线（图4（a））。我们注意到，与Ranvier结、副阳极结和Anastaparanode结中主要成分相关的基因分布在DASE的两个峰值中（E11 d vs E14d和E18 d vs P1 d;图4（b）和附录A中的表S23）[33]。朗氏结及其周围结构域是有颌脊椎动物的一种进化修饰，选择性剪接在发育、分化和进化中发挥着不可或缺的作用[34，35]。因此，DASE的结果表明，选择性剪接可能是一个驱动力的发展的节点Ranvier。我们的研究结果还表明，具有新剪接变体的基因主要分布在E9d，P1d和P12w周围（图1）。 4（c））。我们以这三个时间点为参照，将发育过程分为以下阶段（图4（c））：①早期胚胎阶段（我们还对这些基因进行了WGCNA，并获得了27个共表达模块（附录A中的图S527个模块的功能逐渐从胚胎发育、干细胞维持和细胞分化转变为与髓鞘形成、免疫系统和维持稳态相关的功能（图4（d））。我们的数据提供了一个有价值的资源，为进一步调查新的剪接变异体。与选择性剪接相关的核糖核酸结合蛋白（RBP）也表现出类似的时间模式（图11）。附录A中的S6）。RBP通过相分离对温度变化或其他生物学线索做出反应，并从可溶性蛋白质转变为形成无膜细胞器的不溶性蛋白质[36]。24个RBP具有预测的内在无序区域（IDR，参与相分离[37];附录A中的图S7）在E9 d和E11 d高度表达（附录A中的图S8）。这些RBP可能参与脊髓对微环境信号的反应。CTCF是一种具有多种关键功能的蛋白质，包括选择性剪接和基因转录调控[38，39]，在E9d和E11d也高度表达（附录A）。3.4. 先天免疫力和内在生长能力我们发现，模块M7（图。 1（d））免疫应答基因显著富集（P0.01<）。M7中的基因表达在出生后逐渐上调，并在性成熟时达到峰值（图1（d）中的热图）。我们还发现，许多TF丰富的先天免疫相关的条款后P3d。接下来，由于小胶质细胞是中枢神经系统（CNS）中的主要先天免疫细胞[40]，我们分析了参与小胶质细胞增殖、发育和活化的基因的表达。我们发现这些基因的表达模式，包括小胶质细胞特异性标记物Cx3cr1的表达模式（附录A中的表S25）[41]，与模块M7中的基因相似（附录A中的图S9为了更好地了解脊髓的先天免疫和内在生长能力，我们研究了不同发育阶段可能涉及的因素，并建立了一个近似框架（图5（a）和附录A中的表S26）。 在胚胎期（E9 d-E18 d），与干细胞维持、细胞去分化相关的基因的高表达（图1A和1B）。S1和S10以及附录A中的表S26），未发育的先天免疫系统可能赋予脊髓强大的内在生长能力。然而，与这些过程相关的基因在婴儿期和青春期（P1 d-P4 w）下调（图1A和1B）。 S9和S10）。同时，在稳态（M0）至交替（M2）小胶质细胞极化（抑制细胞凋亡并促进再生[42]）中发挥作用的调节因子（包括TF、miRNA和lncRNA，如Kdm 6 b、miR-124和GAS 5）显示出促进M2小胶质细胞极化的趋势（表S26）。Toll样受体（TLR）水平的小幅增加有助于炎症控制和伤口愈合[43从P1d到P4w，大多数TLR的表达增加（附录中图S11A），这暗示了对内在生长能力的积极影响。因此，在此期间，由免疫系统相关的呼吸爆发（RB）引起的损伤可以被控制到一定程度，从而维持脊髓的内在生长能力。性成熟后（调节因子（如TLR、miR-21和miR-101）的表达变化可能促进M0至M1小胶质细胞极化（表S26）。白细胞介素4受体（IL4R）及其配体IL4将M0小胶质细胞重编程为M2小胶质细胞[42]。我们发现Il4r在P12w中显著下调（附录A中的表S27），这可能会损害重编程。模式识别受体，如TLR4和半胱天冬酶4（CASP4），促进M0至M1小胶质细胞极化[46]，在出生后表现出持续的高表达（附录A中的表S28），这可能导致脊髓中M1小胶质细胞的比例增加。成年哺乳动物在损伤后的CNS中发生大量小胶质细胞再增殖[47]，导致产生大量M1小胶质细胞，可诱导脊髓中的严重RB。此外，轴突修复过程中损伤相关信号逆行转运所必需的双亮氨酸拉链激酶（DLK）[48]在出生后显示出显著下降（附录A中的表S29），这可能会降低脊髓的内在生长能力。总之，这些数据表明，脊髓中复杂的相互作用可能支持J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921599图四、大鼠脊髓发育过程中选择性剪接事件的特定模式（a）大鼠脊髓发育期间的DASE（FDR≤ 0.05）。颜色表示不同的选择性剪接事件：SE（青色），RI（紫色），MXE（绿色），A3SS（深蓝色）和A5SS（红色）。（b）朗维尔结主要组成部分示意图编码这些组分的基因经历DASE，并在部分（a）的峰1和峰2中发现（c）具有新剪接变体的基因同种型的表达模式红星表示参考时间点。颜色图例是表达式的行Z分数。(d)四个阶段的代表性生物过程：早期胚胎阶段（E9 d-E11 d，蓝色）;晚期胚胎和新生儿阶段（E14 d-P1 w，红色）;幼年阶段（P2 w-P4 w，绿色）;和成年阶段（P8 w-P12 w，黄色）。cGMP：鸟苷3 β，5 β-环一磷酸; PKG：蛋白激酶G; CD：分化簇。J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921600图五.组织生长能力逐渐降低，表观遗传修饰和G蛋白偶联受体（GPCR）的特征。(a)示意图显示随着持续发育内在生长能力圆圈代表发育过程中影响内在生长能力的因素绿色和红色分别代表促进内在生长能力的有利和不利条件。内在生长能力随着持续发育而降低（下图）。(b)参与CpH甲基化的基因的表达模式。(c)参与N6-甲基腺苷（m6 A）甲基化的基因的表达趋势。颜色块表示参与不同过程的基因，即：m6A甲基化（绿色），反向m6A甲基化（红色），m6A甲基化的读者（蓝色）和相关的RBP（氖绿色）。（d）10 - 11易位（TET）蛋白家族的表达趋势上图显示Tet 1、Tet 2和Tet 3的表达趋势;下图显示五种TET调节因子（Lin 28 a、Vprbp、Ogt、Parp 1和Sin 3a）的表达趋势红线显示所分析基因的平均表达模式（e）全球项目完成情况报告的时间组群（f）具有嗅觉受体（OR）活性的基因的表达模式上图显示了421个OR的时间聚类;下图显示了27个相对高度表达的OR（FPKM≥0.3）的时间聚类和热图Th 1：辅助性T细胞 1型;REST：RE 1沉默转录因子。（b）、（c）和（f）中的颜色图例是表达式的行Z分数J. 扬湖，澳-地Zhao，S.Yi等人工程7（2021）15921601≥建立恶性循环，导致慢性和进行性神经退行性变，最终导致脊髓再生失败（图1）。 5（a））。3.5. 表观遗传修饰与G蛋白偶联受体表观遗传修饰是神经元发育的重要调控因子我们分析了脊髓发育过程中涉及三种表观遗传修饰的基因：CpH甲基化（非CpG，其中首先，CpH甲基化是CNS成熟的关键调节剂，并由DNA甲基转移酶3 α（DNMT3A）维持[49，50]。RE 1沉默转录因子（REST）可通过miR-29 c-3 p调控Dnmt 3a [51]。Dnmt3a在出生前表现出持续和缓慢下调的表达模式，这可能导致CpH甲基化在发育早期阶段积累（图5（b））[49]。其次，m6A甲基化与可变剪接和mRNA稳定性相关[52，53]。我们的数据显示，参与m6 A甲基化的大多数基因表现出缓慢下调的趋势（图5（c）），这表明m6 A甲基化可能在发育中达到了总体动态平衡[54第三，10 - 11易位（TET）蛋白家族的表达表明，活跃的DNA去甲基化事件可能发生在围产期左右（图1）。5（d）和附录A表S30）。五种主要TET调节剂（Lin 28 a、Vprbp、Ogt 、 Parp 1 和 Sin 3a ） 的 表 达 在 胚 胎 期 降 低 （ 图 57-61 ） 。 5（d））。这些数据表明，活跃的DNA去甲基化可能在围产期的基因调控中发挥关键作用。G蛋白偶联受体（GPCR）包括最大哺乳动物中的受体超家族[62]。共注释了883个具有GPCR活性的基因（GO：0004930）。基于GPCR表达的12个时间点的层次聚类分析显示两个时间聚类： 此外，GPCR 可以分为三个基因簇（图11）。附录A中的S12）。与嗅觉传导相关的通路（KEGG通路：rno 04740 ） 、 神 经 活 性显 著三 组 间 差 异 有 统 计 学 意 义（P0.01）。先前的研究表明，嗅觉受体（OR）在嗅上皮和精子中高度表达，而在其他组织中表达较低[63，64]。我们发现，在新生儿时期左右，421个具有OR活性的基因在脊髓中的表达发生了变化（图1）。 5（f），上图）。此外，27个相对高度表达的OR（FPKM0.3）表现出类似的变化（在图5（f）中由虚线表示，下图）。此外，对所有表达基因的WGCNA分析表明，在模块M6和M10中也发现了OR，并分别在胚胎发育和轴突生长中富集。ORs作为化学感受器，其低表达可能在脊髓发育中发挥多效性作用。4. 结论总之，我们已经产生了广泛的转录组数据，涵盖了脊髓发育期间的12个时间点，包括胚胎期、新生儿期、幼年期和成年期。据我们所知，这是最长期的，基于转录组的大鼠脊髓发育时程分析。通过生物信息学分析，我们探索了发育中脊髓中发生的分子事件，包括差异基因表达，共表达网络、阶段特异性表达模式和可变剪接事件。两个相邻时期的差异表达分析显示，mRNA表达的波动表明E14 d之前的高基因调控活性，E14 d是从神经管向脊髓发育过渡的关键时期[65]。我们发现与先天免疫系统相关的基因随着持续发育（特别是在出生后）而逐渐上调，这表明先天免疫与脊髓内在生长能力之间可能存在负相关[66]。我们的研究结果还表明，与兰氏结相关的基因经历了更多的DASE，这表明同种型多样性的增加可能有助于有颌脊椎动物兰氏结的进化外观[67，68]。本研究为研究脊髓发育过程中的转录变化提供了宝贵的资源。我们的时间转录组数据可以与其他研究产生的空间，单细胞数据相结合。本研究也为再生医学和脊髓相关组织工程技术的改进提供了一个整体的参考框架。致谢本研究得到了国家自然科学基金（31730031）、国家重点研发计划（2017YFA0104700和2016YFC1101603）、江苏省重点医学中心和江苏省高等学校重点学科发展项目（PAPD）的资助。数据和材料大 鼠 脊 髓 不 同 发 育 阶 段 的 RNA-Seq 和 miRNA-Seq 数 据 已 经 以BioProject登录标签PRJNA 505253保存在NCBI数据库中。遵守道德操守准则Jian Yang、Lili Zhao、Sheng Yi、Fei Ding、Yumin Yang、YanLiu、Yongjun Wang、Mei Liu、Chengbin Xue、Lian Xu、LeileiGong、Xinghui Wang、Yu Zhang、Bin Yu、Guo-li Ming和Guo-liGu声明他们没有利益冲突或财务冲突需要披露。附录A.补充数据本文的补充数据可在https://doi.org/10.1016/j.eng.2021.10.001上找到。引用[1] 维斯滕达尔河中国对组织工程的推动。Science 2012;338（6109）：900-2.[2] Yang P，Yang Z.增强内在生长能力促进成体CNS再生。神经科学杂志2012;312（1-2）：1-6。[3] Neumann S，Skinner K，Basbaum AI.维持成年神经元的内在生长能力促进脊髓再生。Proc Natl Acad Sci USA2005;102（46）：16848-52.[4] Schultz MB ， Sinclair DA. 当 干 细 胞 变 老 ： 干 细 胞 衰 老 的 表 型 和 机 制 。Development2016;143（1）：3-14.[5] Wells JM ， Watt FM. 哺 乳 动 物 器 官 内 源 性 再 生 和 修 复 的 多 种 机 制 。Nature2018;557（7705）：322-8.[6] Parkhomchuk D ，Borodina T， Amstislavskiy V， Banaru M， Hallen L，Krobitsch S，et al. Transcriptome analysis by strand-specific sequencing ofcomplementary DNA.Nucleic Acids Res 2009;37（18）：e123.[7] 李文辉， 李文辉，李文辉.鉴定编码小表达RNA 的新基因。Science2001;294（5543）：853-8.[8] Li R ， Yu C ， Li Y ， Lam TW ， Yiu SM ， Kristiansen