Integrative single-cell RNA-seq and ATAC-seq analysis of myogenic differentiation in pig 猪肌原分化的单细胞 RNA-seq 和 ATAC-seq 整合分析
Shufang Cai , Bin Hu , Xiaoyu Wang , Tongni Liu , Zhuhu Lin , Xian Tong , Rong Xu , Meilin Chen ,Tianqi Duo , Qi Zhu , Ziyun Liang , Enru Li , Yaosheng Chen , Jianhao Li , Xiaohong Liu and Delin Mo
Abstract 摘要
Background Skeletal muscle development is a multistep process whose understanding is central in a broad range of fields and applications, from the potential medical value to human society, to its economic value associated with improvement of agricultural animals. Skeletal muscle initiates in the somites, with muscle precursor cells generated in the dermomyotome and dermomyotome-derived myotome before muscle differentiation ensues, a developmentally regulated process that is well characterized in model organisms. However, the regulation of skeletal muscle ontogeny during embryonic development remains poorly defined in farm animals, for instance in pig. Here, we profiled gene expression and chromatin accessibility in developing pig somites and myotomes at single-cell resolution. 背景骨骼肌的发育是一个多步骤的过程，对这一过程的了解在广泛的领域和应用中至关重要，从其对人类社会的潜在医学价值，到其与改良农业动物相关的经济价值，不一而足。骨骼肌始于体节，在肌肉分化之前，肌肉前体细胞在真皮肌节和真皮肌节衍生的肌节中生成，这一发育调控过程在模式生物中具有良好的特征。然而，在农场动物（如猪）中，胚胎发育过程中骨骼肌本体的调控仍不十分明确。在这里，我们以单细胞分辨率分析了发育中的猪体节和肌节的基因表达和染色质可及性。
Results We identified myogenic cells and other cell types and constructed a differentiation trajectory of pig skeletal muscle ontogeny. Along this trajectory, the dynamic changes in gene expression and chromatin accessibility coincided with the activities of distinct cell type-specific transcription factors. Some novel genes upregulated along the differentiation trajectory showed higher expression levels in muscular dystrophy mice than that in healthy mice, suggesting their involvement in myogenesis. Integrative analysis of chromatin accessibility, gene expression data, and in vitro experiments identified EGR1 and RHOB as critical regulators of pig embryonic myogenesis. 结果 我们鉴定了成肌细胞和其他细胞类型，并构建了猪骨骼肌本体发育的分化轨迹。沿着这一轨迹，基因表达和染色质可及性的动态变化与不同细胞类型特异性转录因子的活动相吻合。沿着分化轨迹上调的一些新基因在肌营养不良小鼠中的表达水平高于健康小鼠，这表明它们参与了肌的发生。对染色质可及性、基因表达数据和体外实验的综合分析发现，EGR1和RHOB是猪胚胎肌生成的关键调控因子。
Conclusions Collectively, our results enhance our understanding of the molecular and cellular dynamics in pig embryonic myogenesis and offer a high-quality resource for the further study of pig skeletal muscle development and human muscle disease. 结论 总之，我们的研究结果增进了我们对猪胚胎肌发生过程中分子和细胞动态的了解，为进一步研究猪骨骼肌发育和人类肌肉疾病提供了优质资源。
These authors jointly supervised this work: Jianhao Li, Xiaohong Liu, Delin 这些作者共同指导了这项工作：李建豪、刘晓红、德林Mo. Shufang Cai, Bin Hu, Xiaoyu Wang, Tongni Liu contributed equally to thiswork. 工作。 Correspondence: 通讯：Jianhao Li 李健豪jianhao63@sina.comXiaohong Liu 刘晓红liuxh8@mail.sysu.edu.cnDelin Mo 莫德林modelin@mail.sysu.edu.cn State Key Laboratory of Biocontrol, School of Life Sciences, Sun Yat-Sen 中山大学生命科学学院生物控制国家重点实验室University, Guangzhou 510006, Guangdong, China 中国广东广州 510006 广州大学 Guangdong Key Laboratory of Animal Breeding and Nutrition, State Key 广东省动物育种与营养重点实验室，国家重点实验室Laboratory of Livestock and Poultry Breeding, Institute of Animal Science, 动物科学研究所家畜家禽育种实验室、Guangdong Academy of Agricultural Sciences, Guangzhou 510640, 广东省农业科学院，广州 510640、Guangdong, China 中国广东
Faculty of Forestry, University of British Columbia, Vancouver, BC V6T 不列颠哥伦比亚大学林学院，不列颠哥伦比亚省温哥华，V6T
Background 背景介绍
Skeletal muscle development is a highly complex and tightly coordinated multistep process [1]. In mammals, all skeletal muscles of the body derive from the somites, which are metameric mesodermal structures that form on both sides of the neural tube. The most dorsal portion of the somite remains epithelial and becomes the dermomyotome [2]. Muscle precursor cells labeled by the paired box transcription factors Pax3 and Pax7 are generated in the dermomyotome and dermomyotome-derived myotome and undergo cell-fate commitment, followed by migration along an established pathway to differentiate into skeletal muscles [3]. In pigs, for example, the development of somites occurs between approximately days 14 and 22 of gestation [4]. Pluripotent mesodermal cells originate in somites and are committed to the myogenic lineage, followed by proliferation of myoblasts and a subsequent two waves of myoblast differentiation and fusion to form primary and secondary myofibers. Primary myofibers are generated at days of gestation, followed by the secondary myofibers that form based on the template of the primary myofibers at approximately days of gestation [5]. 骨骼肌的发育是一个高度复杂和紧密协调的多步骤过程[1]。在哺乳动物中，身体的所有骨骼肌都来自体节，体节是在神经管两侧形成的元中胚层结构。体节最背侧的部分仍为上皮细胞，成为皮肌层[2]。肌肉前体细胞由配对盒转录因子 Pax3 和 Pax7 标记，在真皮肌小体和真皮肌小体衍生的肌小体中生成，并进行细胞命运承诺，然后沿着既定途径迁移，分化成骨骼肌 [3]。以猪为例，体节的发育大约发生在妊娠期的第 14 到 22 天之间 [4]。多能中胚层细胞起源于体节，并致力于成肌细胞系，随后成肌细胞增殖，经过两波成肌细胞分化和融合形成初级和次级肌纤维。初生肌纤维在妊娠 天时生成，随后在妊娠约 天时以初生肌纤维为模板形成次生肌纤维[5]。
The myogenic process involves complex gene expression regulatory networks, which mainly exert their function through precise regulation of intercellular signaling and the control of specific gene expression. Pax3 and Pax7 are essential upstream regulators of myogenesis. Pax3/Pax7-positive cells in the dermomyotome provide the self-renewing reserve cell population for muscle formation [6]. Myogenic regulatory factors (MRFs, including Myf5, MyoD1, myogenin, and MRF4), as members of the basic helix-loop-helix family of transcription factors, have well known roles in controlling the determination and differentiation of skeletal myogenic cells during embryonic and postnatal myogenesis . In vivo, molecular and genetic experiments in mice, Drosophila, and chickens have uncovered the genetic and epigenetic regulatory mechanisms necessary for skeletal muscle formation . In vitro, the precise process of myogenesis has been extensively studied using the murine immortalized C2C12 myoblast cell line as a model system because of its high proliferation and differentiation capacities [11, 12]. However, similar progress in pig myogenesis has been limited by the lack of a convenient model system and the high heterogeneity of pig primary myogenic cells. 肌生成过程涉及复杂的基因表达调控网络，主要通过精确调控细胞间信号传导和控制特定基因表达来发挥其功能。Pax3和Pax7是肌形成过程中必不可少的上游调控因子。真皮肌元中的 Pax3/Pax7 阳性细胞为肌肉的形成提供了自我更新的后备细胞群 [6]。肌生成调节因子（MRFs，包括Myf5、MyoD1、myogenin和MRF4）是转录因子碱性螺旋环螺旋家族的成员，在胚胎和出生后的肌生成过程中控制骨骼肌生成细胞的决定和分化 的作用众所周知。在体内，小鼠、果蝇和鸡的分子和遗传实验揭示了骨骼肌形成所必需的遗传和表观遗传调控机制 。在体外，由于小鼠永生化 C2C12 肌母细胞系具有很强的增殖和分化能力，因此人们利用该细胞系作为模型系统，对肌肉发生的精确过程进行了广泛研究 [11，12]。然而，由于缺乏方便的模型系统以及猪原代肌细胞的高度异质性，猪肌生成的类似研究进展一直受到限制。
The skeletal muscle of agricultural animals (typically pig, cattle, and beef) is one of the most significant dietary protein sources for human consumption. The development and growth of skeletal muscle determine meat yield and quality [13]. Thus, a better understanding of porcine skeletal muscle development is needed because of the agricultural importance of pigs and increasing awareness of the benefits of using pigs as model organisms for human development and disease [14, 15]. Tibetan pig, an indigenous pig breed of China, used to live only in the plateau, but now they can also live normally in the plains [16, 17]. The meat yield of Tibetan pig is low, but the meat yield is significantly improved in its hybrid offspring with Duroc, with a higher preslaughter weight ( ) and a larger loin eye area ( ) than those of Tibetan pig (preslaughter weight, ; loin eye area, 16.83 ) [18]. Therefore, using Tibetan pig (ZZ) and Duroc Tibetan pig (DZ) as models to study the regulatory mechanism of skeletal muscle development would be very beneficial for improving pig farming. Whether there are differences in the early embryonic myogenesis of ZZ and DZ , with the same female parent and different male parent, remains to be studied. 农业动物（通常是猪、牛和牛肉）的骨骼肌是人类食用的最重要的膳食蛋白质来源之一。骨骼肌的发育和生长决定了肉的产量和质量 [13]。因此，由于猪在农业上的重要性，以及越来越多的人认识到将猪作为人类发育和疾病模型生物的益处，因此需要更好地了解猪骨骼肌的发育[14, 15]。藏香猪是中国土生土长的猪种，过去只生活在高原地区，现在也能在平原地区正常生活 [16，17]。藏香猪的产肉率较低，但其与杜洛克猪杂交后代的产肉率显著提高，屠宰前体重（ ）和腰眼面积（ ）均高于藏香猪（屠宰前体重， ；腰眼面积，16.83 ）[18]。因此，以藏香猪（ZZ）和杜洛克 藏香猪（DZ）为模型，研究骨骼肌发育的调控机制，对提高养猪业水平大有裨益。ZZ和DZ具有相同的雌性亲本和不同的雄性亲本，其早期胚胎肌肉发生是否存在差异仍有待研究。
Single-cell RNA sequencing (scRNA-seq) allows the comprehensively profiling of the gene expression changes observed in development at a cellular level [19], while single-cell transposase-accessible chromatin sequencing (scATAC-seq) allows the analysis of chromatin accessibility and transcription factor (TF) binding to be profiled at a similar resolution [20]. Recently, scRNA-seq has been applied to studying cellular heterogeneity in skeletal muscle tissues [21, 22]. Scientific researchers outlined the major mononuclear cell types present in mature skeletal muscle of mouse, ranging from 10 to 15 cell types depending on cluster assignments and the granularity of subtyping [21-23]. The major cell types always include the following broad categories: fibro/adipogenic progenitors (FAPs), tenocytes, endothelial cells, smooth muscle cells, immune cells (B cells, T cells, macrophages, neutrophils), neural/glial cells, and satellite cells [22, 23]. In addition, using scRNA-seq, Ziye Xu et al. uncovered cell and lipid dynamics of fat infiltration in skeletal muscle [24]. Haibin Xi et al. profiled human skeletal muscle tissues from embryonic, fetal, and postnatal stages and constructed a "roadmap" of human skeletal muscle ontogeny across development [25]. However, our knowledge about muscle ontogeny in pigs is limited. Here, to investigate the upstream regulatory networks in myogenesis that lead to establishment of myogenic lineage and subsequent differentiation, we performed scRNA-seq and scATACseq of pig somite and myotome cells from Tibetan pigs (ZZ) and Duroc Tibetan pigs (DZ) at several embryonic stages (E16, E18, E21, and E28). We produced a classification of developing myogenic cells, observed dynamic changes in gene expression and chromatin accessibility along the myogenic differentiation trajectory, and defined cell-type-specific regulatory networks. We also investigated key TFs and cell-cell interactions associated with 单细胞 RNA 测序（scRNA-seq）可在细胞水平上全面分析发育过程中观察到的基因表达变化[19]，而单细胞转座酶可进入染色质测序（scATAC-seq）可在类似分辨率下分析染色质可进入性和转录因子（TF）结合情况[20]。最近，scRNA-seq 被应用于研究骨骼肌组织中的细胞异质性[21, 22]。科学研究人员概述了存在于小鼠成熟骨骼肌中的主要单核细胞类型，根据聚类分配和亚型的粒度不同，细胞类型从 10 到 15 种不等 [21-23]。主要细胞类型始终包括以下几大类：纤维/脂肪生成祖细胞（FAPs）、腱细胞、内皮细胞、平滑肌细胞、免疫细胞（B 细胞、T 细胞、巨噬细胞、中性粒细胞）、神经/胶质细胞和卫星细胞 [22，23]。此外，Ziye Xu 等人利用 scRNA-seq 发现了骨骼肌中脂肪浸润的细胞和脂质动态变化 [24]。Haibin Xi 等人分析了胚胎、胎儿和出生后阶段的人类骨骼肌组织，构建了人类骨骼肌在整个发育过程中的本体 "路线图"[25]。然而，我们对猪肌肉本体发育的了解还很有限。在此，为了研究肌肉发生过程中导致成肌系建立和后续分化的上游调控网络，我们对藏香猪（ZZ）和杜洛克 藏香猪（DZ）在几个胚胎阶段（E16、E18、E21和E28）的体细胞和肌细胞进行了scRNA-seq和scATACseq分析。 我们对发育中的成肌细胞进行了分类，观察了成肌细胞分化轨迹上基因表达和染色质可及性的动态变化，并定义了细胞类型特异性调控网络。此外，我们还研究了与 "细胞-细胞相互作用 "相关的关键 TFs 和细胞-细胞相互作用。
embryonic myogenesis. Finally, the molecular and cellular impacts of early growth response 1 (EGR1) and Ras homolog family member B (RHOB) were substantiated by overexpressing the two genes in pig primary myogenic cells (PPMCs) and C2C12 myoblasts, which resulted in a promotion of myogenic differentiation. This extensive analysis enhances our understanding of the molecular and cellular dynamics in pig embryonic myogenesis and provides an invaluable resource for studying animal skeletal muscle ontogeny and human muscle diseases. 胚胎肌形成。最后，通过在猪原代肌原细胞（PPMCs）和C2C12肌母细胞中过表达早期生长应答1（EGR1）和Ras同源家族成员B（RHOB），证实了这两个基因对分子和细胞的影响，从而促进了肌原分化。这项广泛的分析增强了我们对猪胚胎肌生成过程中分子和细胞动态的了解，为研究动物骨骼肌本体和人类肌肉疾病提供了宝贵的资源。
Results 成果
scRNA-seq identified the major cell types in developing pig somites scRNA-seq鉴定了发育中猪体节的主要细胞类型
To gain a comprehensive view of the cell populations present during pig skeletal muscle ontogeny, we used scRNA-seq to evaluate somite and myotome tissues of ZZ and DZ from embryonic day (E) 16 to E28, which covers the transition from progenitors to myocytes. Single cells from 8 samples (E16-ZZ, E16-DZ, E18-ZZ, E18-DZ, E21-ZZ, E21-DZ, E28-ZZ, and E28-DZ) were processed for scRNA-seq using a Chromium system (10× Genomics) (Fig. 1A). Overall, 70,201 cells passed quality control (QC) with an average of 1892 genes per cell and 5276 unique molecular identifiers (UMIs) per cell (Additional file 1: Figure S1A, B). To exclude technical batch effects, the datasets from all samples and tissues were merged using autoencoders (AEs) and applied the batchbalanced k nearest neighbors (BBKNN) approach [26] to the latent space [27] (Additional file 1: Figure S1C). Dimensional reduction and unsupervised clustering for all 70,201 cells identified 31 distinct clusters (Additional file 1: Figure S1D). Based on differential expression analysis and the expression of selected marker genes from the literature, we manually annotated 12 distinct populations (Additional file 2: Table S1 and Fig. 1B, C), including mesenchymal cells, fibroblasts, epithelial cells, neural stem cells, myogenic progenitors/myoblasts, osteogenic cells, neurons, neurogliocytes, endothelial cells, myocytes, chondrocytes, and muscle cells. Bubble plots of marker gene expression demonstrated the accuracy of the cell annotations (Fig. 1D). Gene Ontology (GO) analysis of the differentially expressed genes (DEGs) for each cell type verified the characteristics and functions of different cells (Additional file 3: Table S2 and Additional file 1: Figure S2). We next quantified changes in the cell type percentages during somite development. As shown in Fig. 1E, the patterns of the cell populations vary considerably across development stages, with the proportions of low differentiation cells, such as epithelial cells, endothelial cells, and neural stem cells, decreased with development, while the percentages of highly differentiated cells, such as fibroblasts, osteogenic cells, and myocytes, increased gradually with development, suggesting a rapidly somite development during E16 and E28. 为了全面了解猪骨骼肌本体发育过程中出现的细胞群，我们使用 scRNA-seq 评估了胚胎第 (E) 16 天至 E28 天期间 ZZ 和 DZ 的体节和肌体组织，其中涵盖了从 祖细胞到肌细胞的转变过程。利用 Chromium 系统（10× Genomics）对 8 个样本（E16-ZZ、E16-DZ、E18-ZZ、E18-DZ、E21-ZZ、E21-DZ、E28-ZZ 和 E28-DZ）的单细胞进行了 scRNA-seq 处理（图 1A）。总体而言，70201 个细胞通过了质量控制（QC），平均每个细胞有 1892 个基因，每个细胞有 5276 个唯一分子标识符（UMI）（附加文件 1：图 S1A、B）。为排除技术批次效应，所有样本和组织的数据集均使用自动编码器（AE）进行合并，并将批次平衡 k 近邻（BBKNN）方法[26]应用于潜空间[27]（附加文件 1：图 S1C）。对所有 70 201 个细胞进行降维和无监督聚类，确定了 31 个不同的聚类（附加文件 1：图 S1D）。根据差异表达分析和文献中选定标记基因的表达，我们人工标注了 12 个不同的群体（附加文件 2：表 S1 和图 1B、C），包括间充质细胞、成纤维细胞、上皮细胞、神经干细胞、成肌原/成肌细胞、成骨细胞、神经元、神经胶质细胞、内皮细胞、肌细胞、软骨细胞和肌肉细胞。标记基因表达的气泡图显示了细胞注释的准确性（图 1D）。对每种细胞类型的差异表达基因（DEGs）进行的基因本体（GO）分析验证了不同细胞的特征和功能（附加文件 3：表 S2 和附加文件 1：图 S2）。 接下来，我们对体节发育过程中细胞类型百分比的变化进行了量化。如图 1E 所示，不同发育阶段的细胞群模式差异很大，低分化细胞（如上皮细胞、内皮细胞和神经干细胞）的比例随着发育而降低，而高分化细胞（如成纤维细胞、成骨细胞和肌细胞）的比例则随着发育而逐渐升高，这表明 E16 和 E28 期间体节发育迅速。
To further dissect the cellular heterogeneity and transcriptional landscape of developing myogenic cells, we extracted myogenic progenitors/myoblasts, myocytes, and muscle cells for further clustering. The myogenic cells were further divided into 8 distinct sub-clusters with increased resolution and annotated as progenitors, myogenic progenitors, myoblasts, myocytes, cardiac muscle cells, and other cells (Fig. 1F, G). Among them, progenitors were characterized by the highest expression level of Pax3, whereas myogenic progenitors and myoblasts were characterized by the expression of the muscle stem cell marker Pax7 and the myogenic regulatory factor . Both myocytes and cardiomyocytes expressed the muscle cell marker ACTC1, but the differentially expressed skeletal myogenic cell markers MYOG and MYL1 as well as the cardiomyocyte-specific marker gene MYBPC3 [28] accurately distinguished the two cell types (Fig. 1H). 为了进一步剖析发育中的成肌细胞的细胞异质性和转录格局，我们提取了成肌细胞的祖细胞/成肌细胞、肌细胞和肌肉细胞进行进一步聚类。随着分辨率的提高，肌原细胞被进一步划分为8个不同的亚群，并分别注释为 祖细胞、肌原细胞、肌母细胞、肌细胞、心肌细胞和其他细胞（图1F、G）。其中， 祖细胞的特征是Pax3的表达水平最高，而肌原性祖细胞和肌母细胞的特征是肌肉干细胞标志物Pax7和肌原性调节因子 的表达。肌细胞和心肌细胞都表达肌肉细胞标记 ACTC1，但不同表达的骨骼肌细胞标记 MYOG 和 MYL1 以及心肌细胞特异性标记基因 MYBPC3 [28]能准确区分这两种细胞类型（图 1H）。
Reconstruction of the myogenic differentiation trajectory of progenitors 重建 祖细胞的成肌分化轨迹
To investigate the molecular processes underlying skeletal muscle development, the cells were ordered in a pseudotime manner using Monocle 2 [29]. Pseudotime trajectory analysis revealed seven different cell states 为了研究骨骼肌发育的分子过程，使用 Monocle 2 [29]对细胞进行了伪时间排序。伪时间轨迹分析显示了七种不同的细胞状态
Fig. 1 (See legend on previous page.) 图 1 （见前页图例。）
(states 1 7) and presented the distributions of cell states along with pseudotime flows (Fig. 2A). An organized, branched progression of cells from progenitors to differentiated myocytes was shown by labeling individual cells using the cell population annotations from the unified atlas in Fig. 1E (Fig. 2B, C). Unexpectedly, the first small branch (state 2) was almost entirely enriched by myocytes, and the cells at the end of pseudotime trajectory (state 7) also belonged to myocytes (Fig. 2A, C). The differential gene expression analysis indicated that many muscle development-related genes were highly expressed in both states of myocytes (e.g., MYMK, FNDC5, MEF2C, and TNNI1). However, the myocytes in state 2 expressed the cardiomyocyte-specific marker MYL9, while those in state 7 expressed the skeletal muscle cell-specific markers MYOD1 and MYOG (Fig. 2D). In addition, the DEGs with high levels in state 2 were largely involved in the regulation of biological processes such as "heart process" and "cardiac cell development," while the DEGs with highly levels in state 7 were involved in "skeletal muscle tissue development" and "skeletal muscle cell differentiation" (Fig. 2E). These results indicated that the myocytes in state 2 were actually cardiomyocytes which were excluded from subsequent analysis (Fig. 2F). (状态 1 7），并将细胞状态的分布与伪时间流一起显示出来（图 2A）。通过使用图 1E 中统一图谱的细胞群注释对单个细胞进行标记，显示了细胞从 祖细胞到分化的肌细胞有组织、有分支的发展过程（图 2B, C）。意想不到的是，第一个小分支（状态 2）几乎全部是肌细胞，伪时间轨迹末端（状态 7）的细胞也属于肌细胞（图 2A，C）。差异基因表达分析表明，许多肌肉发育相关基因在两种状态的肌细胞中都有高表达（如 MYMK、FNDC5、MEF2C 和 TNNI1）。然而，状态 2 的肌细胞表达心肌细胞特异性标记物 MYL9，而状态 7 的肌细胞表达骨骼肌细胞特异性标记物 MYOD1 和 MYOG（图 2D）。此外，状态 2 中高水平的 DEGs 主要参与调控 "心脏过程 "和 "心肌细胞发育 "等生物过程，而状态 7 中高水平的 DEGs 则参与 "骨骼肌组织发育 "和 "骨骼肌细胞分化"（图 2E）。这些结果表明，状态 2 中的肌细胞实际上是心肌细胞，因此被排除在后续分析之外（图 2F）。
Then, the percentages of progenitors, myogenic progenitors, myoblasts, and myocytes in each sample were calculated to assess myogenesis progression at different stages. Although the proportion of total myogenic cells in all cells of somite tissue did not change obviously (Fig. 1E), the percentages of four types of myogenic cells in different developmental stages changed significantly during E16 E28. At E16, there were almost no myoblasts and myocytes in the somites, with pax progenitors accounting for about of the four types of cells. At E18, some progenitors were committed to become myoblasts and further differentiated into myocytes. Subsequently, the proportion of pax progenitors decreased markedly, and myogenic progenitors, myoblasts, and myocytes accounted for about of the four types of cells at E21 and E28 (Fig. 2G). Dual immunostaining of Pax7 and MyoD showed that MyoD cells appeared at 然后，计算每个样本中 祖细胞、肌原细胞、肌母细胞和肌细胞的百分比，以评估不同阶段的肌生成进展。虽然体节组织所有细胞中肌原细胞总数的比例没有明显变化（图1E），但在E16 E28期间，四种类型的肌原细胞在不同发育阶段的比例发生了显著变化。E16时，体节中几乎没有肌母细胞和肌细胞，pax 祖细胞约占四种细胞的 。在E18时，一些祖细胞被承诺成为肌母细胞，并进一步分化为肌细胞。随后，pax 祖细胞的比例明显下降，在E21和E28时，成肌祖细胞、肌母细胞和肌细胞约占四种类型细胞的 （图2G）。Pax7和MyoD的双重免疫染色显示，MyoD 细胞在E21和E28出现。
E21 and increased significantly at E28, indicating that E18 28 was a critical period for the establishment of skeletal muscle lineage (Fig. 2H). As expected, cells from E16 and E18 embryos tended to be distributed at the root of the trajectory suggesting versatile progenitor properties, whereas those from E21 and E28 embryos were distributed in the later part of the trajectory indicating the decreased proportion of progenitor cells and increased proportion of differentiating myogenic cells during skeletal muscle development (Fig. 2B). E21和E28胚胎的细胞数量明显增加，表明E18-28是骨骼肌系建立的关键时期（图2H）。正如预期的那样，来自 E16 和 E18 胚胎的细胞往往分布在轨迹的根部，这表明它们具有多功能的祖细胞特性，而来自 E21 和 E28 胚胎的细胞则分布在轨迹的后期，这表明在骨骼肌发育过程中祖细胞的比例下降，而分化的成肌细胞比例上升（图 2B）。
In addition to differentiation, proliferation is also a major biological event during early skeletal muscle development. The cell cycle phase of each cell was predicted to evaluate the proportion of proliferating cells at different states of myogenic progression. We observed a shift in the transcriptomically defined cell cycle state accompanying the change in cell type representation (Fig. 2C and Additional file 1: Figure S3A). In the early part (states 1, 3, and 4) of the pseudotime trajectory, most of the cells are progenitors, so they are primarily predicted to be proliferating (S and G2M phases), with only a small fraction being predicted to be in G1 phase (non-proliferating). However, in the middle stage of the trajectory (State 5), proliferating cells decreased to , while non-proliferating cells increased significantly to . At the end of the trajectory (States 6 and 7), most of the cells have differentiated into myocytes that no longer proliferate, so they are primarily in G1 phase (75.3-85.1%) (Additional file 1: Figure S3B). Consistently, the variation in the expression of cell cycle-related genes and myogenic genes suggested that proliferating cells decreased with the differentiation trajectory (Additional file 1: Figure S3C, D). 除了分化，增殖也是骨骼肌早期发育过程中的一个主要生物学事件。我们预测了每个细胞的细胞周期阶段，以评估增殖细胞在成肌过程的不同状态下所占的比例。我们观察到转录组定义的细胞周期状态伴随着细胞类型代表的变化而发生了转变（图 2C 和附加文件 1：图 S3A）。在伪时间轨迹的早期（状态 1、3 和 4），大部分细胞是祖细胞，因此它们主要处于增殖期（S 期和 G2M 期），只有一小部分细胞处于 G1 期（非增殖期）。然而，在轨迹的中间阶段（状态 5），增殖细胞减少到 ，而非增殖细胞显著增加到 。在轨迹末期（状态 6 和 7），大部分细胞已分化成肌细胞，不再增殖，因此主要处于 G1 期（75.3-85.1%）（附加文件 1：图 S3B）。同样，细胞周期相关基因和成肌基因的表达变化表明，增殖细胞随着分化轨迹而减少（附加文件 1：图 S3C、D）。
Transcriptome dynamics of progenitor differentiation 祖细胞分化的转录组动态变化
To gain insights into the gene expression dynamics along the trajectory, the expression changes of the 1700 top DEGs among the four types of myogenic cells (including progenitors, myogenic progenitors, myoblasts, and myocytes) were analyzed and clustered into five major categories of transcriptional gene clusters (Additional 为了深入了解基因表达沿轨迹的动态变化，我们分析了1700个顶级DEGs在四种类型的肌原细胞（包括 祖细胞、肌原祖细胞、肌母细胞和肌细胞）中的表达变化，并将其聚类为五大类转录基因簇（附加数据
Fig. 2 (See legend on previous page.) 图 2 （见前页图例。）
file 4: Table S3 and Fig. 3A). Genes in cluster 1, with the highest expression in progenitors, were gradually downregulated from the beginning of programming and were largely involved in the regulation of biological processes such as "ATP metabolic process" (e.g., , ATIC, and MEIS1) (Fig. 3B). Subsequently, gene clusters 2 and 3 were transiently upregulated but finally downregulated, representing two temporary transcriptional waves. 文件4：表S3和图3A）。在 祖细胞中表达量最高的基因簇 1 中的基因从编程开始就逐渐下调，主要参与调控 "ATP 代谢过程 "等生物过程（如 、ATIC 和 MEIS1）（图 3B）。随后，基因簇 2 和 3 短暂上调，但最终下调，代表了两个暂时的转录波。
Cluster 2 genes highly expressed in myogenic progenitors were involved in "mitotic cell cycle process" (e.g., PCLAF, CENPA, and HMGB2) indicating the strongest proliferation capacity of myogenic progenitors during myogenesis. Cluster 3 genes highly expressed in myogenic progenitors and myoblasts are involved in "extracellular matrix" and "response to growth factor" (e.g., HES1, HEY1, and SIX1) (Fig. 3B). Concurrent with cluster 3 activation, cluster 4 genes highly expressed in myoblasts and myocytes were upregulated and maintained at high expression levels until the final stage with enrichment of the GO terms "response to endogenous stimulus" and "striated muscle tissue development" (e.g., TCF12, FOS, and EGR1). Finally, cluster 5 genes were activated at the end of the trajectory with predominant involvement in the GO term "muscle cell differentiation" (e.g., , and SOX6) (Fig. 3B). These data illustrated the trajectory of myogenic differentiation and revealed the ordered activation of transcriptional waves throughout this process. 在肌原纤维祖细胞中高表达的第 2 组基因参与了 "有丝分裂细胞周期过程"（如 PCLAF、CENPA 和 HMGB2），表明肌原纤维祖细胞在肌形成过程中的增殖能力最强。在成肌祖细胞和成肌细胞中高表达的第 3 组基因涉及 "细胞外基质 "和 "对生长因子的反应"（如 HES1、HEY1 和 SIX1）（图 3B）。在第 3 组基因激活的同时，在成肌细胞和肌细胞中高表达的第 4 组基因也被上调并维持在高表达水平，直到最后阶段，GO 术语 "对内源性刺激的反应 "和 "横纹肌组织发育"（如 TCF12、FOS 和 EGR1）被富集。最后，第 5 组基因在轨迹末期被激活，主要涉及 GO 术语 "肌肉细胞分化"（如 和 SOX6）（图 3B）。这些数据说明了肌细胞分化的轨迹，并揭示了整个过程中转录波的有序激活。
Single-cell chromatin accessibility profiling of pig skeletal muscle ontogeny 猪骨骼肌本体发育的单细胞染色质可及性分析
To further investigate the regulatory events in developing myogenic cells, the single-cell chromatin accessibility landscape was analyzed (Additional file 1: Figure S4AC). Using a shared nearest neighbor (SNN) modularity optimization-based clustering algorithm, we obtained 15 distinct clusters of differentially accessible peaks (Additional file 5: Table S4 and Additional file 1: Figure S4D). Clusters 4 and 8 were annotated as myogenic cells for their high accessibility of marker genes associated with myogenic lineage (Figure. S4E, F). Then, the myogenic cells were further divided into 7 distinct sub-clusters with increased resolution (Fig. 4A). 为了进一步研究发育中肌原细胞的调控事件，我们分析了单细胞染色质可及性景观（附加文件1：图S4AC）。利用基于共享近邻（SNN）模块化优化的聚类算法，我们得到了 15 个不同的染色质可及性峰簇（附加文件 5：表 S4 和附加文件 1：图 S4D）。第 4 和第 8 聚类因其与成肌系相关的标记基因的高可及性而被注释为成肌细胞（图 S4E、F）。然后，随着分辨率的提高，成肌细胞被进一步划分为 7 个不同的亚群（图 4A）。
To explore the chromatin accessibility profiles across the seven clusters, we examined the accessibility of selected marker genes from our scRNA-seq data (Fig. 4B). In clusters 2 and 4, we observed greater accessibility of marker genes associated with myocytes (e.g., MYOG, MYH3, MYL1, and CKM) and lower accessibility of genes associated with progenitor cells (e.g., PAX3 and PAX7) (Additional file 6: Table S5 and Fig. 4B). In contrast, cells in clusters 0,3 , and 5 showed greater accessibility for marker genes of progenitor cells and lower accessibility for marker genes of myocytes. Clusters 1 and 6 had mixed signatures, with greater accessibility of PAX7, MSC, MYF5, and MYOG. Based on these observations, we manually annotated the seven clusters as Pax3+ progenitors, myogenic progenitors/myoblasts, and myocytes (Fig. 4C). 为了探索七个集群的染色质可及性概况，我们检查了 scRNA-seq 数据中选定标记基因的可及性（图 4B）。在群集 2 和 4 中，我们观察到与肌细胞相关的标记基因（如 MYOG、MYH3、MYL1 和 CKM）的可及性较高，而与祖细胞相关的基因（如 PAX3 和 PAX7）的可及性较低（附加文件 6：表 S5 和图 4B）。相反，0、3 和 5 群组的细胞对祖细胞标记基因的可及性较高，而对肌细胞标记基因的可及性较低。群组 1 和 6 的特征不一，PAX7、MSC、MYF5 和 MYOG 的可及性较高。基于这些观察结果，我们将这七个群组人工标注为 Pax3+ 祖细胞、成肌祖细胞/成肌细胞和成肌细胞（图 4C）。
To characterize different genomic elements captured by scATAC-seq data, the genome was stratified into promoters, exons, and untranslated regions, introns, and distal regions using the GENCODE annotation [30] (Additional file 1: Figure S5A, B). There was little difference in the proportions of different regions between samples or cell types, with exons accounting for about , promoters, introns, and distal regions accounting for about each, and and untranslated regions accounting for about 5% (Fig. 4D, E). To study the open chromatin heterogeneity across cell types and developmental stages, we derived a cell type-specific chromatin accessibility landscape by conducting pairwise Fisher's exact test for each peak between every cluster. In total, we identified 6422 differentially accessible open chromatin peaks (DAPs) across the 3 cell types, which separated the three cell types perfectly (Additional file 7: Table S6 and Fig. 4F). Among these peaks, most were in regions characterized as distal elements or introns, while relatively few were in the promoters or and untranslated regions (Fig. 4G and Additional file 1: Figure S5C), indicating a critical role for enhancer elements in skeletal muscle development. In addition to the cell type-specific peaks, some cell type-independent open chromatin areas (present across progenitors, myogenic progenitors/ myoblasts, and myocytes) also were found, likely consisting of basal housekeeping genes and regulatory elements (Fig. 4H). The overlapping peaks between progenitors and myogenic progenitors/myoblasts, and between myogenic progenitors/myoblasts and myocytes, were more than that between progenitors and Myocytes, which is consistent with their biological similarities and differentiation process (Fig. 4H). The common peaks of the three cell types were much more than those of the other groups, revealing their close lineage relationship. 为了描述scATAC-seq数据捕获的不同基因组元素的特征，利用GENCODE注释[30]将基因组分为启动子、外显子、 和 非翻译区、内含子和远端区（附加文件1：图S5A、B）。不同样本或细胞类型的不同区域所占比例差别不大，外显子约占 ，启动子、内含子和远端区域各约占 ， 和 非翻译区约占5%（图4D、E）。为了研究不同细胞类型和发育阶段的开放染色质异质性，我们通过对每个聚类之间的每个峰进行成对的费雪精确检验，得出了细胞类型特异的染色质可及性景观。我们在 3 种细胞类型中总共发现了 6422 个差异可及的开放染色质峰（DAPs），将 3 种细胞类型完美地区分开来（附加文件 7：表 S6 和图 4F）。在这些峰中，大多数位于远端元件或内含子区域，而相对较少的 位于启动子或 和 非翻译区（图4G和附加文件1：图S5C），这表明增强子元件在骨骼肌发育中起着关键作用。除了细胞类型特异性峰值外，还发现了一些独立于细胞类型的开放染色质区域（存在于 祖细胞、成肌祖细胞/成肌细胞和肌细胞中），这些区域可能由基本的看家基因和调控元件组成（图 4H）。 祖细胞与成肌祖细胞/成肌细胞、成肌祖细胞/成肌细胞与成肌细胞之间的重叠峰均大于 祖细胞与成肌细胞之间的重叠峰，这与它们的生物学相似性和分化过程相一致（图 4H）。三种细胞类型的共同峰值远远高于其他组，显示了它们之间密切的系谱关系。
Cell type-specific gene regulatory landscape of embryonic skeletal muscle in pigs 猪胚胎骨骼肌细胞特异性基因调控图谱
Cell type-specific chromatin opening and closing events associated with TF binding changes establish the cell 与 TF 结合变化相关的细胞特异性染色质开放和关闭事件建立了细胞
A
IIT
GO: 0009719
GO: 0043043
GO: 0014706
对内源性刺激的反应 肽的生物合成过程 横纹肌组织的发育
response to endogenous stimulus
peptide biosynthetic process
striated muscle tissue development
GO: 0030016
myofibril 肌原纤维
muscle cell differentiation 肌肉细胞分化
GO: 0007519
skeletal muscle tissue development 骨骼肌组织发育
GO: 0060537
muscle tissue development 肌肉组织发育
GO: 0031012
extracellular matrix 细胞外基质
GO: 0070848
response to growth factor 对生长因子的反应
GO: 0061061
muscle structure development 肌肉结构发育
GO: 0046034
ATP metabolic process ATP 代谢过程
GO: 0009144
purine nucleoside triphosphate metabolic process 嘌呤核苷三磷酸代谢过程
GO: 0044429
mitochondrial part 线粒体部分
GO: 0070469
respiratory chain 呼吸链
GO: 0000278
mitotic cell cycle 有丝分裂细胞周期
GO: 1903047
mitotic cell cycle process 有丝分裂细胞周期过程
GO: 0000819
sister chromatid segregation 姐妹染色单体分离
GO: 0000280
nuclear division 核分裂
B Gene clusters B 基因群
Gene signatures 基因特征
Highly expressed in myogenic progenitors 在成肌祖细胞中高表达
E
Myocytes 肌细胞
-1 cluster -1组
Fig. 4 Single-cell chromatin accessibility analysis of pig myogenic cells. A The myogenic cells in the scATAC-seq dataset are shown in the Uniform Manifold Approximation and Projection (UMAP) space, colored by cluster. B Top: Bar plot showing the average accessibility of 13 selected marker genes from our scRNA-seq data considering all myogenic cells. Bottom: dot plot of the standardized accessibility of the marker genes (gene body ) in each of the seven clusters. For each gene, the minimum accessibility value is subtracted, and the result is divided by its maximum accessibility value. The dot size indicates the percentage of cells in each cluster in which the gene of interest is accessible. The standardized accessibility level is indicated by color intensity. C UMAP visualization of the myogenic cells in the scATAC-seq dataset, colored by cell identity. D Percentage distribution of open chromatin elements in each scATAC-seq sample. E Percentage distribution of open chromatin elements in scATAC-seq myogenic cell types. F Heatmap showing cell type-specific differentially accessible peaks (DAPs) (yellow: open chromatin, purple: closed chromatin). G Distribution of open chromatin elements among DAPs in myogenic cell types. Number of shared and unique peaks among snATAC-seq cell types 图 4 猪成肌细胞的单细胞染色质可及性分析。A scATAC-seq 数据集中的成肌细胞显示在均匀簇逼近和投影（UMAP）空间中，并按簇着色。B 顶部：条形图显示了从我们的 scRNA-seq 数据中选取的 13 个标记基因的平均可及性，考虑了所有肌原细胞。下图：七个聚类中每个聚类的标记基因（基因体 ）的标准化可及性点图。对于每个基因，减去最小可及性值，然后将结果除以最大可及性值。点的大小表示每个群组中相关基因可及的细胞百分比。标准化可及性水平用颜色强度表示。C scATAC-seq 数据集中肌原细胞的 UMAP 可视化，按细胞身份着色。D 每个 scATAC-seq 样本中开放染色质元素的百分比分布。E 开放染色质元素在 scATAC-seq 成肌细胞类型中的百分比分布。F 热图显示特定细胞类型的差异可及峰（DAPs）（黄色：开放染色质，紫色：封闭染色质）。G 成肌细胞类型 DAPs 中开放染色质元素的分布。 snATAC-seq细胞类型中共享峰和独特峰的数量
type-specific regulatory landscape, resulting in celltype specification and development. Therefore, the motif enrichment analysis was performed on the cell type-specific open chromatin regions using Genomics. The full list of cell type-specific TF binding motifs is shown in Additional file 8: Table S7. We next correlated 细胞类型特异性调控格局，导致细胞类型的分化和发育。因此，利用 Genomics对细胞类型特异性开放染色质区域进行了基序富集分析。细胞类型特异性TF结合基序的完整列表见附加文件8：表S7。接下来，我们将
the motif enrichment with scRNA-seq TF expression (Additional file 9: Table S8). Using this combined motif enrichment and gene expression approach, the pig skeletal muscle cell type-specific TF landscape was defined (Fig. 5A, B). Correlation of RNA expression and chromatin accessibility in individual single cells revealed two characteristic patterns of ATAC-RNA pairs: (i) RNA expression of TFs directly matches accessibility of corresponding TF bindings sites as exemplified for ARID3A, MEIS2, MEIS1, HOXB4, and HOXD4 in progenitors, and MYOG, MYOD1, KLF2, and SOX6 in myocytes, suggesting that these TFs actively regulate their respective target genes at the specific developmental stage; (ii) RNA expression of TFs precedes the increase in accessibility of the corresponding TF binding sites. This scenario was apparent for MYF6, MYF5, SNAI1, and TGIF1, which reached their highest expression levels in myogenic progenitors/ myoblasts but showed the strongest motif enrichment in myocytes, suggesting that additional epigenetic regulation could occur before TFs take action (Fig. 5A, B). 通过scRNA-seq TF的表达来确定基因主题富集（附加文件9：表S8）。利用这种结合了基序富集和基因表达的方法，确定了猪骨骼肌细胞特异性 TF 图谱（图 5A、B）。单个单细胞中 RNA 表达和染色质可及性的相关性揭示了 ATAC-RNA 对的两种特征模式：(i) TF 的 RNA 表达与相应 TF 结合位点的可及性直接匹配，例如 祖细胞中的 ARID3A、MEIS2、MEIS1、HOXB4 和 HOXD4 以及肌细胞中的 MYOG、MYOD1、KLF2 和 SOX6。这种情况在 MYF6、MYF5、SNAI1 和 TGIF1 中很明显，它们在成肌祖细胞/成肌细胞中的表达水平最高，但在肌细胞中却表现出最强的基序富集，这表明在 TF 采取行动之前可能会发生额外的表观遗传调控（图 5A，B）。
To study the putative target genes of TFs, single-cell regulatory network inference and clustering (SCENIC) was performed to examine TF regulon activity [31]. The activity of each regulon in each cell was quantified and then binarized to "on" or "off" based on activity distribution across cells. The SCENIC results indicated strong enrichment of HOXB9 and MEIS1 regulon activity in Pax progenitors; RAD21, EZH2, and CTCF activity in myogenic progenitors; TCF12, EGR1, and FOSB activity in myoblasts; and MYOD1, MEF2C, MYOG, and SOX6 activity in myocytes (Additional file 10: Table S9 and Fig. 5C). Although the activity of and , which belong to the same family as , was elevated in myocytes, it was not as significant as MEF2C, indicating that they did not play a leading role in myogenesis. SCENIC also successfully inferred multiple downstream target genes. The complete list of regulons and their respective predicted target genes can be found in Additional file 11: Table S10. The scaled and binarized regulon activity is available in Additional file 12: Table S11. Examples of regulon activity, corresponding TF expression, and predicted target gene expression are depicted in Fig. 5D and Additional file 1: Figure S6A. TFs such as SOX6, TEAD4, and FOXO1 were predicted as targets of skeletal muscle-specific TF MYOG, and FOS was predicted as a target of , indicating a critical transcriptional hierarchy of skeletal muscle development. Corresponding chromatin accessibility in scATAC data for these TFs and predicted target genes are shown in Figure. S6B. 为了研究 TF 的推定靶基因，研究人员进行了单细胞调控网络推断和聚类（SCENIC），以检查 TF 调控子的活性[31]。每个细胞中每个调控子的活性都被量化，然后根据细胞间的活性分布二值化为 "开 "或 "关"。SCENIC的结果表明，HOXB9和MEIS1调控子的活性在Pax 祖细胞中强富集；RAD21、EZH2和CTCF的活性在成肌细胞中强富集；TCF12、EGR1和FOSB的活性在成肌细胞中强富集；MYOD1、MEF2C、MYOG和SOX6的活性在成肌细胞中强富集（附加文件10：表S9和图5C）。虽然与 同属一个家族的 和 在肌细胞中的活性升高，但不如MEF2C那么显著，这表明它们在肌形成过程中没有发挥主导作用。SCENIC 还成功地推断出了多个下游靶基因。调控子及其各自预测的靶基因的完整列表见附加文件 11：表 S10。按比例和二值化的调节子活性见附加文件 12：表 S11。图 5D 和附加文件 1 描述了调节子活性、相应 TF 表达和预测靶基因表达的例子：图 S6A。SOX6、TEAD4和FOXO1等TF被预测为骨骼肌特异性TF MYOG的靶标，FOS被预测为 的靶标，这表明骨骼肌发育的转录层次至关重要。图S6B显示了这些TF和预测靶基因在scATAC数据中相应的染色质可及性。S6B。
Chromatin dynamics of progenitor differentiation 祖细胞分化的染色质动力学
The pseudotime trajectory in the scATAC-seq dataset was evaluated, which resulted in a similar cellular differentiation trajectory to scRNA-seq dataset (Fig. 6A). We integrated the scRNA-seq and scATAC-seq datasets to correlate and cross-validate gene expression profiles and the chromatin accessibility landscape in the myogenic cells using the Harmony algorithm [32]. The coembedded Uniform Manifold Approximation and Projection (UMAP) plots with cell type assignment for the scATACseq and scRNA-seq data suggested that the changes in chromatin accessibility and corresponding gene transcript expression in most myocytes occurred in a synchronous manner, whereas in other cell types, it was not fully synchronized, suggesting that other regulation is involved (Fig. 6B, C). Correlations between cell types of scATAC-seq and scRNA-seq cell types were computed with scRNA-seq variable genes (Fig. 6D). These results obtained from two independent approaches demonstrate that our two datasets are highly concordant and cross-validated. 我们对 scATAC-seq 数据集的伪时间轨迹进行了评估，其结果与 scRNA-seq 数据集的细胞分化轨迹相似（图 6A）。我们整合了 scRNA-seq 和 scATAC-seq 数据集，利用和谐算法（Harmony algorithm）[32]关联和交叉验证了基因表达谱和成肌细胞染色质可及性图谱。scATACseq和scRNA-seq数据与细胞类型分配的共嵌统一图谱逼近和投影（UMAP）图表明，大多数肌细胞中染色质可及性和相应基因转录本表达的变化是同步发生的，而在其他细胞类型中，变化并不完全同步，这表明还涉及其他调控（图6B、C）。利用 scRNA-seq 可变基因计算了 scATAC-seq 细胞类型与 scRNA-seq 细胞类型之间的相关性（图 6D）。通过两种独立方法得出的这些结果表明，我们的两个数据集高度一致并经过交叉验证。
We next performed pseudotime ordering of the chromatin accessibility-based TF motif enrichment of individual cells and correlated changes of the TF motif patterns with TF expression (Fig. 7A, B). To this end, we also investigated TFs and target genes differentially expressed over scRNA-seq pseudotime. We noticed good concordance of time-dependent changes in TF and predicted target gene expression along with motif enrichment, suggesting that a set of TFs cooperatively regulates myogenic differentiation (Fig. 7C and Additional file 1: Figure S7). 接下来，我们对单个细胞基于染色质可及性的TF基序富集进行了伪时间排序，并将TF基序模式的变化与TF表达相关联（图7A, B）。为此，我们还研究了在 scRNA-seq 伪时间中差异表达的 TF 和靶基因。我们注意到 TF 和预测的靶基因表达随时间的变化与主题富集有很好的一致性，这表明一组 TF 协同调控了肌原分化（图 7C 和附加文件 1：图 S7）。
EGR1 and RHOB play critical roles in myogenesis EGR1 和 RHOB 在肌肉生成过程中发挥关键作用
Although the roles of several identified TFs in myoblast specification and differentiation have been established, 尽管已经确定了几种已识别的 TFs 在成肌细胞规格化和分化过程中的作用、
Fig. 5 (See legend on previous page.) 图 5 （见前页图例。）
Fig. 6 Integrated analysis of scATAC-seq and scRNA-seq data. A The pseudotime trajectory in the scATAC-seq dataset. B UMAP representation of scATAC-scRNA integration results. Cells are colored by technology (scATAC red, scRNA=blue). C UMAP representation of scATAC-scRNA integration results. Cells are colored by cell type assignment. D Dot plot showing scATAC-scRNA integration cell type assignment using confusion matrix. Each column represents the original cell type assignment of scRNA-seq data, and each row represents the cell type assignment predicted after integration with scATAC-seq data. The size of the dots represents the absolute value of the correlation, and the red and gray dots represent the positive and negative correlations, respectively 图 6 scATAC-seq 和 scRNA-seq 数据的综合分析。A scATAC-seq 数据集的伪时间轨迹。B scATAC-scRNA 整合结果的 UMAP 表示。细胞按技术着色（scATAC 红色，scRNA=蓝色）。C scATAC-scRNA 整合结果的 UMAP 表示。细胞按细胞类型着色。D 点阵图显示使用混淆矩阵进行的 scATAC-scRNA 整合细胞类型分配。每列代表 scRNA-seq 数据的原始细胞类型分配，每行代表与 scATAC-seq 数据整合后预测的细胞类型分配。点的大小代表相关性的绝对值，红点和灰点分别代表正相关性和负相关性
the expression and functions of many other genes that present specific expression profiles over scRNA-seq pseudotime have not been well studied. In consideration of the possibility that pig is used as a model of muscle disease in the future, and in order to identify the conserved functional genes related to muscle development among mammals, we compared the expression of these genes in skeletal muscle from wild-type and Duchenne muscular dystrophy (DMD) mice using an RNA-Seq dataset (GSE162455). Consistent with the patterns of classical myogenic genes, most of the genes upregulated along the pseudotime trajectory were expressed at higher levels in DMD mice than that in wild-type mice (Additional file 4: Table S3 and Additional file 1: Figure S8A, B), suggesting that these genes may be induced by muscle regeneration in muscular dystrophy mice and play an important role in myogenesis. The expression levels of genes downregulated along the pseudotime trajectory were also downregulated in DMD mice, suggesting that they are involved in muscular dystrophy and may not contribute 在scRNA-seq伪时间中呈现特定表达谱的许多其他基因的表达和功能尚未得到很好的研究。考虑到猪将来可能被用作肌肉疾病的模型，为了确定哺乳动物中与肌肉发育相关的保守功能基因，我们使用 RNA-Seq 数据集（GSE162455）比较了这些基因在野生型小鼠和杜氏肌营养不良症（DMD）小鼠骨骼肌中的表达情况。与经典肌生成基因的模式一致，大多数沿伪时间轨迹上调的基因在 DMD 小鼠中的表达水平高于野生型小鼠（附加文件 4：表 S3 和附加文件 1：图 S8A、B），这表明这些基因可能由肌肉萎缩症小鼠的肌肉再生诱导，并在肌生成中发挥重要作用。沿着伪时间轨迹下调的基因的表达水平在 DMD 小鼠中也下调了，这表明这些基因参与了肌肉萎缩症的发生，可能对肌肉生成没有贡献。
(See figure on next page.) (见下页图）
Fig. 7 Activity and RNA expression dynamics of TFs along the pseudotime trajectory. A Heatmap showing the activity of TFs along the differentiation trajectory. B TF expression heatmap corresponding to the motif enrichment along the differentiation trajectory. Pseudotime-dependent chromatin accessibility and gene expression changes along the myogenic lineages. The first column shows the dynamics of the Genomics TF enrichment score. The second column shows the dynamics of TF gene expression values, and the third and fourth columns represent the dynamics of SCENIC-reported target gene expression values of corresponding TFs. Error bars denote confidence intervals of local polynomial regression fitting. Additional examples are given in Additional file 1: Figure S7 图 7 TFs 沿伪时间轨迹的活性和 RNA 表达动态。A 热图显示 TFs 沿分化轨迹的活性。B TF表达热图，与沿分化轨迹的主题富集相对应。 依赖于伪时间的染色质可及性和基因表达沿成肌系的变化。第一列显示了 Genomics TF富集得分的动态变化。第二列显示了TF基因表达值的动态变化，第三列和第四列代表了SCENIC报告的相应TF靶基因表达值的动态变化。误差条表示局部多项式回归拟合的 置信区间。更多例子见附加文件 1：图 S7
Fig. 7 (See legend on previous page.) 图 7 （见前页图例。）
to terminal myoblast differentiation (Additional file 4: Table S3 and Additional file 1: Figure S8C). The scRNA analysis highlighted that the expression of EGR1 and its predicted target gene gradually increased along the pseudotime trajectory, and scATAC analysis showed that EGR1 reached its strongest motif enrichment in myocytes (Fig. 7C). The peaks in the vicinity of EGR1 and , which were located at +3979 bp and of the EGR1 and RHOB transcriptional start sites (TSSs) respectively, were the most accessible in myocytes (Additional file 1: Figure S9). Based on gene expression correlation and binding motif analysis, we identified a total of 215 high-confidence annotated TF-target pairs between 12 TFs and 75 target genes, and constructed a putative gene regulatory network associated with skeletal muscle development (Fig. 8A). Among them, the number of links between transcription factor EGR1 and myogenesis-related target genes is second only to the classical myogenic transcription factor MYOD1 (Fig. 8B). RHOB is an important regulator of cell and tissue morphology and function, acting mainly through the cellular cytoskeleton [33, 34]. A few studies have revealed that RHOB is a key mediator during diverse cellular and physiological processes like cell division, cell migration in smooth muscle cells [35, 36]. Totally, we speculated that EGR1 and RHOB are likely to play positive roles during myogenic differentiation in pigs, so they were selected for in vitro functional validation. 附加文件 4：表 S3 和附加文件 1：图 S8C）。scRNA 分析显示，EGR1 及其预测靶基因 的表达量沿着伪时间轨迹逐渐增加，scATAC 分析显示 EGR1 在肌细胞中达到了最强的图案富集（图 7C）。EGR1和RHOB转录起始位点（TSSs）的+3979 bp和 附近的峰是肌细胞中最容易获得的（附加文件1：图S9）。基于基因表达相关性和结合基序分析，我们在12个TFs和75个靶基因之间共鉴定出215对高置信度注释的TF-靶基因对，并构建了一个与骨骼肌发育相关的假定基因调控网络（图8A）。其中，转录因子EGR1与肌生成相关靶基因之间的联系数量仅次于经典的肌生成转录因子MYOD1（图8B）。RHOB 是细胞和组织形态与功能的重要调节因子，主要通过细胞骨架发挥作用 [33，34]。一些研究表明，RHOB 是平滑肌细胞分裂、细胞迁移等多种细胞和生理过程的关键介质[35, 36]。总之，我们推测 EGR1 和 RHOB 可能在猪的肌原分化过程中发挥积极作用，因此选择它们进行体外功能验证。
To confirm the dynamic expression of EGR1 and during myogenic differentiation, quantitative PCR (qPCR) was performed on porcine primary myogenic cells (PPMCs) at several differentiation points (day 0 , day 2 , day 4 , day 6 , and day 8 ). It was found that EGR1 and RHOB had the same expression pattern as the well-known myogenic differentiation makers (Fig. 8C). To validate the effect of EGR1 and RHOB on myogenesis progression, they were overexpressed in PPMCs and mouse C2C12 myoblasts. Ethynyl-2'-deoxyuridine (EdU) incorporation and immunofluorescence assays showed that EGR1 overexpression did not influence cell proliferation but promoted myogenic differentiation, inducing a significant increase in the fusion index (Fig. 8D, E, and Additional file 1: Figure S10). In line with this, the expression of myogenic differentiation markers increased when EGR1 was overexpressed (Fig. 8F). Consistent with EGR1, RHOB overexpression did not change the proliferation ability of myoblasts, but the expression level of MYOD1 was significantly upregulated (Additional file 1: Figure S11A-C, and Fig. 8G). This prompted us to further verify whether regulates myoblast fate commitment. Immunofluorescence co-staining of PAX7 and MYOD showed that accelerated the transformation of myogenic progenitors into MYOD myoblasts (Fig. 8H, I). After induction of differentiation, the overexpression group formed more myotubes with a significantly increased the fusion index (Fig. 8J, and Additional file 1: Figure S11D). These findings revealed that EGR1 and RHOB are critical regulators of pig embryonic myogenesis. 为了证实EGR1和 在成肌细胞分化过程中的动态表达，对猪原代成肌细胞（PPMCs）在几个分化点（第0天、第2天、第4天、第6天和第8天）进行了定量PCR（qPCR）检测。结果发现，EGR1 和 RHOB 与众所周知的肌原分化制造者具有相同的表达模式（图 8C）。为了验证 EGR1 和 RHOB 对成肌过程的影响，我们在 PPMCs 和小鼠 C2C12 成肌细胞中过表达了它们。乙炔基-2'-脱氧尿苷（EdU）掺入和免疫荧光检测表明，EGR1的过表达不会影响细胞增殖，但会促进成肌分化，诱导融合指数显著增加（图8D、E和附加文件1：图S10）。与此相一致，当 EGR1 被过表达时，成肌分化标记物的表达也增加了（图 8F）。与 EGR1 一致，RHOB 的过表达并没有改变成肌细胞的增殖能力，但 MYOD1 的表达水平却显著上调（附加文件 1：图 S11A-C 和图 8G）。这促使我们进一步验证 是否调控了成肌细胞的命运承诺。PAX7和MYOD的免疫荧光共染显示， 加速了 肌原细胞向 MYOD 肌母细胞的转化（图8H，I）。诱导分化后， 过表达组形成了更多的肌管，融合指数显著增加（图 8J；附加文件 1：图 S11D）。这些发现揭示了EGR1和RHOB是猪胚胎肌形成的关键调控因子。
Cell-cell communications 小区-小区通信
To predict the cellular communications involved in pig skeletal muscle ontogeny, we evaluated the potential cell-cell interactions by using CellPhoneDB [37]. First, the cellular interactions between myogenic cells were analyzed, and it was found that the interactions involving Pax progenitor were predicted to be more significant than those involving myogenic progenitor, myoblast, and myocyte (Fig. 9A). This result illustrated that the stronger the characteristics of stem cells, the greater the possibility of their interaction with other cells. Then, all somite cell populations, including mesenchymal cells, fibroblasts, epithelial cells, neural stem cells, osteogenic cells, neurons, neurogliacytes, endothelial cells, chondrocytes, and myogenic cells, were included in the analysis. Interestingly, the interaction between cells of the same cell 为了预测参与猪骨骼肌本体形成的细胞通讯，我们使用 CellPhoneDB [37] 评估了潜在的细胞-细胞相互作用。首先，我们分析了肌原细胞之间的细胞相互作用，结果发现，涉及Pax 祖细胞的相互作用比涉及肌原祖细胞、肌母细胞和肌细胞的相互作用更显著（图9A）。这一结果说明，干细胞的特征越强，与其他细胞相互作用的可能性就越大。然后，所有体细胞群，包括间充质细胞、成纤维细胞、上皮细胞、神经干细胞、成骨细胞、神经元、神经胶质细胞、内皮细胞、软骨细胞和成肌细胞，都被纳入分析。有趣的是，同一细胞的细胞之间的相互作用
Fig. 8 (See legend on previous page.) 图 8 （见前页图例。）
type was weaker than that between different cell types. Myogenic cells tended to communicate with fibroblasts, osteogenic cells, chondrocytes, neurons, and endothelial cells. In myogenic cells, except for the stronger interaction between progenitors and other cells, the communication between myogenic progenitors, myoblasts, and myocytes was less (Fig. 9B). Among the identified interactions, many were related to IGF2 and the corresponding receptors IGF1R and IGF2R (Fig. 9C), which is consistent with the well-known roles of IGF signaling in skeletal muscle development [38]. The myogenic cells showed a significant ligand-receptor interaction between FGF9 and its receptors FGFR4 and FGFR1 (Fig. 9C), underscoring its critical role in myogenic differentiation [39]. The cellular interactions between myogenic cells and non-myogenic cells were mainly related to CD74, IGF2, PTN, ERBB3, FN1, and CADM1 (Fig. 9D). ERBB3-NRG1 signal plays a critical role in sustainable myogenesis by restraining myogenic progenitors from precocious differentiation [40]. Because relatively little is known about some of these genes, the significance of these putative interactions requires further investigation. For example, in skeletal muscle, PTN is upregulated during myogenesis and postsynaptic induction [41], but little is known regarding its effects on muscle formation. Interestingly, fibroblasts, osteoblasts, and chondrocytes communicate with progenitors, myogenic progenitors, and myoblasts via PTN-PTPRS interaction pair, but this cellular communication no longer exists in differentiated myocytes (Fig. 9D). 与不同类型细胞之间的通讯相比，肌原性细胞与成纤维细胞、成骨细胞、软骨细胞和内皮细胞之间的通讯较弱。成肌细胞倾向于与成纤维细胞、成骨细胞、软骨细胞、神经元和内皮细胞交流。在成肌细胞中，除了 祖细胞与其他细胞之间的相互作用较强之外，成肌细胞祖细胞、成肌细胞和成肌细胞之间的交流较少（图9B）。在已发现的相互作用中，许多与 IGF2 以及相应的受体 IGF1R 和 IGF2R 有关（图 9C），这与众所周知的 IGF 信号在骨骼肌发育中的作用是一致的[38]。成肌细胞显示 FGF9 与其受体 FGFR4 和 FGFR1 之间存在明显的配体-受体相互作用（图 9C），这突出表明了 FGF9 在成肌细胞分化中的关键作用[39]。成肌细胞与非成肌细胞之间的细胞相互作用主要与 CD74、IGF2、PTN、ERBB3、FN1 和 CADM1 有关（图 9D）。ERBB3-NRG1信号通过抑制成肌祖细胞的过早分化，在可持续成肌过程中发挥着关键作用[40]。由于对其中一些基因的了解相对较少，因此需要进一步研究这些假定相互作用的意义。例如，在骨骼肌中，PTN 在肌生成和突触后诱导过程中上调[41]，但人们对其对肌肉形成的影响知之甚少。有趣的是，成纤维细胞、成骨细胞和软骨细胞通过 PTN-PTPRS 相互作用对与 祖细胞、成肌祖细胞和成肌细胞交流，但这种细胞交流在分化的肌细胞中已不复存在（图 9D）。
Discussion 讨论
In the present study, we depicted the first gene expression and open chromatin maps in developing pig somites and myotomes at single-cell resolution. Using this dataset, we explored the regulatory dynamics along the myogenic differentiation trajectory and defined cell-type-specific regulatory networks. We also investigated key TFs for embryonic myogenesis and cell-cell interactions associated with skeletal muscle development. These results shed light on the upstream regulatory networks of pig skeletal muscle ontogeny. 在本研究中，我们首次以单细胞分辨率描绘了发育中猪体节和肌节的基因表达和开放染色质图谱。利用该数据集，我们探索了肌形成分化轨迹的调控动态，并定义了细胞类型特异性调控网络。我们还研究了胚胎成肌过程中的关键 TFs 以及与骨骼肌发育相关的细胞-细胞相互作用。这些结果揭示了猪骨骼肌本体发育的上游调控网络。
The widespread changes and period specificity of gene expression precisely capture the transcriptional characteristics of myogenic differentiation. For example, genes involved in vital metabolic pathways, such as "ATP metabolic process" and "purine nucleoside triphosphate metabolic process", had their highest expression levels at the beginning of the differentiation trajectory, revealing higher metabolic activity in progenitor cells [42, 43]. The high expression of genes associated with "mitotic cell cycle" and "sister chromatid segregation" in myogenic progenitors indicated that these cells had the strongest proliferation activity during myogenesis, and they expand the pool of myogenic cells for the subsequent myofiber formation [44, 45]. The inactivation of cell cycle-related genes in myoblasts and myocytes coincided with cells exiting the cell cycle and beginning to differentiate [46]. Moreover, genes associated with "myofibril" and "muscle tissue development" were over-represented at the end of the differentiation trajectory, representing the skeletal muscle terminal differentiation processes [9, 47]. These observations suggest that distinct transcriptional programs and biological processes are activated or depressed at specific stages, which provides valuable clues for understanding their functions in skeletal muscle development. 基因表达的广泛变化和时期特异性恰好捕捉到了肌原分化的转录特征。例如，参与重要代谢途径的基因，如 "ATP 代谢过程 "和 "嘌呤核苷三磷酸代谢过程"，在分化轨迹的起始阶段表达水平最高，表明祖细胞的代谢活性较高[42, 43]。与 "有丝分裂细胞周期 "和 "姐妹染色单体分离 "相关的基因在成肌祖细胞中的高表达表明，这些细胞在成肌过程中具有最强的增殖活性，它们为随后的肌纤维形成扩充了成肌细胞池[44, 45]。肌母细胞和肌细胞中细胞周期相关基因的失活与细胞退出细胞周期并开始分化相吻合[46]。此外，与 "肌原纤维 "和 "肌肉组织发育 "相关的基因在分化轨迹的末端比例过高，代表了骨骼肌的末端分化过程[9, 47]。这些观察结果表明，不同的转录程序和生物过程在特定阶段被激活或抑制，这为了解它们在骨骼肌发育中的功能提供了宝贵的线索。
During mammalian development, differences in chromatin state coincide with cellular differentiation and reflect changes in the gene regulatory landscape [48]. Our studies revealed a similar pattern between open chromatin information and gene expression when tracked with myogenic cell differentiation. The expression levels and motif enrichment of HOXB8 and HOXA1, which are known regulators of smooth muscle or extraocular muscle differentiation [49, 50], gradually decreased along the pseudotime trajectory. Myogenic differentiation correlated with increased expression of SOX6 and SOX8, among which, SOX6 has been known to play important roles in myoblast differentiation [51]. This is the first report of SOX8 in skeletal muscle development. 在哺乳动物的发育过程中，染色质状态的差异与细胞分化相吻合，反映了基因调控格局的变化 [48]。我们的研究揭示了开放染色质信息与基因表达之间的类似模式，并将其与成肌细胞分化联系起来。HOXB8和HOXA1是已知的平滑肌或眼外肌分化的调控因子[49, 50]，它们的表达水平和基因位点富集度沿着伪时间轨迹逐渐降低。成肌分化与 SOX6 和 SOX8 的表达增加有关，其中 SOX6 在成肌细胞分化中发挥重要作用[51]。这是首次报道 SOX8 在骨骼肌发育过程中的作用。
Fig. 9 (See legend on previous page.) 图 9 （见前页图例。）
Numerous studies have reported that the MYC family of basic helix-loop-helix zipper (bHLHZ) transcription factors and their binding protein MAX control multiple cellular functions and are widely involved in oncogenesis [52, 53]. C-MYC inhibits myogenic differentiation by suppressing MYOD1 expression [54]. Here, we revealed 大量研究表明，MYC 家族的碱性螺旋环-螺旋拉链（bHLHZ）转录因子及其结合蛋白 MAX 控制着多种细胞功能，并广泛参与肿瘤的发生 [52，53]。C-MYC 通过抑制 MYOD1 的表达来抑制肌源性分化 [54]。在此，我们发现
that the expression of and remained high in the late-pseudotime cells, and their motif enrichment increased with differentiation, suggesting that these genes also play an essential role in the formation of skeletal muscle during pig embryo development. 结果表明， 和 在伪时后期细胞中的表达量仍然很高，而且它们的主题富集度随着分化而增加，这表明这些基因在猪胚胎发育过程中对骨骼肌的形成也起着至关重要的作用。
Integrative single-cell RNA-seq and ATAC-seq analysis also demonstrated a cell type-specific gene regulatory landscape in pig embryonic skeletal muscle. For example, motif enrichment analysis showed that the TF binding sites of TCF12 reached their highest level of accessibility in myocytes, which is consistent with the reported function of TCF12 in regulating myocyte differentiation [55]. However, the analysis of TF regulon activity with SCENIC indicated the most vigorous TCF12 regulon activity in myoblasts. This scenario was also apparent for EGR1. The phenomenon that RNA expression of TFs precedes accessibility of related TF binding sites suggests that additional epigenetic regulation might occur before TFs take action [56]. Furthermore, although individual TFs have bound to cognate motifs, this might not be a sufficient condition to initiate target gene transcription, which would lead to mismatches between RNA-seq and ATAC-seq results [57, 58]. Of course, such conditions do not apply to all genes. TFs such as MYOD1, MYOG, SOX6, and MEF2D, critical mediators of myogenesis, were indicated to reach their highest expression level, strongest motif enrichment, and most vigorous regulon activity in myocytes, which is likely responsible for their specific and vital regulatory roles in terminally differentiation. 单细胞RNA-seq和ATAC-seq整合分析还显示了猪胚胎骨骼肌中细胞类型特异性基因调控格局。例如，图案富集分析表明，TCF12 的 TF 结合位点在肌细胞中的可及性达到最高水平，这与报道的 TCF12 调控肌细胞分化的功能相一致 [55]。然而，用 SCENIC 对 TF 调节子活性的分析表明，在成肌细胞中 TCF12 调节子的活性最强。这种情况在 EGR1 中也很明显。TFs的RNA表达先于相关TF结合位点的可及性这一现象表明，在TFs发挥作用之前可能会发生额外的表观遗传调控[56]。此外，虽然单个 TF 与同源基团结合，但这可能不是启动靶基因转录的充分条件，这将导致 RNA-seq 和 ATAC-seq 结果不匹配 [57，58]。当然，这些条件并不适用于所有基因。MYOD1、MYOG、SOX6和MEF2D等TFs是肌形成的关键介导因子，它们在肌细胞中的表达水平最高、基因主题富集最强、调控子活性最活跃，这可能是它们在终末分化中发挥特殊和重要调控作用的原因。
Distal regulatory regions (DRR) mediated chromosomal interactions play an essential role in differentiating cells [59]. Distal promoter contacts are highly cell-type specific [60] and form networks of co-regulated genes correlated with their biological functions. The cis-regulatory contact upon lineage commitment is a dynamic process that includes acquiring and losing specific promoter interactions [61]. Our results revealed that most of the differential accessible peaks were in regions characterized as distal elements or introns. Genome-wide studies showed that MRFs predominantly occupied the extragenic regions, marked by acetylated histones [62, 63]. The DRR of MYOD1 were transcribed to RNA which contributed to establishing a cell-type-specific transcriptional circuitry by directing chromatin-remodeling events [64]. Thus, we hypothesized that intronic enhancers [65] and DRRs play an important role in regulating myogenic differentiation, consistent with mouse kidney development [58]. However, because cis-acting regulatory elements can be located in kilo-bases away from their target genes, it is challenging to identify the true functional targets of regulatory elements [66]. 远端调控区（DRR）介导的染色体相互作用在细胞分化过程中发挥着至关重要的作用 [59]。远端启动子接触具有高度的细胞类型特异性 [60]，并形成与其生物功能相关的共调基因网络。系承时的顺式调控接触是一个动态过程，包括特定启动子相互作用的获得和丧失 [61]。我们的研究结果表明，大多数差异可及峰位于远端元件或内含子区域。全基因组研究表明，MRFs 主要占据以乙酰化组蛋白为标志的基因外区域[62, 63]。MYOD1 的 DRR 转录为 RNA，通过引导染色质重塑事件，有助于建立细胞类型特异性转录回路 [64]。因此，我们推测内含子增强子[65]和 DRRs 在调节肌原分化中起着重要作用，这与小鼠肾脏发育[58]一致。然而，由于顺式作用的调控元件可能位于离其靶基因几千个碱基的地方，因此要确定调控元件的真正功能靶点是很有挑战性的[66]。
Comparative analysis of transcriptome dynamics during porcine myogenic differentiation with transcriptome characteristics in the muscular dystrophy mouse model suggested that the molecular regulatory network of porcine skeletal muscle development is similar to that of muscular dystrophy. Therefore, in addition to partially revealing the molecular mechanism of skeletal muscle development in porcine embryo, this study also provides clues for studying the occurrence and development of muscular dystrophy. Furthermore, our study not only provided a large amount of high-quality data but also identified a lot of new genes involved in myogenesis. EGR1, an immediate-early response zinc-finger transcription factor with various functions in numerous contexts, including the regulation of growth and differentiation , was validated to promote myoblast differentiation. RHOB, a member of the small Rho GTPase family [69], was demonstrated to promote the fate commitment of myogenic progenitors into myoblasts, thus actively regulating myoblast differentiation. 猪肌原分化过程中转录组动态与肌肉营养不良症小鼠模型转录组特征的比较分析表明，猪骨骼肌发育的分子调控网络与肌肉营养不良症的分子调控网络相似。因此，本研究除了部分揭示猪胚胎骨骼肌发育的分子机制外，还为研究肌肉萎缩症的发生和发展提供了线索。此外，我们的研究不仅提供了大量高质量的数据，还发现了许多参与肌肉发生的新基因。EGR1是一种即时早期反应锌指转录因子，在多种情况下具有多种功能，包括调节生长和分化 ，研究证实它能促进肌母细胞分化。RHOB 是小 Rho GTPase 家族的成员 [69]，已被证实能促进成肌细胞祖细胞向成肌细胞的命运承诺，从而积极调节成肌细胞的分化。
There are significant differences in skeletal muscle development between pig breeds differing in growth rate, muscle production, muscle fiber diameter, and meat quality [70]. In this study, the results of immunohistochemistry, scRNA, and scATAC analysis did not show significant differences between ZZ and DZ in early embryonic development, despite their considerable difference in meat production. We speculate that the early embryonic skeletal muscle development of ZZ and DZ is conserved, which may be due to their close genetic relationship. Whether there are differences between purebred Tibetan and Duroc pigs in early embryonic skeletal muscle development and which stage of skeletal muscle development is responsible for the difference in meat production between ZZ and Duroc remain to be studied in the future. 不同猪种在生长速度、肌肉产量、肌纤维直径和肉质等方面的骨骼肌发育存在显著差异 [70]。在本研究中，免疫组化、scRNA 和 scATAC 分析结果表明，尽管 ZZ 和 DZ 在产肉量上存在显著差异，但它们在早期胚胎发育上并无明显差异。我们推测，ZZ 和 DZ 的早期胚胎骨骼肌发育是一致的，这可能与它们的遗传关系密切有关。纯种藏香猪和杜洛克猪在早期胚胎骨骼肌发育方面是否存在差异，以及骨骼肌发育的哪个阶段是造成藏香猪和杜洛克猪产肉量差异的原因，还有待今后的研究。
Conclusions 结论
In summary, our study provided critical insight into the cell type-specific gene regulatory network and cell differentiation program of pig embryonic myogenesis, and screened essential genes that regulate the development of skeletal muscle (Fig. 9E), which will provide a theoretical basis for the study of pig breeding and human muscle diseases. 总之，我们的研究对猪胚胎肌发生的细胞类型特异性基因调控网络和细胞分化程序提供了重要的见解，筛选出了调控骨骼肌发育的重要基因（图 9E），这将为猪育种和人类肌肉疾病的研究提供理论依据。
Methods 方法
Preparation of single-cell suspensions 制备单细胞悬浮液
Somite tissues [a mixture of embryos from the same sow] of ZZ and DZ at E16, E18, E21, and myotome tissues [a mixture of embryos from the same sow] at E28 were isolated and transferred to RPMI 1640 分离E16、E18和E21期ZZ和DZ的体节组织[来自同一母猪的胚胎 的混合物]和E28期的肌体组织[来自同一母猪的胚胎 的混合物]，并将其转移到RPMI 1640中。
medium (GIBCO, 11875093) containing 10% fetal bovine serum (HyClone, SH30070.03) on ice. Then, the tissues were washed with PBS three times and transferred to prewarmed digestion medium containing DNase I (Sigma-Aldrich, DN25) and Collagenase I (Sigma, C2674) in RPMI 1640. The digestion reactions were shaken vigorously for 30 s and further incubated at for 30 to 40 min in an incubator with general shaking every 6 min to release cells. The released cells were passed through a cell strainer (BD, 352350) and were collected in tubes. Sample viability was assessed via Trypan Blue (Thermo Fisher) and an automatic cell counter (Countstar). Cells derived from embryos of different stages and from different pigs were sorted and processed independently in all experiments. Thus, each cell can be traced back to its specific embryo and tissue of origin. 将含有 10%胎牛血清（HyClone，SH30070.03）的培养基（GIBCO，11875093）置于冰上。然后，用 PBS 冲洗组织三次，并将其转移到预热的消化培养基中，该培养基含有 RPMI 1640 中的 DNase I（Sigma-Aldrich，DN25）和 胶原酶 I（Sigma，C2674）。将消化反应剧烈振荡 30 秒，然后在培养箱中以 的温度培养 30 至 40 分钟，每 6 分钟摇动一次，以释放细胞。释放出的细胞通过 细胞过滤器（BD，352350），收集到 管中。通过胰蓝（赛默飞世尔）和自动细胞计数器（Countstar）评估样品的存活率。在所有实验中，来自不同阶段胚胎和不同猪的细胞都是独立分拣和处理的。因此，每个细胞都可以追溯到其特定的胚胎和组织来源。
Droplet-based scRNA-seq datasets were produced using a Chromium system (10× Genomics, PN120263) following the manufacturer's instructions. After droplet generation, samples were transferred into prechilled 8-well tube strips (Eppendorf), and reverse transcription was performed using a Veriti 96 -well thermal cycler (Thermo Fisher). After reverse transcription, cDNA was recovered using the Recovery Agent provided by Genomics, followed by clean-up with Silane Dynabeads (Thermo Fisher) as outlined in the user guide. Purified cDNA was amplified for 12 cycles before being cleaned up using SPRIselect beads (Beckman). Samples were diluted at a ratio of 1:4 and run on a Bioanalyzer (Agilent Technologies) to determine the cDNA concentration. cDNA libraries were prepared as outlined by the Single Cell Reagent Kit v3 user guide with appropriate modifications to the PCR cycles based on the calculated cDNA concentration (as recommended by Genomics). 基于液滴的 scRNA-seq 数据集是使用 Chromium 系统（10× Genomics，PN120263）按照制造商的说明制作的。液滴生成后，将样本转移到预冷的 8 孔试管条（Eppendorf）中，并使用 Veriti 96 孔热循环仪（Thermo Fisher）进行反转录。反转录后，使用 Genomics公司提供的Recovery Agent回收cDNA，然后用Silane Dynabeads（赛默飞世尔）进行净化，具体操作见用户指南。纯化的 cDNA 扩增 12 个循环，然后使用 SPRIselect beads（贝克曼）进行净化。按照单细胞 试剂盒 v3 用户指南的要求制备 cDNA 文库，并根据计算出的 cDNA 浓度对 PCR 循环进行适当修改（根据 Genomics 的建议）。
scRNA-seq data analysis
Quality control and doublet removal scRNA-seq数据分析 质量控制和双倍体去除
Raw sequencing data were aligned to the pig gene expression reference Sscrofa11 (the reference database was prepared with the gene annotation gff3 dataset v11.1.98 and genomic fasta dataset v11.1.98, which were downloaded from Ensemble). The CellRanger (10× Genomics) analysis pipeline was used to generate a digital gene expression matrix from these data. The raw digital gene expression matrix (UMI counts per gene per cell) was filtered, normalized, and clustered using R [v3.5.2, https://www.R-project.org/]. Cell and gene filtering was performed as follows: cells with fewer than 200 genes detected, or more than 30,000 UMIs and cells containing greater than of reads from mitochondrial genes were removed. To eliminate doublets, cells with more than 40,000 transcripts were identified using DoubletFinder [71] and then removed from the dataset. Genes detected (UMI count ) in fewer than three cells were removed. A total of 70,201 cells were included in downstream analysis, with an average of 1892.47 genes per cell and 5276.32 UMIs per cell. 原始测序数据与猪基因表达参考文献Sscrofa11（ 参考文献数据库是用基因注释gff3数据集v11.1.98和基因组fasta数据集v11.1.98准备的，这两个数据集是从Ensemble下载的）进行比对。CellRanger (10× Genomics)分析管道用于从这些数据中生成数字基因表达矩阵。原始数字基因表达矩阵（每个细胞每个基因的 UMI 计数）使用 R [v3.5.2, https://www.R-project.org/] 进行过滤、归一化和聚类。细胞和基因过滤方法如下：检测到的基因少于 200 个或 UMI 超过 30,000 个的细胞以及线粒体基因读数超过 的细胞被剔除。为了剔除双转录本，使用 DoubletFinder [71]识别出转录本超过 40,000 个的细胞，然后从数据集中剔除。在少于三个细胞中检测到的基因（UMI计数 ）也被剔除。共有 70 201 个细胞被纳入下游分析，平均每个细胞有 1892.47 个基因，每个细胞有 5276.32 个 UMI。
Cell clustering and cell-type annotation 细胞聚类和细胞类型注释
Normalization was performed in the Seurat R package (v3.1.1) using the default parameters [72]. Multiple samples were integrated by the CCA-based anchor method [73]. After principal component analysis (PCA), the first 30 principal components were selected for clustering the cells using standard package procedures. The Louvain algorithm with a resolution of 0.6 was used to cluster cells, which resulted in 31(29) distinct cell clusters. A gene was considered to be differentially expressed if it was detected in at least of one group and with at least 0.25 log fold change between two groups, and a Benjamini-Hochberg ( BH adjusted value in Wilcoxon rank-sum test was considered to indicate significance. To investigate the heterogeneity within the myogenic cells, we conducted sub-clustering analysis of myogenic cells. Using recently reported marker genes, we identified 4 myogenic sub cell types, including Pax3+ progenitors (Pax3P), myogenic progenitors (MyoP), myoblasts (Myob), and myocytes (Myoc). 归一化在 Seurat R 软件包（v3.1.1）中使用默认参数进行[72]。采用基于 CCA 的锚方法对多个样本进行整合[73]。主成分分析（PCA）后，使用标准软件包程序选择前 30 个主成分对细胞进行聚类。使用分辨率为 0.6 的卢万算法对细胞进行聚类，得出 31（29）个不同的细胞群。如果一个基因在一个组中至少 被检测到，并且在两个组之间至少有0.25对折的变化，则该基因被认为是差异表达的，本杰明-霍奇伯格（Benjamini-Hochberg）（BH 调整 值 在Wilcoxon秩和检验中被认为是显著的。为了研究肌原细胞内部的异质性，我们对肌原细胞进行了亚聚类分析。利用最近报道的标记基因，我们确定了4种肌原亚细胞类型，包括Pax3+祖细胞（Pax3P）、肌原祖细胞（MyoP）、肌母细胞（Myob）和肌细胞（Myoc）。
scRNA-seq trajectory analysis scRNA-seq 轨迹分析
Monocle 2 (version 2.10.1) [74] was used to construct the pseudotemporal path of myogenic cell differentiation. The myogenic cells were ordered in pseudotime along a trajectory using reduceDimension with the DDRTree method and orderCells functions. We detected genes that followed similar kinetic trends along the myogenic cell trajectory from the starting state. Hierarchical clustering was applied to cluster genes into five subgroups according to the expression patterns. A gene expression heatmap following pseudotime was plotted by the plot-pseudotimeheatmap function. Gene expression curves following pseudotime were calculated and plotted according to the mean value of gene expression of each pseudotime unit. The -axis (pseudotime, scaled) is the conversion value mapped from the pseudotime value of the cell in the time trajectory analysis result to the interval . Monocle 2（2.10.1 版）[74] 被用来构建成肌细胞分化的伪时间路径。使用 DDRTree 方法和 orderCells 函数 reduceDimension，沿着轨迹对肌原细胞进行伪时间排序。我们检测了从起始状态开始沿着成肌细胞轨迹遵循相似动力学趋势的基因。根据表达模式，我们对基因进行了分层聚类，将其分为五个亚组。使用 plot-pseudotimeheatmap 函数绘制了伪时间后的基因表达热图。根据每个伪时间单位的基因表达平均值计算并绘制伪时间后的基因表达曲线。 轴（伪时间，缩放）是从时间轨迹分析结果中细胞的伪时间值映射到 区间的换算值。
Gene regulatory network inference 基因调控网络推断
To identify TFs and characterize cell states, we employed cis-regulatory analysis using the R package SCENIC v1.1.2-2 with human orthologous genes; this program infers gene regulatory networks based 为了识别TFs并描述细胞状态，我们使用R软件包SCENIC v1.1.2-2对人类同源基因进行顺式调控分析。
on coexpression and DNA motif analysis. The network activity is then analyzed in each cell to identify recurrent cellular states. In short, TFs were identified using GENIE3 (v1.4.3) and compiled into modules (regulons), which were subsequently subjected to cis-regulatory motif analysis using RcisTarget (v1.2.1) with two gene-motif rankings: 10 kb around the TSS and 500 bp upstream. The regulon activity in every cell was then scored using AUCell (v1.4.1). Finally, binarized regulon activity was projected into t-distributed stochastic neighbor embedding (t-SNE) plots. 共表达和 DNA 主题分析。然后对每个细胞中的网络活动进行分析，以确定反复出现的细胞状态。简而言之，使用 GENIE3（v1.4.3）识别 TFs 并将其编译成模块（调控子），随后使用 RcisTarget（v1.2.1）对其进行顺式调控基因图谱分析，并进行两次基因图谱排序：在 TSS 周围 10 kb 和上游 500 bp 处进行。然后使用 AUCell（v1.4.1）对每个细胞中的调控子活性进行评分。最后，将二值化的调节子活性投射到 t 分布随机邻接嵌入（t-SNE）图中。
Ligand-receptor interactions 配体与受体的相互作用
To assess cellular crosstalk between different cell types, we used human orthologous genes converted UMI count data and the Python package CellPhoneDB v2.1.2 with the human database v2.0.0 to infer cell-cell communication networks from scRNA-seq data as per the authors' instructions. Only interactions that met the values threshold ( value ) were considered. 为了评估不同细胞类型之间的细胞串扰，我们按照作者的说明，使用人类同源基因转换的 UMI 计数数据和 Python 软件包 CellPhoneDB v2.1.2 及人类数据库 v2.0.0，从 scRNA-seq 数据中推断细胞-细胞通讯网络。只有符合 值阈值（ 值 ）的相互作用才会被考虑。
scATAC-seq was performed using the Chromium platform. All protocols to generate scATAC-seq data, including sample preparation, library preparation, and instrument and sequencing settings, are described below and are available at https://support.10xgenomics.com/single-cell-atac. scATAC-seq 使用 Chromium 平台进行。生成 scATAC-seq 数据的所有规程，包括样本制备、文库制备、仪器和测序设置，均在下文描述，可在 https://support.10xgenomics.com/single-cell-atac 上获取。
Cell filtering and cell clustering 细胞过滤和细胞聚类
Raw sequencing data aligned to the pig genome ATAC reference Sscrofa11 ( reference database was prepared with gene annotation gff3 data v11.1.98 and genomic fasta data v11.1.98 downloaded from the Ensemble). The digital peaks matrix was quantified using the Cell Ranger ATAC ( Genomics, v1.2.0) analysis pipeline. The following analysis was performed with the Seurat R package (v3.1.1) and Signac R package (v0.2.2). We kept only those valid barcodes that passed the following thresholds: (1) number of fragments ranging from 2000 to 30,000 , (2) mitochondria ratio less than , (3) transcription start site (TSS) enrichment score (by TSSEnrichment function) between 2 and 40, (4) strength of the nucleosome signal per cell (by NucleosomeSignal function) less than 5, and (5) fraction of reads overlapping peaks per cell greater than . After this stringent quality control, we obtained 48,514 single cells. Normalization and latent semantic indexing (LSI) dimensionality reduction were performed for each sample independently with the RunTFIDF and RunSVD functions. Batch correction by samples was applied by Signac on the first 50 LSI components. We applied the UMAP and TSNE algorithms to the first 50 原始测序数据与猪基因组ATAC参考Sscrofa11（ 参考数据库是用从Ensemble下载的基因注释gff3数据v11.1.98和基因组fasta数据v11.1.98）进行比对）。使用 Cell Ranger ATAC ( Genomics, v1.2.0) 分析管道对数字峰矩阵进行量化。以下分析使用 Seurat R 软件包（v3.1.1）和 Signac R 软件包（v0.2.2）进行。我们只保留通过以下阈值的有效条形码：（1）片段数在 2000 到 30000 之间；（2）线粒体比率小于 ；（3）转录起始位点（TSS）富集得分（通过 TSSEnrichment 函数）在 2 到 40 之间；（4）每个细胞的核小体信号强度（通过 NucleosomeSignal 函数）小于 5；（5）每个细胞的峰重叠读数分数大于 。经过如此严格的质量控制，我们获得了 48,514 个单细胞。使用 RunTFIDF 和 RunSVD 函数对每个样本独立进行了归一化和潜在语义索引（LSI）降维处理。Signac 对前 50 个 LSI 成分按样本进行了批量校正。我们对前 50 个样本应用了 UMAP 和 TSNE 算法。
LSI components corrected by Signac. We identified 15 distinct clusters using the Seurat function FindClusters (resolution equal to 0.3 ). LSI 成分经 Signac 校正。我们使用 Seurat 函数 FindClusters（分辨率等于 0.3）确定了 15 个不同的聚类。
Genomic element stratification 基因组元素分层
Pig ssc11 genome annotation files were downloaded from Ensembl (http://ftp.ensembl.org/pub/release-104/ gtf/sus_scrofa/Sus_scrofa. -Sscrofa11.1.104.gtf.gz). The 2-kb regions upstream of each TSS were included as promoter regions. We considered only open chromatin peaks detected in more than of myogenic cells. Since one open chromatin region could overlap with multiple genomic elements, we defined the priority order for the genomic elements as follows: exon promoter -UTR UTR intron distal elements. For example, if one peak overlapped with both an exon and a -UTR, the algorithm counted it as an exon-region peak. 猪 ssc11 基因组注释文件从 Ensembl 下载 ( http://ftp.ensembl.org/pub/release-104/ gtf/sus_scrofa/Sus_scrofa。-Sscrofa11.1.104.gtf.gz）。每个 TSS 上游的 2 kb 区域被列为启动子区域。我们只考虑在 以上的肌原细胞中检测到的开放染色质峰。由于一个开放染色质区域可能与多个基因组元素重叠，我们将基因组元素的优先顺序定义如下：外显子 启动子 -UTR UTR 内含子 远端元素。例如，如果一个峰同时与一个外显子和一个 -UTR 重叠，算法就把它算作一个外显子区域峰。
Identification of differentially accessible peaks 差异可及峰的识别
Differentially accessible peaks were defined as peaks detected in at least of cells in one group and exhibited a log fold change of at least 0.25 between two groups. The significance level was set as a Benjamini-Hochberg (BH) adjusted value in the logistic regression framework by the FindMarkers function with the parameter latent.vars = 'peak_region_fragments'. 差异可及峰的定义是：在一个组中至少有 个细胞检测到峰，且两组间的对数折叠变化至少为 0.25。在逻辑回归框架中，通过参数 latent.vars = 'peak_region_fragments' 的 FindMarkers 函数将显著性水平设置为 Benjamini-Hochberg (BH) 调整后的 值 。
scATAC trajectory analysis scATAC 轨迹分析
The imputed expression data of scATAC-seq were generated by integrating scRNA-seq and scATAC-seq data with the TransferData function. Monocle 2 was used to construct a pseudotemporal path from the imputed expression data. 利用 TransferData 函数整合 scRNA-seq 和 scATAC-seq 数据，生成了 scATAC-seq 的估算表达数据。使用 Monocle 2 从估算的表达数据中构建伪时间路径。
Predict cis-regulatory elements 预测顺式调节元件
We implemented two methods to study cis-regulatory elements in the snATAC-seq data. The first method was direct prediction. Peaks that overlapped with the TSSs of specific genes were assigned to cis-regulatory elements to those genes. The second method was the supplementary method, which followed Miao et al. There was coenrichment of the snATAC-seq cell-type-specific peaks and scRNA-seq cell-type-specific genes in the genome. This method links a gene with a peak if (1) the gene and peak were both specific to the same cell type, (2) the gene and peak were in cis, meaning that the peak was within region of the TSS of the corresponding gene, and (3) the peak did not directly overlap with the TSS of the gene. 我们采用了两种方法来研究 snATAC-seq 数据中的顺式调控元件。第一种方法是直接预测。将与特定基因的 TSS 重叠的峰归入这些基因的顺式调控元件。snATAC-seq细胞特异性峰与scRNA-seq细胞特异性基因在基因组中存在共富集。这种方法将基因与峰值联系起来，条件是：（1）基因和峰值都是同一细胞类型的特异基因；（2）基因和峰值是顺式的，即峰值在相应基因的 TSS 区域内；（3）峰值与基因的 TSS 没有直接重叠。
Correlation of scATAC-seq and scRNA-seq data scATAC-seq 和 scRNA-seq 数据的相关性
The consistency of myogenic cell subtypes defined in scRNA and scATAC was evaluated with correlation 通过相关性评估了 scRNA 和 scATAC 中定义的肌原细胞亚型的一致性。
analysis. Briefly, we summed the peaks intersecting the gene body and the region 2 kb upstream to calculate a gene activity score for each gene in each cell. This procedure is implemented in the Create Gene Activity Matrix in Seurat v3. Subsequently, the average gene activity score for scATAC and the average gene expression for scRNA were calculated for each cell type. The Pearson correlations of these two average values of scRNA-seq variable genes were computed for each cell type between the two data types. 分析。简而言之，我们将与基因体和上游 2 kb 区域相交的峰值相加，计算出每个细胞中每个基因的基因活性得分。Seurat v3 中的 "创建基因活性矩阵 "实现了这一程序。随后，计算出每种细胞中 scATAC 的平均基因活性得分和 scRNA 的平均基因表达量。为每种细胞类型计算这两种数据类型之间 scRNA-seq 可变基因的这两个平均值的皮尔逊相关性。
Integrated analysis of scATAC-seq and scRNA-seq data 综合分析 scATAC-seq 和 scRNA-seq 数据
scRNA-seq and scATAC-seq data from myogenic cells were integrated using the transfer procedure implemented in Signac v0.2.5 and Seurat v3.2.0. Briefly, the expression levels were imputed from the scATAC-seq data by defining a genomic region for each gene, which included the gene body and 2 kb upstream of the transcription start site, and taking the sum of scATAC-seq fragments within that region (genes from chromosomes , and M were excluded). Anchors were identified for scRNA- and scATAC-seq using the Find Transfer Anchors function with the canonical correlation analysis (CCA) reduction method with the scRNA-seq dataset as a reference. Then, scATAC-imputed gene expression values were projected on the scRNA expression values and weighted by LSI reduction by TransferData. Merged scATAC-seq imputed expression data and scRNA-seq data were then subjected to scaling, PCA reduction, and UMAP dimensional reduction to completely coembed the two types of data. 使用 Signac v0.2.5 和 Seurat v3.2.0 中的转移程序整合了肌原细胞的 scRNA-seq 和 scATAC-seq 数据。简言之，通过为每个基因定义一个基因组区域（包括基因体和转录起始位点上游 2 kb），并取该区域内的 scATAC-seq 片段之和（染色体 和 M 的基因不包括在内），从 scATAC-seq 数据中推算表达水平。以 scRNA-seq 数据集为参照，使用查找转移锚（Find Transfer Anchors）功能和典型相关分析（CCA）还原法确定了 scRNA 和 scATAC-seq 的锚。然后，将 scATAC 计算出的基因表达值投影到 scRNA 表达值上，并通过 TransferData 的 LSI 缩减进行加权。然后对合并的 scATAC-seq 估算表达数据和 scRNA-seq 数据进行缩放、PCA 缩减和 UMAP 降维处理，以完全合并两类数据。
Cell culture, proliferation, and differentiation 细胞培养、增殖和分化
Pig primary myogenic cells (PPMCs) were maintained in our laboratory. C 2 C 12 cells were purchased from the American Type Culture Collection (ATCC). The cells were cultured in Dulbecco's modified Eagle's medium (DMEM) with (v/v) fetal bovine serum (growth medium, GM). When the cells were cultured in GM at sub-confluent densities, it was defined as day 0 (d 0 ). To induce differentiation, cells were switched to DMEM with horse serum (differentiation medium, DM) after reaching confluence and then maintained in the culture medium for another , or 8 days (day 2 , day 3 , day 4 , day 5 , day 6 , or day 8 ). All cells were cultured in a incubator with . We examined cell proliferation using a CellLight EdU DNA Cell Proliferation Kit (RiboBio, China), as described previously [75]. 猪原代成肌细胞（PPMCs）由本实验室保存。C 2 C 12细胞购自美国类型培养物保藏中心（ATCC）。细胞在含有 （v/v）胎牛血清的杜氏改良老鹰培养基（DMEM）（生长培养基，GM）中培养。当细胞在 GM 中培养到亚融合密度时，将其定义为第 0 天（d 0）。为了诱导分化，细胞在达到 融合度后转入含有 马血清的DMEM（分化培养基，DM），然后在培养基中再维持 或8天（第2天、第3天、第4天、第5天、第6天或第8天）。所有细胞都是在 培养箱中用 培养的。我们使用 CellLight EdU DNA 细胞增殖试剂盒（RiboBio，中国）检测了细胞增殖情况，如前所述 [75]。
Gene overexpression 基因过表达
For the EGR1 and RHOB expression vectors, the coding sequences (CDSs) of the pig EGR1 and RHOB genes were inserted into the pcDNA3.1 vector (Invitrogen, Carlsbad, USA). Pig primary myogenic cells and C 2 C 12 cells were seeded into 6 - or 12 -well plates at 12 h before treatment and then transfected with expression plasmids using Lipofectamine 3000 (Invitrogen, Carlsbad, USA). Transfections were performed in at least triplicate for each experiment. 对于 EGR1 和 RHOB 表达载体，将猪 EGR1 和 RHOB 基因的编码序列（CDS）插入 pcDNA3.1 载体（Invitrogen，Carlsbad，USA）。在处理前 12 小时将猪原代成肌细胞和 C 2 C 12 细胞播种到 6 孔或 12 孔平板中，然后用 Lipofectamine 3000（Invitrogen 公司，美国卡尔斯巴德）转染表达质粒。每个实验至少转染三份。
Total RNA was extracted from cultured cells according to the instructions of TRIzol Reagent (Invitrogen, Shanghai, China), and cDNA was synthesized from of total RNA using a reverse-transcription kit (Genestar, Beijing, China). Real-time quantitative PCR (qPCR) was performed using a SYBR Green qPCR Kit (Genestar, Beijing, China) and detected on a LightCycler 480 II system (Roche, Basel, Switzerland). The primers used for qPCR are given in Additional file 13: Table S12. Gapdh was used as an internal control, and all reactions were run in triplicate. 按照TRIzol 试剂（Invitrogen公司，中国上海）的说明从培养细胞中提取总RNA，并用反转录试剂盒（Genestar公司，中国北京）从 总RNA中合成cDNA。使用 SYBR Green qPCR 试剂盒（Genestar，中国北京）进行实时定量 PCR（qPCR），并在 LightCycler 480 II 系统（罗氏，瑞士巴塞尔）上检测。用于 qPCR 的引物见附加文件 13：表 S12。Gapdh 用作内部对照，所有反应均一式三份。
Immunofluorescence
Cells cultured in 12 -well plates were fixed in paraformaldehyde for 10 min , followed by permeabilization in Triton X-100 for min. The cells were blocked with BSA in Tris-buffered saline with Tween (TBST) for 1 h . Then, the cells were incubated with primary antibodies overnight at . Afterwards, the cells were washed in PBS three times and incubated with secondary antibodies for 1 h at room temperature. Finally, the cells were washed three times in PBS, and the nuclei were counterstained with -diamidino2-phenylindole (DAPI; 1:1000 in PBS). The antibodies used were as follows: anti-MyoD antibody (Abcam, ab212662), anti-Pax7 antibody (Thermo Fisher Scientific, PA1117), anti-Fast Myosin Skeletal Heavy chain antibody (Abcam, ab51263), anti-mouse (H+L), F (ab')2 fragment (Alexa Fluor 555 conjugate) (Cell Signaling, #4409), and anti-rabbit (ab')2 fragment (Alexa Fluor 488 conjugate) (Cell Signaling, #4412). Immunostaining images were obtained via fluorescence microscopy on an inverted microscope (Nikon, Tokyo, Japan). 在 12 孔板中培养的细胞在 多聚甲醛中固定 10 分钟，然后在 Triton X-100 中渗透 分钟。细胞在含吐温的三缓冲盐水（TBST）中用 BSA 封闭 1 小时。然后，在 条件下用一抗孵育细胞过夜。然后，用 PBS 冲洗细胞三次，并与二抗在室温下孵育 1 小时。最后，用 PBS 冲洗细胞三次，并用 -二脒基 2-苯基吲哚（DAPI；1:1000，PBS 中）对细胞核进行反染色。使用的抗体如下抗MyoD抗体（Abcam，ab212662）、抗Pax7抗体（Thermo Fisher Scientific，PA1117）、抗快速肌球蛋白骨骼重链抗体（Abcam，ab51263）、抗小鼠 （H+L）、F (ab')2 片段（Alexa Fluor 555 连接物）（Cell Signaling，#4409）和抗兔 (ab')2 片段（Alexa Fluor 488 连接物）（Cell Signaling，#4412）。通过倒置显微镜（尼康，日本东京）上的荧光显微镜获得免疫染色图像。
Immunohistochemistry 免疫组化
Somite tissues of ZZ and DZ at E21 and E28 were fixed in paraformaldehyde for 19 h at , then dehydrated using gradient alcohol and embedded with paraffin. Paraffin-embedded samples were cut into somite cross sections. The paraffin sections were placed in an oven at for 30 min and immediately moved to xylene for dewaxing. The sections were rehydrated in gradient 将ZZ和DZ在E21和E28时的体节组织在 多聚甲醛中固定19小时，然后用梯度酒精脱水并用石蜡包埋。将石蜡包埋的样本切成 体节横切面。将石蜡切片放入 烘箱中烘烤30分钟，然后立即转移到二甲苯中脱蜡。切片在梯度
alcohol, and antigen retrieval was performed using citrate antigen retrieval solution. Finally, immunofluorescence staining was performed using IHC kit (Abcam, Cambridge, UK) according to the manufacturer's instructions. The antibodies used were as follows: anti-MyoD antibody (Cell Signaling Technology, 13812S), anti-Pax7antibody (Abcam, ab199010), Anti-rabbit IgG (H+L), F (ab')2 Fragment (Alexa Fluor Conjugate) (CellSignaling Technology, #4412), and Anti-mouse IgG (H+L), F (ab')2 Fragment (Alexa Fluor Conjugate) (CellSignaling Technology, #4409). Immunostaining images were obtained with a fluorescence microscope (Nikon, Tokyo, Japan). 然后用柠檬酸抗原回收液进行抗原回收。最后，根据生产商的说明使用 IHC 试剂盒（Abcam，英国剑桥）进行免疫荧光染色。所用抗体如下抗 MyoD 抗体（Cell Signaling Technology，13812S）、抗 Pax7 抗体（Abcam，ab199010）、抗兔 IgG（H+L），F（ab'）2 片段（Alexa Fluor 结合物）（CellSignaling Technology、#4412）和抗小鼠 IgG (H+L)、F (ab')2 片段（Alexa Fluor 结合物）（CellSignaling Technology，#4409）。用荧光显微镜（尼康，日本东京）获得免疫染色图像。
Statistical information 统计信息
The differential accessibility analysis was conducted by logistic regression. The differential expression analysis was conducted using the Wilcoxon rank sum test. Motif enrichment was determined based on the modified score (https://support.10xgenomics.com/singlecell-atac/software/pipelines/latest/algorithms/overview# zscore). All statistical tests were corrected for multiple testing using the Bonferroni method, and a significance level of 0.05 was used throughout the manuscript. 差异可及性分析采用逻辑回归法。差异表达分析采用 Wilcoxon 秩和检验。根据修正的 得分（https://support.10xgenomics.com/singlecell-atac/software/pipelines/latest/algorithms/overview# zscore）确定分子富集度。所有统计检验均采用 Bonferroni 方法进行多重检验校正，稿件中的显著性水平均为 0.05。
Supplementary Information 补充信息
The online version contains supplementary material available at https://doi. org/10.1186/s12915-023-01519-z. 在线版本包含补充材料，可查阅 https://doi. org/10.1186/s12915-023-01519-z。
Additional file 1: Figure S1. Quality control and batch effect correction in scRNA-Seq, related to Figure 1 A. Violin plots showing the number of expressed genes, the number of reads uniquely mapped against the reference genome, and the fraction of mitochondrial genes compared to all genes per cell in scRNA-Seq data. B. Box plot showing the number of genes (left) and the number of uniquely mapped reads (right) per cell in each identified cell type in scRNA-Seq data. C. tSNE plot visualization of the sample source for all 70,201 cells. Each dot is a cell. Different colors represent different samples. D. tSNE plot visualization of unsupervised clustering analysis for all 70,201 cells based on scRNA-Seq data after quality control, which gave rise to 31 distinct clusters. Figure S2. Gene Ontology (GO) analysis of the DEGs for each cell type was performed and the representative enriched GO terms are presented, related to Figure 1. Figure S3. Expression of selected marker genes along the differentiation trajectory, related to Figure 2 A. tSNE plot demonstrating cell cycle regression (left). Visualization of myogenic differentiation trajectory by cell cycle phases (G1, S, and G2/M) (right). B. Donut plots showing the percentages of cells in G1, S, and G2M phase at different cell states. C. Expression levels of cell cycle-related genes in the myogenic cells organized into the Monocle trajectory. D. Expression levels of muscle related genes in the myogenic cells organized into the Monocle trajectory. Figure S4. Unsupervised clustering analysis for all cells in scATAC-Seq data and myogenicspecific scATAC-seq peaks, related to Figure 4 A-C. tSNE plot visualization of the sample source for all 48514 cells in scATAC-Seq. Each dot is a cell. Different colors represent different pigs (A), different embryonic stages (B), or different samples (C). D. tSNE plot visualization of unsupervised clustering analysis for all 48514 cells after quality control in scATAC-Seq data, which gave rise to 15 distinct clusters. E. tSNE plot visualization of myogenic cells and other cells. Clusters 4 and 8 in Figure S4D were annotated as myogenic cells due to their high levels of accessibility of marker genes associated with myogenic lineage. F. Genome browser view of myogenic-specific peaks at the TSS of MyoG and Myf5 for myogenic cells and other cells in the scATAC-seq dataset. Figure S5. Percentage distribution of open chromatin elements in scATAC-Seq data, related to Figure 4 A. Distribution of open chromatin elements in each snATAC-seq sample. B. Distribution of open chromatin elements in snATAC-seq of myogenic cell types. C. Percentage distribution of open chromatin elements among DAPs in myogenic cell types. Figure S6. Integrative analysis of transcription factors and target genes, related to Figure 5 A. tSNE depiction of regulon activity ("on-blue", "off-gray"), TF gene expression (red scale), and expression of predicted target genes (purple scale) of MyoG, FOSB, and TCF12. B. Corresponding chromatin accessibility in sCATAC data for TFs and predicted target genes are depicted. Figure S7. Pseudotimedependent chromatin accessibility and gene expression changes, related to Figure 7. The first column shows the dynamics of the Genomics TF enrichment score. The second column shows the dynamics of TF gene expression values, and the third and fourth columns represent the dynamics of the SCENIC-reported target gene expression values of corresponding TFs, respectively. Figure S8. Myogenesis related gene expression in DMD (Duchenne muscular dystrophy) mice. Comparison of RNA-seq data of flexor digitorum short (FDB), extensor digitorum long (EDL), and soleus (SOL) in DMD and wild-type mice including 2-month and 5-month age. A. The expression levels of myogenesis related genes (Myod1, Myog, Myf5, Pax7). B. The expression levels of related genes that were upregulated during porcine embryonic myogenesis (EGR1, RHOB, KLF4, SOX8, NGFR, MAX, RBFOX2, ANXA6, HES6, RASSF4, PLS3, SPG21). C. The expression levels of related genes that were downregulated during porcine embryonic myogenesis COX5A, HOMER2, BNIP3, CNCS). Data were obtained from the GEO database (GSE162455;WT, ). Figure S9. Genome browser view of differentially accessible peaks at the TSS of EGR1 and between myogenic cells in the sCATAC-seq dataset, related to Figure 8. Figure S10. Functional analysis of EGR1 in myogenesis, related to Figure 8 A-B. EdU assays for the proliferation of pig primary myogenic cells (A) and C2C12 myoblasts following EGR1 overexpression. C. qPCR analysis of the mRNA levels of cell cycle regulators in C2C12 cells following EGR1 附加文件 1：图 S1.与图 1 A 有关的 scRNA-Seq 质量控制和批次效应校正。Violin 图显示表达基因的数量、与参考基因组唯一映射的读数数量，以及与 scRNA-Seq 数据中每个细胞的所有基因相比线粒体基因的比例。B. 方框图显示了 scRNA-Seq 数据中每种已识别细胞类型中每个细胞的基因数量（左）和唯一映射读数数量（右）。C. 所有 70 201 个细胞样本源的 tSNE 图。每个点代表一个细胞。不同颜色代表不同样本。D. 基于质量控制后的 scRNA-Seq 数据对所有 70 201 个细胞进行无监督聚类分析的 tSNE 图，结果显示有 31 个不同的聚类。图 S2.对每种细胞类型的 DEGs 进行了基因本体（GO）分析，并展示了与图 1 相关的具有代表性的富集 GO 术语。图 S3.分化轨迹上选定标记基因的表达，与图 2 有关 A. 显示细胞周期回归的 tSNE 图（左）。按细胞周期阶段（G1、S 和 G2/M）显示的成肌细胞分化轨迹（右图）。B. 多纳图显示不同细胞状态下 G1、S 和 G2M 期细胞的百分比。C. 成肌细胞中细胞周期相关基因的表达水平按 Monocle 轨迹排列。D. 按 Monocle 轨迹排列的肌原细胞中肌肉相关基因的表达水平。图 S4.对 scATAC-Seq 数据中所有细胞和肌原特异性 scATAC-seq 峰的无监督聚类分析，与图 4 A-C 相关。每个点代表一个细胞。 不同颜色代表不同的猪（A）、不同的胚胎阶段（B）或不同的样本（C）。D. 对 scATAC-Seq 数据进行质量控制后，对所有 48514 个细胞进行无监督聚类分析的 tSNE 图，结果显示有 15 个不同的聚类。E. 肌原细胞和其他细胞的 tSNE 图。图 S4D 中的第 4 和第 8 个簇被注释为肌原细胞，因为它们与肌原细胞系相关的标记基因的可及性很高。F. 在基因组浏览器中查看成肌细胞和 scATAC-seq 数据集中 MyoG 和 Myf5 的 TSS 处的成肌特异性峰值。图 S5.与图 4 A 有关的 scATAC-Seq 数据中开放染色质元素的百分比分布。每个 snATAC-seq 样本中开放染色质元素的分布。B. 成肌细胞类型的 snATAC-seq 中开放染色质元素的分布。C. 成肌细胞类型 DAPs 中开放染色质元素的百分比分布。图 S6.与图 5 有关的转录因子和靶基因的综合分析 A. tSNE 描述了 MyoG、FOSB 和 TCF12 的调控子活性（"蓝色"、"灰色"）、TF 基因表达（红色标度）和预测靶基因表达（紫色标度）。B. 在 sCATAC 数据中描绘了 TF 和预测靶基因的相应染色质可及性。图 S7.与图 7 相关的伪时间依赖性染色质可及性和基因表达变化。第一列显示了 Genomics TF 富集得分的动态。第二列显示了 TF 基因表达值的动态变化，第三列和第四列分别代表了 SCENIC 报告的相应 TF 靶基因表达值的动态变化。图 S8. DMD（杜氏肌营养不良症）小鼠肌生成相关基因的表达。比较 DMD 小鼠和野生型小鼠（包括 2 个月和 5 个月大）短屈肌（FDB）、长伸肌（EDL）和比目鱼肌（SOL）的 RNA-seq 数据。A. 肌生成相关基因（Myod1、Myog、Myf5、Pax7）的表达水平。B. 猪胚胎肌生成过程中上调的相关基因（EGR1、RHOB、KLF4、SOX8、NGFR、MAX、RBFOX2、ANXA6、HES6、RASSF4、PLS3、SPG21）的表达水平。C. 猪胚胎成肌过程中下调的相关基因的表达水平（COX5A、HOMER2、BNIP3、CNCS）。数据来自 GEO 数据库（GSE162455;WT, ）。图 S9.基因组浏览器查看 sCATAC-seq 数据集中肌原细胞之间 EGR1 和 的 TSS 处不同的可访问峰，与图 8 有关。图 S10.EGR1 在成肌过程中的功能分析，见图 8 A-B。过表达 EGR1 后猪原代成肌细胞（A）和 C2C12 成肌细胞增殖的 EdU 检测。C.对 EGR1 使用后 C2C12 细胞中细胞周期调节因子 mRNA 水平的 qPCR 分析
overexpression. D. Immunofluorescence staining for MyHC in C2C12 cells following EGR1 overexpression and differentiation for 3 d . Then, the fusion index was calculated. Figure S11. Functional analysis of RHOB in myogenesis, related to Figure 8 A-B. EdU assays for proliferation of pig primary myogenic cells (A) and C2C12 myoblasts following RHOB overexpression. C. qPCR analysis of the mRNA levels of cell-cycle regulators in C2C12 cells following RHOB overexpression. D. Immunofluorescence staining for MyHC in C2C12 cells following RHOB overexpression and differentiation for 3 d . Then, the fusion index was calculated. 过表达。D. EGR1 过表达并分化 3 d 后 C2C12 细胞中 MyHC 的免疫荧光染色。然后计算融合指数。图 S11.RHOB 在成肌过程中的功能分析，与图 8 A-B 相关。猪原代成肌细胞（A）和 C2C12 成肌细胞在过表达 RHOB 后的增殖的 EdU 检测。C. RHOB 过表达后 C2C12 细胞中细胞周期调节因子 mRNA 水平的 qPCR 分析。D. RHOB 过表达并分化 3 d 后 C2C12 细胞中 MyHC 的免疫荧光染色。然后计算融合指数。
Additional file 2: Supplementary Table S1. DEG between 31 distinct cell clusters in scRNA-seq data after ambient RNA cleaning. Related to Figure 1. 附加文件 2：补充表 S1。环境 RNA 清除后 scRNA-seq 数据中 31 个不同细胞群之间的 DEG。与图 1 有关。
Additional file 3: Supplementary Table S2. Gene Ontology analysis of the DEGs for each cell type. Related to Figure S2. 附加文件 3：补充表 S2。各细胞类型 DEGs 的基因本体分析。与图 S2 有关。
Additional file 4: Supplementary Table S3. The 1,700 top DEGs among 4 myogenic cell types were analyzed and clustered into five major categories of transcriptional gene clusters according to their expression changes. Related to Figure 3. 附加文件 4：补充表 S3。分析了 4 种肌原细胞类型中的 1,700 个顶级 DEGs，并根据其表达变化将其聚类为五大类转录基因簇。与图 3 有关。
Additional file 5: Supplementary Table S4. Differentially accessible peaks between 15 distinct clusters in scATAC-seq data. Related to Figure 4 附加文件 5：补充表 S4。scATAC-seq 数据中 15 个不同聚类之间的不同可访问峰。与图 4 相关
Additional file 6: Supplementary Table S5. Differentially accessible peaks between 7 distinct sub-clusters of myogenic cells in ScATAC-seq data. Related to Figure 4. 附加文件 6：补充表 S5。ScATAC-seq 数据中 7 个不同的成肌细胞亚群之间可访问的不同峰值。与图 4 有关。
Additional file 7: Supplementary Table S6. Cell type-specific open chromatin derived from scATAC-seq analysis. Related to Figure 4. 附加文件 7：补充表 S6。通过 scATAC-seq 分析获得的细胞类型特异性开放染色质。与图 4 有关。
Additional file 8: Supplementary Table S7. Full list of cell type-specific motif enrichment. Related to Figure 5. 附加文件 8：补充表 S7。细胞类型特异性图案富集的完整列表。与图 5 有关。
Additional file 9: Supplementary Table S8. Expression of transcription factors associated with cell type-specific motif enrichment. Related to Figure 5 . 附加文件 9：补充表 S8。与细胞类型特异性图案富集相关的转录因子的表达。与图 5 相关。
Additional file 10: Supplementary Table S9. The top 10 regulons and respective target genes inferred by SCENIC. Related to Figure 5. 附加文件 10：补充表 S9。SCENIC 推断的前 10 个调控子和各自的靶基因。与图 5 有关。
Additional file 11: Supplementary Table S10. The complete list of regulons and their respective predicted target genes inferred by SCENIC. Related to Figure 5. 附加文件 11：补充表 S10。由 SCENIC 推断的调控子及其各自预测的靶基因的完整列表。与图 5 有关。
Additional file 12: Supplementary Table S11. Scaled and binarized regulon activities in each cell type inferred by SCENIC. Related to Figure 5. 附加文件 12：补充表 S11。SCENIC 推断的每种细胞类型中按比例和二值化的调控子活动。与图 5 有关。
Additional file 13: Supplementary Table S12. The primers for qPCR. Related to Figure 8, Figure S10, and Figure S11. 附加文件 13：补充表 S12。用于 qPCR 的引物。与图 8、图 S10 和图 S11 有关。
Acknowledgements 致谢
Not applicable 不适用
Authors' contributions 作者的贡献
Shufang Cai, Bin Hu, Xiaoyu Wang, Tongni Liu, contributed equally to this work. D.M. designed and conceived the project. B.H. provided experimental pigs and arranged embryo collection. S.C., X.W., and T.L. analyzed the scRNAseq and scATAC-seq data. Z.L., R.X., and M.C. performed RNA extraction, qPCR and immunofluorescence assays. X.T., T.D., and Q.Z. helped in sample collection. Z.L. and E.L. helped in cell culture. Y.C., J.L., and X.L. contributed to data processing, discussions, and advice. S.C. and D.M. wrote the manuscript. All authors read and approved the final manuscript. 蔡淑芳、胡斌、王晓宇、刘彤妮对本研究做出了同等贡献。D.M. 设计并构思了该项目。B.H. 提供实验猪并安排胚胎采集。S.C.、X.W.和T.L.分析了scRNAseq和scATAC-seq数据。Z.L.、R.X.和M.C.进行了RNA提取、qPCR和免疫荧光检测。X.T.、T.D.和 Q.Z. 协助收集样本。Z.L. 和 E.L. 协助细胞培养。Y.C.、J.L.和X.L.参与了数据处理、讨论和建议。S.C.和D.M.撰写了手稿。所有作者阅读并批准了最终手稿。
Funding 资金筹措
This work was supported by the National Natural Science Foundation of China (32072697, 31772565), Guangdong Basic and Applied Basic Research Foundation (2019B1515210013), China Agriculture Research System (CASR-35). 本研究得到了国家自然科学基金（32072697、31772565）、广东省基础与应用基础研究基金（2019B1515210013）、中国农业科研系统（CASR-35）的资助。
Availability of data and materials 数据和材料的可用性
All data generated or analyzed during this study are included in this published article, its supplementary information files and publicly available repositories. Raw data and processed data of scRNA-seq and scATAC-seq have been deposited in GEO (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE206914) with the accession codes GSE206914 [76]. In this study, we downloaded public 本研究中生成或分析的所有数据都包含在这篇已发表的文章、其补充信息文件和可公开获取的资料库中。scRNA-seq和scATAC-seq的原始数据和处理数据已存入GEO（https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE206914），加入代码为GSE206914[76]。在本研究中，我们下载了公开的
RNA-Seq dataset of skeletal muscle from wild-type and Duchenne muscular dystrophy (DMD) mice with an accession number GSE162455 (https://wwwncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE162455) [77]. All codes used for data analysis and graph generation in this research could be obtained from https://github.com/wangxiaoyu0687/Integrative-single-cell-RNA-seq-and-ATACseq-analysis-demonstrates-cellular-myogenic-differentiation. 野生型和杜氏肌营养不良症（DMD）小鼠骨骼肌的 RNA-Seq 数据集，登录号为 GSE162455 ( https://www ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE162455)[77]。本研究中用于数据分析和图表生成的所有代码均可从 https://github.com/wangxiaoyu0687/Integrative-single-cell-RNA-seq-and-ATACseq-analysis-demonstrates-cellular-myogenic-differentiation 获取。
Declarations 声明
Ethics approval and consent to participate 伦理批准和参与同意书
The animal experimental procedures used in this experiment were approved by the Animal Care and Use Committee of Guangdong Province, China. Approval ID or permit numbers are SCXK (Guangdong) 2011-0029 and SYXK (Guangdong) 2011-0112. 本实验中使用的动物实验程序已获得中国广东省动物保育和使用委员会的批准。批准文号为SCXK（粤）2011-0029和SYXK（粤）2011-0112。
Consent for publication 同意出版
Not applicable. 不适用。
Competing interests 竞争利益
The authors declare that they have no competing interests. 作者声明他们没有利益冲突。
Received: 11 October 2022 Accepted: 18 January 2023 收到：接受： 2022 年 10 月 11 日2023 年 1 月 18 日
Published online: 01 February 2023 在线出版日期：2023 年 2 月 1 日
References 参考资料
Buckingham M, Rigby PWJ. Gene regulatory networks and transcriptional mechanisms that control myogenesis. Dev Cell. 2014;28:225-38. Buckingham M, Rigby PWJ.控制肌肉生成的基因调控网络和转录机制。Dev Cell.2014;28:225-38.
Bentzinger CF, Wang YX, Rudnicki MA. Building muscle: molecular regulation of myogenesis. Cold Spring Harb Perspect Biol. 2012;4(2):a008342. Bentzinger CF, Wang YX, Rudnicki MA.构建肌肉：肌肉生成的分子调控。Cold Spring Harb Perspect Biol.
Buckingham M, Montarras D, Relaix F, Cumano A, Mansouri A, Morgan J, et al. Pax3 and Pax7 mark a major population of muscle progenitor cells that contribute to skeletal muscle formation and regeneration. Neuromuscul Disord. 2006;16:S48-8 Buckingham M, Montarras D, Relaix F, Cumano A, Mansouri A, Morgan J, et al. Pax3 和 Pax7 标志着有助于骨骼肌形成和再生的主要肌肉祖细胞群。Neuromuscul Disord.2006;16:S48-8
Murani E, Muraniova M, Ponsuksili S, Schellander K, Wimmers K. Identification of genes differentially expressed during prenatal development of skeletal muscle in two pig breeds differing in muscularity. BMC Dev Biol. 2007;7:109. Murani E, Muraniova M, Ponsuksili S, Schellander K, Wimmers K. Identification of genes differentially expressed during prenatal development of skeletal muscle in two pig breed different in muscularity.BMC Dev Biol.
Picard B, Berri C, Lefaucheur L, Molette C, Sayd T, Terlouw C. Skeletal muscle proteomics in livestock production. Brief Funct Genomics. 2010;9:259-78. Picard B、Berri C、Lefaucheur L、Molette C、Sayd T、Terlouw C. 牲口生产中的骨骼肌蛋白质组学。Brief Funct Genomics.2010;9:259-78.
Buckingham M, Relaix F. PAX3 and PAX7 as upstream regulators of myogenesis. Semin Cell Dev Biol. 2015;44:115-25. Buckingham M, Relaix F.PAX3 和 PAX7 作为肌生成的上游调节因子。Semin Cell Dev Biol.
Hernandez-Hernandez M, Garcia-Gonzalez EG, Brun CE, Rudnicki MA. The myogenic regulatory factors, determinants of muscle development, cell identity and regeneration. Semin Cell Dev Biol. 2017;72:10-8. Hernandez-Hernandez M, Garcia-Gonzalez EG, Brun CE, Rudnicki MA.肌肉发育、细胞特性和再生的决定因素--肌肉生成调节因子。Semin Cell Dev Biol.
Chal J, Pourquie O. Making muscle: skeletal myogenesis in vivo and in vitro. Development. 2017;144:2104-22. Chal J, Pourquie O. 制造肌肉：体内和体外的骨骼肌生成。发育。2017;144:2104-22.
Zammit PS. Function of the myogenic regulatory factors Myf5, MyoD, Myogenin and MRF4 in skeletal muscle, satellite cells and regenerative myogenesis. Semin Cell Dev Biol. 2017;72:19-32. Zammit PS.骨骼肌、卫星细胞和再生肌生成中肌生成调控因子Myf5、MyoD、Myogenin和MRF4的功能。Semin Cell Dev Biol.
Cao Y, Yao Z, Sarkar D, Lawrence M, Sanchez GJ, Parker MH, et al. Genomewide MyoD binding in skeletal muscle cells: a potential for broad cellular reprogramming. Dev Cell. 2010;18:662-74. Cao Y, Yao Z, Sarkar D, Lawrence M, Sanchez GJ, Parker MH, et al. 骨骼肌细胞中的全基因组 MyoD 结合：广泛细胞重编程的潜力。Dev Cell.2010;18:662-74.
Asp P, Blum R, Vethantham V, Parisi F, Micsinai M, Cheng J, et al. Genomewide remodeling of the epigenetic landscape during myogenic differentiation. Proc Natl Acad Sci U S A. 2011;108:E149-58 Asp P, Blum R, Vethantham V, Parisi F, Micsinai M, Cheng J, et al. Myogenic 分化过程中表观遗传景观的全基因组重塑。Proc Natl Acad Sci U S A. 2011;108:E149-58
Fu Y, Shang P, Zhang B, Tian XL, Nie RX, Zhang R, et al. Function of the porcine TRPC1 gene in myogenesis and muscle growth. Cells. 2021;10(1):147 Fu Y, Shang P, Zhang B, Tian XL, Nie RX, Zhang R, et al. 猪 TRPC1 基因在肌生成和肌肉生长中的功能。细胞。2021;10(1):147
Groenen MAM, Archibald AL, Uenishi H, Tuggle CK, Takeuchi Y, Rothschild , et al. Analyses of pig genomes provide insight into porcine demography and evolution. Nature. 2012;491:393-8. Groenen MAM, Archibald AL, Uenishi H, Tuggle CK, Takeuchi Y, Rothschild , et al. 对猪基因组的分析为了解猪的人口和进化提供了线索。自然》。2012;491:393-8.
Huang Z, Li Q, Li M, Li C. Transcriptome analysis reveals the long intergenic noncoding RNAs contributed to skeletal muscle differences between Yorkshire and Tibetan pig. Sci Rep. 2021;11:2622. Huang Z, Li Q, Li M, Li C. 转录组分析揭示长基因间非编码RNA对约克夏猪和藏香猪骨骼肌差异的贡献科学报告，2021；11:2622.
Zhang B, Ban D, Gou X, Zhang Y, Yang L, Chamba Y, et al. Genome-wide DNA methylation profiles in Tibetan and Yorkshire pigs under highaltitude hypoxia. J Anim Sci Biotechnol. 2019;10:25. Zhang B, Ban D, Gou X, Zhang Y, Yang L, Chamba Y, et al. 高海拔缺氧条件下藏香猪和约克夏猪的全基因组DNA甲基化图谱。J Anim Sci Biotechnol.2019;10:25.
Yang Y, Huang Q, Xin H, Li B, Du Z, Meng F, et al. Analysis of differences in growth, meat quality traits and expression of genes related to fat deposition between Tibetan and Duzang pigs. J China South Agric. 2022;04:926-34 Yang Y, Huang Q, Xin H, Li B, Du Z, Meng F, et al. 藏香猪与堆藏猪生长、肉质性状及脂肪沉积相关基因表达的差异分析。J China South Agric.2022;04:926-34
Hwang B, Lee JH, Bang D. Single-cell RNA sequencing technologies and bioinformatics pipelines (vol 50, pg 1, 2018). Exp Mol Med. 2021;53:1005-5. Hwang B、Lee JH、Bang D. 单细胞RNA测序技术和生物信息学管道（第50卷，第1页，2018年）。Exp Mol Med.2021;53:1005-5.
Chen X, Miragaia RJ, Natarajan KN, Teichmann SA. A rapid and robust method for single cell chromatin accessibility profiling. Nat Commun. 2018;9:5345. Chen X, Miragaia RJ, Natarajan KN, Teichmann SA.快速、稳健的单细胞染色质可及性分析方法。Nat Commun.2018;9:5345.
Giordani L, He GJ, Negroni E, Sakai H, Law JYC, Siu MM, et al. High-dimensional single-cell cartography reveals novel skeletal muscle-resident cell populations. Mol Cell. 2019;74:609-621 e606 Giordani L, He GJ, Negroni E, Sakai H, Law JYC, Siu MM, et al. 高维单细胞制图揭示了新的骨骼肌驻留细胞群。Mol Cell.2019;74:609-621 e606
Dell'Orso S, Juan AH, Ko KD, Naz F, Perovanovic J, Gutierrez-Cruz G , et al. Single cell analysis of adult mouse skeletal muscle stem cells in homeostatic and regenerative conditions. Development. 2019;146(12):dev174177. Dell'Orso S, Juan AH, Ko KD, Naz F, Perovanovic J, Gutierrez-Cruz G , et al. 单细胞分析小鼠骨骼肌成体干细胞的稳态和再生状态。发育。2019;146(12):dev174177.
De Micheli AJ, Laurilliard EJ, Heinke CL, Ravichandran H, Fraczek P, SoueidBaumgarten S , et al. Single-cell analysis of the muscle stem cell hierarchy identifies heterotypic communication signals involved in skeletal muscle regeneration. Cell Rep. 2020;30:3583-3595 e3585. De Micheli AJ, Laurilliard EJ, Heinke CL, Ravichandran H, Fraczek P, SoueidBaumgarten S , et al. 肌肉干细胞层次结构的单细胞分析确定了参与骨骼肌再生的异型通讯信号。Cell Rep. 2020;30:3583-3595 e3585.
Xu Z, You W, Chen W, Zhou Y, Nong Q, Valencak TG, et al. Single-cell RNA sequencing and lipidomics reveal cell and lipid dynamics of fat infiltration in skeletal muscle. J Cachexia Sarcopenia Muscle. 2021;12:109-29. Xu Z, You W, Chen W, Zhou Y, Nong Q, Valencak TG, et al. 单细胞 RNA 测序和脂质组学揭示骨骼肌中脂肪浸润的细胞和脂质动态。J Cachexia Sarcopenia Muscle.2021;12:109-29.
Xi HB, Langerman J, Sabri S, Chien P, Young CS, Younesi S, et al. A human skeletal muscle atlas identifies the trajectories of stem and progenitor cells across development and from human pluripotent stem cells (vol 27, pg 158, 2020). Cell Stem Cell. 2020;27:181-5. Xi HB、Langerman J、Sabri S、Chien P、Young CS、Younesi S 等：《人类骨骼肌图谱》（A human skeletal muscle atlas identifies the trajectories of stem and progenitor cells across development and from human pluripotent stem cells）（第 27 卷，第 158 页，2020 年）。细胞干细胞》。2020;27:181-5.
Polanski K, Young MD, Miao Z, Meyer KB, Teichmann SA, Park JE. BBKNN: fast batch alignment of single cell transcriptomes. Bioinformatics. 2020;36:964-5.
Tangherloni A, Spolaor S, Rundo L, Nobile MS, Cazzaniga P, Mauri G, et al. GenHap: a novel computational method based on genetic algorithms for haplotype assembly. BMC Bioinformatics. 2019;20:172. Tangherloni A, Spolaor S, Rundo L, Nobile MS, Cazzaniga P, Mauri G, et al. GenHap：基于遗传算法的单倍型组装新计算方法。BMC Bioinformatics.2019;20:172.
Suay-Corredera C, Pricolo MR, Herrero-Galan E, Velazquez-Carreras D, Sanchez-Ortiz D, Garcia-Giustiniani D, et al. Protein haploinsufficiency drivers identify MYBPC3 variants that cause hypertrophic cardiomyopathy. J Biol Chem. 2021;297(1):100854. Suay-Corredera C、Pricolo MR、Herrero-Galan E、Velazquez-Carreras D、Sanchez-Ortiz D、Garcia-Giustiniani D 等蛋白单倍缺失驱动因子识别出导致肥厚型心肌病的 MYBPC3 变体。J Biol Chem.2021;297(1):100854.
Qiu XJ, Mao Q, Tang Y, Wang L, Chawla R, Pliner HA, et al. Reversed graph embedding resolves complex single-cell trajectories. Nat Methods. 2017;14:979-82. Qiu XJ, Mao Q, Tang Y, Wang L, Chawla R, Pliner HA, et al.Nat Methods.2017;14:979-82.
Frankish A, Diekhans M, Ferreira AM, Johnson R, Jungreis I, Loveland , et al. GENCODE reference annotation for the human and mouse genomes. Nucleic Acids Res. 2019;47:D766-73. Frankish A, Diekhans M, Ferreira AM, Johnson R, Jungreis I, Loveland , et al. GENCODE 人类和小鼠基因组参考注释。Nucleic Acids Res. 2019;47:D766-73.
Aibar S, Gonzalez-Blas CB, Moerman T, Huynh-Thu VA, Imrichova H, Hulselmans G, et al. SCENIC: single-cell regulatory network inference and clustering. Nat Methods. 2017;14:1083-6. Aibar S, Gonzalez-Blas CB, Moerman T, Huynh-Thu VA, Imrichova H, Hulselmans G, et al. SCENIC：单细胞调控网络推断与聚类。Nat Methods.2017;14:1083-6.
Korsunsky I, Millard N, Fan J, Slowikowski K, Zhang F, Wei K, et al. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat Methods. 2019;16:1289-96. Korsunsky I, Millard N, Fan J, Slowikowski K, Zhang F, Wei K, et al. 用 Harmony 快速、灵敏、准确地整合单细胞数据。Nat Methods.2019;16:1289-96.
Thumkeo D, Watanabe S, Narumiya S. Physiological roles of Rho and Rho effectors in mammals. Eur J Cell Biol. 2013;92:303-15. Thumkeo D, Watanabe S, Narumiya S. 哺乳动物中 Rho 和 Rho 效应体的生理作用。Eur J Cell Biol.
Sun GX, Song H, Wu SY. miR-19a promotes vascular smooth muscle cell proliferation, migration and invasion through regulation of Ras homolog family member B. Int J Mol Med. 2019;44:1991-2002. Sun GX, Song H, Wu SY. miR-19a通过调控Ras同源家族成员B促进血管平滑肌细胞增殖、迁移和侵袭 Int J Mol Med.2019;44:1991-2002.
Chan KC, Lin MC, Huang CN, Chang WC, Wang CJ. Mulberry 1-deoxynojirimycin pleiotropically inhibits glucose-stimulated vascular smooth muscle cell migration by activation of AMPK/RhoB and down-regulation of FAK. J Agric Food Chem. 2013;61:9867-75. Chan KC, Lin MC, Huang CN, Chang WC, Wang CJ.桑椹 1-脱氧野尻霉素通过激活 AMPK/RhoB 和下调 FAK 多效抑制葡萄糖刺激的血管平滑肌细胞迁移。J Agric Food Chem.2013;61:9867-75.
Vento-Tormo R, Efremova M, Botting RA, Turco MY, Vento-Tormo M, Meyer KB , et al. Single-cell reconstruction of the early maternal-fetal interface in humans. Nature. 2018;563:347-53 Vento-Tormo R、Efremova M、Botting RA、Turco MY、Vento-Tormo M、Meyer KB 等人早期母胎界面的单细胞重建。自然2018;563:347-53
Zanou N, Gailly P. Skeletal muscle hypertrophy and regeneration: interplay between the myogenic regulatory factors (MRFs) and insulin-like growth factors (IGFs) pathways. Cell Mol Life Sci. 2013;70:4117-30. Zanou N, Gailly P. 骨骼肌肥大与再生：肌生成调节因子（MRFs）与胰岛素样生长因子（IGFs）通路之间的相互作用。Cell Mol Life Sci.
Marics I, Padilla F, Guillemot JF, Scaal M, Marcelle C. FGFR4 signaling is a necessary step in limb muscle differentiation. Development. 2002;129:4559-69. Marics I, Padilla F, Guillemot JF, Scaal M, Marcelle C. FGFR4 信号传导是四肢肌肉分化的必要步骤。发育。2002;129:4559-69.
Van Ho AT, Hayashi S, Brohl D, Aurade F, Rattenbach R, Relaix F. Neural crest cell lineage restricts skeletal muscle progenitor cell differentiation through Neuregulin1-ErbB3 signaling. Dev Cell. 2011;21:273-87. Van Ho AT, Hayashi S, Brohl D, Aurade F, Rattenbach R, Relaix F.神经嵴细胞系通过 Neuregulin1-ErbB3 信号限制骨骼肌祖细胞分化。Dev Cell.2011;21:273-87.
Camerino GM, Pierno S, Liantonio A, De Bellis M, Cannone M, Sblendorio , et al. Effects of pleiotrophin overexpression on mouse skeletal muscles in normal loading and in actual and simulated microgravity. PLoS One. 2013;8:e72028 Camerino GM, Pierno S, Liantonio A, De Bellis M, Cannone M, Sblendorio , et al. 多养蛋白过表达对小鼠骨骼肌在正常负荷以及实际和模拟微重力条件下的影响。PLoS One.2013;8:e72028
Chen . Mitochondrial dynamics in regulating the unique phenotypes of cancer and stem cells. Cell Metab. 2017;26:39-48. 陈 .调控癌细胞和干细胞独特表型的线粒体动力学。细胞代谢。2017;26:39-48.
Shyh-Chang N, Daley GQ, Cantley LC. Stem cell metabolism in tissue development and aging. Development. 2013;140:2535-47. Shyh-Chang N, Daley GQ, Cantley LC.组织发育和衰老过程中的干细胞代谢。发育。2013;140:2535-47.
Agarwal M, Sharma A, Kumar P, Kumar A, Bharadwaj A, Saini M, et al. Myosin heavy chain-embryonic regulates skeletal muscle differentiation during mammalian development. Development. 2020;147(7):dev184507.
Schuster-Gossler K, Cordes R, Gossler A. Premature myogenic differentiation and depletion of progenitor cells cause severe muscle hypotrophy in Delta1 mutants. Proc Natl Acad Sci U S A. 2007;104:537-42. Schuster-Gossler K, Cordes R, Gossler A. 过早的肌原分化和祖细胞耗竭导致 Delta1 突变体肌肉严重萎缩。Proc Natl Acad Sci U S A. 2007;104:537-42.
Ruijtenberg S, van den Heuvel S. G1/S inhibitors and the SWI/SNF complex control cell-cycle exit during muscle differentiation. Cell. 2015;162:300-13. Ruijtenberg S, van den Heuvel S. G1/S 抑制剂和 SWI/SNF 复合物在肌肉分化过程中控制细胞周期的退出。细胞》。2015;162:300-13.
Spletter ML, Barz C, Yeroslaviz A, Schonbauer C, Ferreira IR, Sarov M, et al. The RNA-binding protein Arrest (Bruno) regulates alternative splicing to enable myofibril maturation in Drosophila flight muscle. EMBO Rep. 2015;16:178-91. Spletter ML, Barz C, Yeroslaviz A, Schonbauer C, Ferreira IR, Sarov M, et al. RNA 结合蛋白 Arrest (Bruno) 调节替代剪接以实现果蝇飞行肌肉中肌原纤维的成熟。EMBO Rep. 2015;16:178-91.
Nord AS, Blow MJ, Attanasio C, Akiyama JA, Holt A, Hosseini R, et al. Rapid and pervasive changes in genome-wide enhancer usage during mammalian development. Cell. 2013;155:1521-31. Nord AS, Blow MJ, Attanasio C, Akiyama JA, Holt A, Hosseini R, et al. 哺乳动物发育过程中全基因组增强子使用的快速和普遍变化。细胞。2013;155:1521-31.
El-Mounayri O, Triplett JW, Yates CW, Herring BP. Regulation of smooth muscle-specific gene expression by homeodomain proteins, Hoxa10 and Hoxb8. J Biol Chem. 2005;280:25854-63. El-Mounayri O, Triplett JW, Yates CW, Herring BP.同源结构域蛋白 Hoxa10 和 Hoxb8 对平滑肌特异性基因表达的调控。J Biol Chem.2005;280:25854-63.
Graeber CP, Hunter DG, Engle EC. The genetic basis of incomitant strabismus: consolidation of the current knowledge of the genetic foundations of disease. Semin Ophthalmol. 2013;28:427-37 Graeber CP, Hunter DG, Engle EC.不对称斜视的遗传基础：疾病遗传基础的现有知识整合。Semin Ophthalmol.2013;28:427-37
Taglietti V, Maroli G, Cermenati S, Monteverde S, Ferrante A, Rossi G, et al. Nfix induces a switch in Sox6 transcriptional activity to regulate MyHC-I expression in fetal muscle. Cell Rep. 2016;17:2354-66. Taglietti V, Maroli G, Cermenati S, Monteverde S, Ferrante A, Rossi G, et al. Nfix诱导Sox6转录活性转换以调节胎儿肌肉中MyHC-I的表达。Cell Rep. 2016;17:2354-66.
Mathsyaraja H, Freie B, Cheng P, Babaeva E, Eisenman RN. Max deletion destabilizes MYC protein and abrogates E-Myc lymphomagenesis. Genes Dev. 2019;33(17-18):1252-64. Mathsyaraja H、Freie B、Cheng P、Babaeva E、Eisenman RN。Max缺失会破坏MYC蛋白的稳定性，并抑制E-Myc淋巴瘤的发生。Genes Dev.2019;33(17-18):1252-64.
Carroll PA, Freie BW, Mathsyaraja H, Eisenman RN. The MYC transcription factor network: balancing metabolism, proliferation and oncogenesis, Front Med. 2018;12(4):412-25 Carroll PA, Freie BW, Mathsyaraja H, Eisenman RN.MYC转录因子网络：平衡新陈代谢、增殖和肿瘤发生》，Front Med.2018;12(4):412-25
Conacci-Sorrell M, Ngouenet C, Eisenman RN. Myc-Nick: a cytoplasmic cleavage product of Myc that promotes alpha-tubulin acetylation and cell differentiation. Cell. 2010;142:480-93. Conacci-Sorrell M, Ngouenet C, Eisenman RN.Myc-Nick：一种促进α-微管蛋白乙酰化和细胞分化的 Myc 细胞质裂解产物。细胞》。2010;142:480-93.
Gao MJ, Li X, Yang ZJ, Zhao S, Ling XX, Li JJ, et al. circHIPK3 regulates proliferation and differentiation of myoblast through the miR-7/TCF12 pathway. J Cell Physiol. 2021;236:6793-805.
McKinney-Freeman S, Hall T. Evaluation of Ranzoni et al.: integrative single-cell RNA-Seq and ATAC-Seq analysis of human developmental hematopoiesis. Cell Stem Cell. 2021;28:357-8. McKinney-Freeman S、Hall T. 对 Ranzoni 等人的评估：人类造血发育的单细胞 RNA-Seq 和 ATAC-Seq 整合分析。细胞干细胞》。2021;28:357-8.
Jia GS, Preussner J, Chen X, Guenther S, Yuan XJ, Yekelchyk M, et al. Single cell RNA-seq and ATAC-seq analysis of cardiac progenitor cell transition states and lineage settlement. Nat Commun. 2018;9(1):4877. Jia GS, Preussner J, Chen X, Guenther S, Yuan XJ, Yekelchyk M, et al. 单细胞RNA-seq和ATAC-seq分析心脏祖细胞的过渡状态和系谱定居。Nat Commun.2018;9(1):4877.
Miao Z, Balzer MS, Ma ZY, Liu HB, Wu JN, Shrestha R, et al. Single cell regulatory landscape of the mouse kidney highlights cellular differentiation programs and disease targets. Nat Commun. 2021;12(1):2277. Miao Z, Balzer MS, Ma ZY, Liu HB, Wu JN, Shrestha R, et al. 小鼠肾脏的单细胞调控图谱突显了细胞分化程序和疾病靶点。Nat Commun.2021;12(1):2277.
Schoenfelder, an Se, Furlan-Magaril, Mayra, Mifsud, Borbala, et al. The pluripotent regulatory circuitry connecting promoters to their longrange interacting elements. Genome Res. 2015; 25(4):582-597. Schoenfelder, an Se, Furlan-Magaril, Mayra, Mifsud, Borbala, et al.Genome Res. 2015; 25(4):582-597.
Javierre BM, Burren OS, Wilder SP, Kreuzhuber R, Hill SM, Sewitz S, et al. Lineage-specific genome architecture links enhancers and non-coding disease variants to target gene promoters. Cell. 2016;167:1369-1384. e1319. Javierre BM, Burren OS, Wilder SP, Kreuzhuber R, Hill SM, Sewitz S, et al. Lineage-specific genome architecture links enhancers and non-coding disease variants to target gene promoters.Cell.2016;167:1369-1384. e1319.
Freire-Pritchett P, Schoenfelder S, Várnai C, Wingett SW, Cairns J, Collier AJ, et al. Global reorganisation of cis-regulatory units upon lineage commitment of human embryonic stem cells. Elife. 2017;6:e21926. Freire-Pritchett P、Schoenfelder S、Várnai C、Wingett SW、Cairns J、Collier AJ 等人胚胎干细胞系承时顺式调控单元的全球重组。Elife.2017;6:e21926.
Blum R, Vethantham V, Bowman C, Rudnicki M, Dynlacht BD. Genomewide identification of enhancers in skeletal muscle: the role of MyoD1. Genes Dev. 2012;26:2763-79. Blum R, Vethantham V, Bowman C, Rudnicki M, Dynlacht BD.骨骼肌增强子的全基因组鉴定：MyoD1 的作用。Genes Dev.2012;26:2763-79.
Soleimani VD, Yin H, Jahani-AsI A, Ming H, Kockx CEM, van Ijcken WFJ, et al. Snail Regulates MyoD Binding-Site Occupancy to Direct Enhancer Switching and Differentiation-Specific Transcription in Myogenesis. Mol Cell. 2012;47:457-68 Soleimani VD, Yin H, Jahani-AsI A, Ming H, Kockx CEM, van Ijcken WFJ, et al. 蜗牛调控 MyoD 结合位点占位，引导增强子切换和分化特异性转录在肌肉发生中的作用。Mol Cell.2012;47:457-68
Mousavi K, Zare H, Dell'Orso S, Grontved L, Gutierrez-Cruz G, Derfoul A, et al. eRNAs promote transcription by establishing chromatin accessibility at defined genomic Loci. Mol Cell. 2013;51:606-17. Mousavi K, Zare H, Dell'Orso S, Grontved L, Gutierrez-Cruz G, Derfoul A, et al. eRNAs 通过在确定的基因组位点建立染色质可及性来促进转录。Mol Cell.2013;51:606-17.
Yoshida T, Delafontaine P. An intronic enhancer element regulates angiotensin II type 2 receptor expression during satellite cell differentiation, and its activity is suppressed in congestive heart failure. J Biol Chem. 2016;291:25578-90 Yoshida T, Delafontaine P. 在卫星细胞分化过程中，一个内含增强子元件调节血管紧张素 II 2 型受体的表达，其活性在充血性心力衰竭中受到抑制。J Biol Chem.2016;291:25578-90
Smemo S, Tena JJ, Kim KH, Gamazon ER, Sakabe NJ, Gomez-Marin C, et al. Obesity-associated variants within FTO form long-range functional connections with IRX3. Nature. 2014;507:371-5. Smemo S, Tena JJ, Kim KH, Gamazon ER, Sakabe NJ, Gomez-Marin C, et al. FTO 中的肥胖相关变体与 IRX3 形成长程功能连接。自然》杂志。2014;507:371-5.
Pagel JI, Deindl E. Early growth response 1-A transcription factor in the crossfire of signal transduction cascades. Indian J Biochem Biophys. 2011;48:226-35. 早期生长应答 1--信号转导级联交火中的转录因子。Indian J Biochem Biophys.
Zhang WW, Tong HL, Zhang ZH, Shao SL, Liu D, Li SF, et al. Transcription factor EGR1 promotes differentiation of bovine skeletal muscle satellite cells by regulating MyoG gene expression. J Cell Physiol. 2018;233:350-62. Zhang WW, Tong HL, Zhang ZH, Shao SL, Liu D, Li SF, et al. 转录因子EGR1通过调控MyoG基因表达促进牛骨骼肌卫星细胞分化J Cell Physiol.
Vega FM, Ridley AJ. The RhoB small GTPase in physiology and disease. Small GTPases. 2018;9:384-93 Vega FM, Ridley AJ.生理和疾病中的 RhoB 小 GTP 酶。Small GTPases.2018;9:384-93
Yang YL, Yan JY, Fan XH, Chen JX, Wang ZS, Liu XQ, et al. The genome variation and developmental transcriptome maps reveal genetic differentiation of skeletal muscle in pigs. PLoS Genet. 2021:17(11):e1009910. Yang YL, Yan JY, Fan XH, Chen JX, Wang ZS, Liu XQ, et al. 基因组变异和发育转录组图谱揭示了猪骨骼肌的遗传分化。PLoS Genet.2021:17(11):e1009910.
McGinnis CS, Murrow LM, Gartner ZJ. DoubletFinder: doublet detection in single-cell RNA sequencing data using artificial nearest neighbors. Cell Syst. 2019;8:329-37. McGinnis CS, Murrow LM, Gartner ZJ.DoubletFinder：使用人工近邻法检测单细胞 RNA 测序数据中的双点。Cell Syst.2019;8:329-37.
Butler A, Hoffman P, Smibert P, Papalexi E, Satija R. Integrating single-cel transcriptomic data across different conditions, technologies, and species. Nat Biotechnol. 2018;36:411-20 Butler A、Hoffman P、Smibert P、Papalexi E、Satija R. 跨不同条件、技术和物种整合单细胞转录组数据。Nat Biotechnol.2018;36:411-20
Stuart T, Butler A, Hoffman P, Hafemeister C, Papalexi E, Mauck WM, et al. Comprehensive integration of single-cell data. Cell. 2019;177:1888-902. Stuart T, Butler A, Hoffman P, Hafemeister C, Papalexi E, Mauck WM, et al. 单细胞数据的全面整合。细胞。2019;177:1888-902.
Trapnell C, Cacchiarelli D, Grimsby J, Pokharel P, Li S, Morse M, et al. The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells. Nat Biotechnol. 2014;32:381-6. Trapnell C, Cacchiarelli D, Grimsby J, Pokharel P, Li S, Morse M, et al. 单细胞伪时序排列揭示了细胞命运决定的动态和调节因素。Nat Biotechnol.2014;32:381-6.
Cai S, Zhu Q, Guo C, Yuan R, Zhang X, Nie Y, et al. MLL1 promotes myo genesis by epigenetically regulating Myf5. Cell Prolif. 2020;53:e12744. Cai S, Zhu Q, Guo C, Yuan R, Zhang X, Nie Y, et al. MLL1通过表观遗传调控Myf5促进肌生成 Cell Prolif.
Somites and myotomes at embryonic stages E18, E21, and E28 from Tibetan pigs (ZZ) and Duroc x Tibetan pigs (DZ). 2022. NCBI Sequence Read Archive https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc= GSE206914 藏猪（ZZ）和杜洛克×藏猪（DZ）胚胎 E18、E21 和 E28 阶段的体节和肌节。2022.NCBI 序列读取档案 https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc= GSE206914
Ralston E, Gutierrez-Cruz G, Kenea A, Brooks SR. Transcriptomic analysis of mouse muscles reveals a signature of early human Duchenne muscular dystrophy. bioRxiv. 2021:07.16.452553. https://doi.org/10.1101/ 2021.07.16.452553. NCBI Sequence Read Archive https://www.ncbi.nIm. nih.gov/geo/query/acc.cgi?acc=GSE162455. Ralston E, Gutierrez-Cruz G, Kenea A, Brooks SR. 小鼠肌肉转录组分析揭示了早期人类杜氏肌营养不良症的特征。bioRxiv. 2021:07.16.452553. https://doi.org/10.1101/ 2021.07.16.452553.NCBI Sequence Read Archive https://www.ncbi.nIm. nih.gov/geo/query/acc.cgi?acc=GSE162455.
Publisher's Note 出版商说明
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. 《施普林格-自然》对出版地图和机构隶属关系中的管辖权主张保持中立。
Ready to submit your research? Choose BMC and benefit from: 准备好提交您的研究成果了吗？选择 BMC，您将受益于
fast, convenient online submission 快速、便捷的在线提交
thorough peer review by experienced researchers in your field 由您所在领域经验丰富的研究人员进行全面的同行评审
rapid publication on acceptance 接受后迅速出版
support for research data, including large and complex data types 支持研究数据，包括大型复杂数据类型
gold Open Access which fosters wider collaboration and increased citations 金色的开放存取，促进更广泛的合作和更多的引用
maximum visibility for your research: over 100M website views per year 最大限度地提高研究成果的知名度：每年网站浏览量超过 1 亿次
At BMC, research is always in progress 在 BMC，研究始终在进行中
(c) The Author(s) 2023. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data. (c) 作者 2023。开放存取 本文采用知识共享署名 4.0 国际许可协议进行许可，该协议允许以任何媒介或格式使用、共享、改编、分发和复制本文，但需适当注明原作者和来源，提供知识共享许可协议的链接，并说明是否进行了修改。本文中的图片或其他第三方材料均包含在文章的知识共享许可协议中，除非在材料的署名栏中另有说明。如果材料未包含在文章的知识共享许可协议中，且您打算使用的材料不符合法律规定或超出许可使用范围，则您需要直接从版权所有者处获得许可。如需查看该许可的副本，请访问 http://creativecommons.org/licenses/by/4.0/。除非在数据的署名行另有说明，否则知识共享公共领域专用免责声明 ( http://creativeco mmons.org/publicdomain/zero/1.0/)适用于本文提供的数据。
(See figure on next page.) (见下页图）
Fig. 1 scRNA-seq identified major cell types in developing pig somites. A Experimental workflow schematic. Somite tissues [a mixture of embryos ( ) from the same sow] of ZZ and DZ at E16, E18, E21, and myotome tissues [a mixture of embryos ( ) from the same sow] at E28 were isolated. Tissue samples were dissociated into a single-cell solution and then single-cell transcriptomes were captured and analyzed using Genomics. The minimum scale of the ruler in the embryo photograph is 1 mm . t-Stochastic neighbor embedding (tSNE) plots showing the distribution of the main cell populations using the scRNA-seq. Using marker genes, cells were annotated as mesenchymal cells, fibroblasts, epithelial cells, neural stem cells, myogenic progenitors/myoblasts, osteogenic cells, neurons, neurogliocytes, endothelial cells, myocytes, chondrocytes, or muscle cells. Colors indicate cell types. Each dot represents one cell. C Heatmap showing the top 20 markers for each of the 12 cell populations. D Dot plot of the mean expression of canonical marker genes for 12 cell populations. Bar plot showing the percentage of different cell types within each sample. F Visualization of myogenic cell (including myogenic progenitors/myoblasts, myocytes, and muscle cells in Fig. 1B) sub-clusters via t-SNE by developmental stage (left) and sub-cluster number (right). G tSNE plots showing the cell identities of myogenic cell sub-clusters. Violin plots showing feature gene expression in each cell sub-cluster. Colors represent sub-clusters described in Fig. 1G 图 1 scRNA-seq 鉴定了发育中猪体节的主要细胞类型。A 实验工作流程示意图。分离 E16、E18 和 E21 期 ZZ 和 DZ 的体节组织[来自同一母猪的胚胎 ( ) 的混合物]和 E28 期的肌体组织[来自同一母猪的胚胎 ( ) 的混合物]。将组织样本解离成单细胞溶液，然后用 基因组学软件捕获和分析单细胞转录组。胚胎照片中标尺的最小刻度为 1 毫米。 t-Stochastic neighbor embedding (tSNE) 图显示了使用 scRNA-seq 的主要细胞群的分布情况。利用标记基因，细胞被注释为间充质细胞、成纤维细胞、上皮细胞、神经干细胞、成肌原/成肌细胞、成骨细胞、神经元、神经胶质细胞、内皮细胞、肌细胞、软骨细胞或肌肉细胞。颜色表示细胞类型。每个点代表一个细胞。C 热图显示 12 个细胞群中每个细胞群的前 20 个标记物。12 个细胞群的典型标记基因平均表达量的点阵图。 柱状图显示每个样本中不同细胞类型的百分比。F 根据发育阶段（左）和亚群数量（右），通过 t-SNE 显示肌原细胞（包括图 1B 中的肌原祖细胞/肌母细胞、肌细胞和肌肉细胞）亚群。G tSNE图显示成肌细胞亚簇的细胞特征。 显示各细胞亚簇特征基因表达的 Violin 图。颜色代表图 1G 中描述的亚簇
(See figure on next page.) (见下页图）
Fig. 2 Reconstruction of the myogenic differentiation trajectory in a pseudotime manner. A Pseudotime analysis of myogenic cells (including Pax3+ progenitors, myogenic progenitors, myoblasts, and myocytes in Fig. 1E) was performed by Monocle 2 and revealed seven different cell states (states 1 7). The distributions of cell states were presented along with pseudotime flows. Each dot is a cell. B Visualization of myogenic differentiation trajectory by cell origins (left) and developmental stages (right). C Visualization of myogenic differentiation trajectory by cell identity. D Violin plots showing feature gene expression in each cell cluster. MYMK, FNDC5, MEF2C, and TNNI1 are muscle development-related genes. MYL9 is a cardiomyocyte-specific marker. MYOD1 and MYOG are skeletal muscle cell-specific markers. Gene Ontology (GO) analysis of the differentially expressed genes with high levels in myocytes (state 2) and myocytes (state 7). F Visualization of myogenic differentiation trajectory by cell state, with cardiac cells distinguished from differentiating skeletal muscle cells. G Bar plot showing the percentage of cells within each sample assigned to the annotated myogenic cell types. H Immunofluorescence staining for Pax7 and MyoD on somite cross sections of ZZ and DZ at E21 and E28. Scale bar 图 2 以假时间方式重建肌原细胞分化轨迹。Monocle 2对成肌细胞（包括图1E中的Pax3+祖细胞、成肌细胞、肌母细胞和肌细胞）进行了伪时间分析，发现了七种不同的细胞状态（状态1-7）。细胞状态的分布与伪时间流一起呈现。每个点代表一个细胞。B 按细胞起源（左）和发育阶段（右）显示的成肌细胞分化轨迹。C 按细胞特征可视化肌原分化轨迹。D 显示各细胞群特征基因表达的小提琴图。MYMK、FNDC5、MEF2C 和 TNNI1 是肌肉发育相关基因。MYL9 是心肌细胞特异性标记。MYOD1和MYOG是骨骼肌细胞特异性标记。 肌细胞（状态 2）和肌细胞（状态 7）中高水平差异表达基因的基因本体（GO）分析。F 按细胞状态显示肌细胞分化轨迹，心肌细胞与分化中的骨骼肌细胞区分开来。G 柱状图显示每个样本中分配到注释的成肌细胞类型的细胞百分比。H E21 和 E28 时 ZZ 和 DZ 体节横切面上 Pax7 和 MyoD 的免疫荧光染色。比例尺 。
(See figure on next page.) (见下页图）
Fig. 3 Transcriptome dynamics of the myogenic differentiation. A Heatmap showing the expression changes of the 1700 top differentially expressed genes (DEGs) in a pseudotemporal order, with the DEGs, were cataloged into 5 five major clusters in characterized patterns (right). The GO analysis was performed for each gene cluster, and the representative enriched biological process (BP) terms are presented (left). B The expression dynamics of DEGs in gene clusters. Thick lines indicate the average gene expression patterns in each cluster (left). Gene signatures and expression dynamics of representative genes in each gene cluster (right) 图 3 成肌分化的转录组动态。热图显示了 1700 个顶级差异表达基因（DEGs）的表达变化，以假时序排列，DEGs 按特征模式被编入 5 个五大基因簇（右图）。对每个基因簇进行了 GO 分析，并给出了具有代表性的富集生物过程（BP）术语（左图）。B 基因簇中 DEGs 的表达动态。粗线表示每个基因簇的平均基因表达模式（左图）。每个基因簇中代表性基因的基因特征和表达动态（右图）
(See figure on next page.) (见下页图）
Fig. 5 Cell type-specific gene regulatory landscape of pig embryonic skeletal muscle. A Heatmap showing motif enrichment analysis on the cell type-specific open chromatin regions using Genomics (full results are shown in Additional file 7: Table S6). B TF expression z score heatmap that corresponding to the motif enrichment in each cell type. C Heatmap of cell type-specific regulons, as inferred by the SCENIC algorithm. Regulon activity was binarized to "on" (black) or "off" (white). D tSNE depiction of regulon activity ("on-blue","off-gray"), TF gene expression (red scale), and expression of predicted target genes (purple scale) of exemplary regulons for Pax3+ progenitors (ME/S1), myogenic progenitors (EZH2 and HDAC2), myoblasts (EGR1), and myocytes (MYOG). Examples of target gene expression of the TFs (PAX3, PCNA, SRSF7, RHOB, and SPG21) are shown in purple scale. Additional examples are given in Figure S6. The full list of regulons and their respective predicted target genes can be found in Additional file 11:Table S10 图 5 猪胚胎骨骼肌细胞特异性基因调控图谱。A 热图显示了使用 Genomics 对细胞类型特异性开放染色质区域进行的基因主题富集分析（全部结果见附加文件 7：表 S6）。B TF表达z得分热图，与每种细胞类型中的主题富集相对应。C SCENIC 算法推断出的细胞类型特异性调控子热图。调控子活性被二值化为 "开"（黑色）或 "关"（白色）。D tSNE 描述了 Pax3+ 祖细胞（ME/S1）、成肌细胞祖细胞（EZH2 和 HDAC2）、成肌细胞（EGR1）和成肌细胞（MYOG）的调控子活性（"开-蓝"、"关-灰"）、TF 基因表达（红色标度）和预测靶基因表达（紫色标度）。紫色标尺显示了 TFs（PAX3、PCNA、SRSF7、RHOB 和 SPG21）的靶基因表达实例。其他例子见图 S6。调控子及其各自预测的靶基因的完整列表见附加文件 11：表 S10。
(See figure on next page.) (见下页图）
Fig. 8 Functional analysis of EGR1 and RHOB in myogenesis. A To explore the connection network between TFs and targets in skeletal muscle development, 172 genes associated with muscle development (http://wiki.geneontology.org/index.php/Muscle biology) were extracted as target genes. Then, 215 high-confidence annotation TF-target pairs were selected from regulon activity network to construct the regulatory network. B Target genes counts of TFs. C The mRNA levels of EGR1, RHOB, and myogenic markers during the differentiation of pig primary myogenic cells (PPMCs) at several indicated time points. When the cells were cultured in growth medium (GM) at sub-confluent density, it was defined as day 0 (day 0); when the cells reached confluence, GM was changed to differentiation medium (DM). D The mRNA levels of EGR1, Myf5, and MyoD in proliferating C2C12 cells at 36 h after transfection with the pcDNA3.1-EGR1 vector. Immunofluorescence staining for MyHC in PPMCs after transfection with the pcDNA3.1-EGR1 vector and differentiation induction for 5 days. The fusion index (the percentage of nuclei in fused myotubes out of the total nuclei) was calculated. F The mRNA levels of myogenic differentiation markers in C2C12 cells after transfection with the pcDNA3.1-EGR1 vector and induction of differentiation for 3 days. G The mRNA levels of RHOB, Myf5, and MyoD1 in proliferating C2C12 cells at 36 h after transfection with the pcDNA3.1-RHOB vector. H Immunofluorescence co-staining for Pax7 and MyoD in C2C12 cells transfected with the pcDNA3.1-EGR1 vector and cultured in growth medium for 36 h . I Statistical analysis was performed to quantify the percentages of the three myogenic cell populations. J Immunofluorescence staining for MyHC in PPMCs after transfection with the pcDNA3.1-RHOB vector and induced differentiation for 5 days. The fusion index was calculated. , , 图 8 EGR1 和 RHOB 在肌肉发生过程中的功能分析A 为了探索骨骼肌发育过程中 TF 与靶标之间的连接网络，提取了 172 个与肌肉发育相关的基因（http://wiki.geneontology.org/index.php/Muscle biology）作为靶标基因。然后，从调控子活动网络中筛选出 215 个高置信度注释的 TF-靶标对，构建调控网络。B TF 的靶基因数量。C 猪原代成肌细胞（PPMCs）分化过程中 EGR1、RHOB 和成肌标志物的 mRNA 水平。当细胞以亚融合密度在生长培养基（GM）中培养时，定义为第0天（day 0）；当细胞达到 融合时，将生长培养基换成分化培养基（DM）。D 转染pcDNA3.1-EGR1载体36小时后，增殖的C2C12细胞中EGR1、Myf5和MyoD的mRNA水平。 转染pcDNA3.1-EGR1载体并诱导分化5天后，PPMCs中MyHC的免疫荧光染色。计算融合指数（融合肌管中的细胞核占总细胞核的百分比）。F 转染 pcDNA3.1-EGR1 载体并诱导分化 3 天后 C2C12 细胞中成肌分化标志物的 mRNA 水平。G 转染 pcDNA3.1-RHOB 载体 36 h 后，增殖的 C2C12 细胞中 RHOB、Myf5 和 MyoD1 的 mRNA 水平。H 用 pcDNA3.1-EGR1 载体转染并在生长培养基中培养 36 h 的 C2C12 细胞中 Pax7 和 MyoD 的免疫荧光共染。I 进行统计分析以量化三种成肌细胞群的百分比。 J 转染 pcDNA3.1-RHOB 载体并诱导分化 5 天后 PPMCs 中 MyHC 的免疫荧光染色。计算融合指数。 , ,
(See figure on next page.) (见下页图）
Fig. 9 Cell-cell communication analysis in the developing pig somite. Heatmaps showing the number of cell-cell interactions in the scRNA-seq dataset of myogenic cells (A) and all somite cells (B), as inferred by CellPhoneDB. Dark blue and dark red colors denote low and high numbers of cell-cell interactions, respectively. C CellPhoneDB-derived measures of cell-cell interaction scores and values among myogenic cells. Each row shows a ligand-receptor pair, and each column shows the two interacting cell types binned by cell type. Columns are sub-ordered by first interacting cell type into Pax3+ progenitors, myogenic progenitors, myoblasts, and myocytes. The color scale denotes the mean values for all the interacting partners, where the mean value refers to the total mean of the individual partner average expression values in the interacting cell-type pairs. D CellPhoneDB-derived measures of cell-cell interaction scores and values between myogenic cells and non-myogenic cells. Columns are sub-ordered by interacting myogenic cell type into Pax3+ progenitors, myogenic progenitors, myoblasts, and myocytes. E A diagram demonstrating the regulation of pig skeletal muscle ontogeny during embryonic development. In this model, during E16 to E 28, chromatin accessibility regulates the expression of classical and newly identified myogenic-related genes, which in turn promotes myogenic lineage fate determination and further myogenic differentiation 图 9 发育中猪体节的细胞-细胞通讯分析。 热图显示了由 CellPhoneDB 推断的肌原细胞（A）和所有体节细胞（B）的 scRNA-seq 数据集中细胞-细胞间相互作用的数量。深蓝色和深红色分别表示细胞间相互作用的低数量和高数量。C CellPhoneDB 导出的肌原细胞间细胞-细胞相互作用评分和 值。每行显示一对配体-受体，每列显示按细胞类型分类的两种相互作用的细胞类型。各列按第一个相互作用的细胞类型细分为 Pax3+祖细胞、肌原细胞、肌母细胞和肌细胞。色标表示所有相互作用伙伴的平均值，其中平均值是指相互作用细胞类型对中各个伙伴平均表达值的总平均值。D CellPhoneDB 导出的成肌细胞与非成肌细胞之间的细胞-细胞相互作用评分和 值。各列按相互作用的成肌细胞类型分为 Pax3+ 祖细胞、成肌细胞祖细胞、成肌细胞和成肌细胞。E 展示胚胎发育过程中猪骨骼肌本体发育调控的示意图。在该模型中，在 E16 至 E 28 期间，染色质可及性调控经典和新发现的肌生成相关基因的表达，进而促进肌生成系命运的确定和进一步的肌生成分化。
