神经轴突再生过程中微小RNA(miRNA)、转录因子和靶基因的网络功能分析

Network analysis of microRNAs, transcription factors, and target genes involved in axon regeneration # # 神经轴突再生过程中微小RNA(miRNA)、转录因子和靶基因的网络功能分析 1 first-author 1 2 3 1 1 corresp ‡Corresponding Author ‡Corresponding Author †E-mail:whp123456@sina.com †E-mail:whp123456@sina.com whp123456@sina.com Department of Basic Medicine, Hebei North University, Zhangjiakou 075029, China Department of Basic Medicine, Hebei North University, Zhangjiakou 075029, China Department of Pharmacy, Hebei North University, Zhangjiakou 075029, China Department of Pharmacy, Hebei North University, Zhangjiakou 075029, China Shenzhen RealOmics (Biotech) Co., Ltd., Shenzhen 518081, China Shenzhen RealOmics (Biotech) Co., Ltd., Shenzhen 518081, China Axon regeneration is crucial for recovery from neurological diseases. Numerous studies have identified several genes, microRNAs (miRNAs), and transcription factors (TFs) that influence axon regeneration. However, the regulatory networks involved have not been fully elucidated. In the present study, we analyzed a regulatory network of 51 miRNAs, 27 TFs, and 59 target genes, which is involved in axon regeneration. We identified 359 pairs of feed-forward loops (FFLs), seven important genes ( Nap1l1 , Arhgef12 , Sema6d , Akt3 , Trim2 , Rab11fip2 , and Rps6ka3 ), six important miRNAs (hsa-miR-204-5p, hsa-miR-124-3p, hsa-miR-26a-5p, hsa-miR-16-5p, hsa-miR-17-5p, and has-miR-15b-5p), and eight important TFs (Smada2, Fli1, Wt1, Sp6, Sp3, Smad4, Smad5, and Creb1), which appear to play an important role in axon regeneration. Functional enrichment analysis revealed that axon-associated genes are involved mainly in the regulation of cellular component organization, axonogenesis, and cell morphogenesis during neuronal differentiation. However, these findings need to be validated by further studies. 创新点 2 本研究比较详细地阐明了与神经轴突再生相关的微小RNA(miRNA)、转录因子和靶基因的相互作用关系,为神经系统疾病的恢复奠定基础. 方法 2 本研究从基因表达综合数据库中获得与轴突再生相关的转录因子和基因数据,利用文献报道及生物信息学数据库预测的方法筛选与轴突再生相关的基因靶向miRNA.利用生物信息学方法分析了三者之间的网络作用关系,预测了在相互作用网络中发挥重要作用的节点.最后对目标基因的基因本体(GO)功能进行了富集. 结论 2 通过分析,初步筛选了与神经轴突再生相关的51个miRNA、27个转录因子和59个靶标基因.进一步分析得到359对前馈环路,在此基础上推测了神经轴突再生过程中发挥重要作用的7个核心基因( Nap1l1 、 Arhgef12 、 Sema6d 、 Akt3 、 Trim2 、 Rab11fip2 和 Rps6ka3 ),6个miRNA(hsa-miR-204-5p、hsa-miR-124-3p、hsa-miR-26a-5p、hsa-miR-16-5p、hsa-miR-17-5p和hsa-miR-15b-5p)和8个转录因子(Smada2、Fli1、Wt1、Sp6、Sp3、 Smad4、Smad5和Creb1). Transcription factors miRNAs Target genes Axon Network analysis 转录因子 miRNA 靶基因 神经轴突 网络分析 Introduction 1 Axons are slender neuronal projections that transmit information from the cell body of a neuron around the body. Axons are also known as nerve fibers. Guidepost cells, typically other (sometimes immature) neurons, aid in the process of axon growth (Kunik, 2011 2011 ). Axon damage can result in the prevention of signal transmission, thereby contributing to the pathology of many neurological diseases (Debanne et al., 2011 2011 ). Previous studies have shown that an injured adult mammalian central nervous system, rich in inhibitory proteins and glycoproteins, has no capacity to initiate axonal regeneration spontaneously (Geoffroy and Zheng, 2014 2014 ). Therefore, promoting axon regeneration after injury is crucial. However, researchers have discovered that as long as the cell body of a neuron is not damaged, damaged axons can regenerate and recover function with the help of guidepost cells (Kunik, 2011 2011 ). The process of axon regrowth after injury may be fundamentally similar to the process of axon growth in embryonic stages. In the past decade, many studies have shown that axon regeneration can be induced by the knockout of negative regulators or overexpression of positive factors in neurons (Ferguson and Son, 2011 2011 ; Geoffroy and Zheng, 2014 2014 ; Baldwin and Giger, 2015 2015 ; McKerracher and Rosen, 2015 2015 ; He and Jin, 2016). For example, knockout of phosphatase and the tensin homolog (Pten) resulted in the promotion of axon regeneration in adult corticospinal neurons (Huang et al., 2017 2017 ). Overexpression of the inhibitor of DNA binding 2 (Id2) promoted axonal growth after injury (Yu et al., 2011 2011 ). Transcription factors (TFs) bind to target DNA sequences, either alone or as part of a complex, to increase or decrease gene transcription (Lee and Young, 2000 2000 ). TFs play important roles in developmental processes (Lobe, 1992 1992 ), signaling cascades (Osborne et al., 2001 2001 ), cell cycle control (Evan et al., 1994 1994 ), disease pathogenesis (Boch and Bonas, 2010 2010 ), and responses to environmental stimuli (Pullamsetti et al., 2016 2016 ). In recent years, studies have shown that TFs can stimulate axon regeneration via different molecular mechanisms (Venkatesh and Blackmore, 2017 2017 ). MicroRNAs (miRNAs) are small, non-coding RNA molecules that degrade target mRNAs or inhibit their translation by binding to the 3' untranslated region (3' UTR) (Tan et al., 2009 2009 ). miRNAs can regulate axon regeneration by controlling target gene expression (Jiang et al., 2015 2015 ). Therefore, understanding how TFs and miRNAs regulate the expression of axon-associated genes is critical for the development of therapies that promote axon regeneration. The expression of miRNAs can be mediated by TFs. TF and miRNA co-regulatory networks contain multiple feed-forward loops (FFLs) and feedback loops (Lin et al., 2015 2015 ). Interactions between these loops control gene expression to ensure the most suitable response to external stimuli (Wang et al., 2016 2016 ). FFLs include TFs that directly and indirectly regulate miRNA-target genes. In other words, miRNAs and TFs can co-regulate the same target genes. In this study, we investigated axon-associated TFs, miRNAs, and target genes using bioinformatics to construct a TF-miRNA co-regulatory network. In addition, we identified that FFLs were involved in axon regeneration. By analyzing these networks, we will better understand the cooperative TF-miRNA regulatory mechanisms and the role of FFLs in gene regulation, and be able to identify and describe regulatory mechanisms that could be used for promoting axon regeneration after being damaged. Methods 1 Identification of differentially expressed genes and hierarchical cluster analysis 2 GSE84975 microarray data were downloaded from the Gene Expression Omnibus database ( http://www.ncbi.nlm.nih.gov/geo http://www.ncbi.nlm.nih.gov/geo ). The GSE84975 dataset submitted on July 29, 2016 (Bigler et al., 2016 2016 ) included three neuron samples (GSM2254874, GSM2254875, GSM2254876) derived from human embryonic stem cells (hESC-Neuron samples) and three axon samples (GSM2254871, GSM2254872, GSM2254873) derived from hESC-Neurons (Axon hESC-Neuron samples). The dataset was based on GPL16686 (Affymetrix Human Gene 2.0 ST Array) platform information. The data were background-corrected and standardized using Expression Console™ software (Affymetrix) (Irizarry et al., 2003 2003 ). Genes that were differentially expressed between hESC-Neuron samples and Axon hESC-Neuron samples were identified by one-way between-subjects analysis of variance using Transcriptome Analysis Console 3.1 software (Affymetrix) (Torres et al., 2015 2015 ). Genes with thresholds of log 2 fold-changes > 2, P < 0.05, and false discovery rates (FDRs) < 0.05 were selected as differentially expressed genes (DEGs) (Table S1 S1 ). Non-annotated probes were filtered out manually. Hierarchical cluster analysis of DEGs was carried out using Transcriptome Analysis Console 3.1 software (Affymetrix). Identification of axon-related DEGs 2 We searched for axon-related genes in PubMed using “axon” as a search term. Genes related to axon growth, axon regeneration, and axon development were selected. Next, we identified those genes that overlapped between DEGs and axon-related genes. Finally, overlapping genes and genes in the same family as the overlapping genes were defined as Genes Cluster 1 (Table S2 S2 ). Axon-associated miRNAs and predictions of miRNA targets 2 We identified axon-related miRNAs by searching PubMed using “axon AND miRNA” as search terms. We selected miRNAs related to axon growth, axon regeneration, and axon development. We identified 51 mature miRNAs related to human axons (Table S3 S3 ). miRNA-target genes were predicted using the miRNA-target gene prediction databases miRanda and TargetScan. MiRanda software was developed to predict miRNA-target genes by bioinformatics methods. MiRanda screens the 3' UTR of mRNAs based on sequence matching, the thermal stability of the miRNA and mRNA double strand, and conservation of the target site. TargetScan was used to predict mammalian miRNA-target genes. It predicts conserved miRNA binding sites between different species based on a combination of thermodynamic models and RNA sequence analysis. The target genes identified by TargetScan and miRanda were analyzed and a total of 1852 targets were identified by both databases. These target genes were defined as Genes Cluster 2 (Table S4 S4 ). Prediction of axon-associated TFs 2 To select TFs that regulated axon-related DEGs, we downloaded human TFs from the transcriptional regulatory element database (TRED) ( http://rulai.cshl http://rulai.cshl . edu/TRED/TFlist.htm). Then, we screened for genes that overlapped between axon-related DEGs and human TFs. These were defined as TFs Cluster 1. Genes Cluster 3 included 59 genes that overlapped between Genes Cluster 1 and Genes Cluster 2. We predicted the interactions between the TFs in TFs Cluster 1 and 51 axon-related miRNAs or 59 axon-related genes. We searched for gene or miRNA sequences near the transcriptional start site area, and then used the TRANSFACT database (Biobase) to predict TF binding sites based on the Match™ algorithm. The search algorithm uses two score values, the matrix similarity score (MSS) and the core similarity score (CSS), to evaluate the result. In this study, 1 was the cutoff value for the MSS and CSS. FFLs and miRNA-TF-target gene networks 2 After identifying the regulatory relationships between miRNAs, TFs, and target genes, we used this information to construct FFLs and a miRNA-TF-target gene co-expression network. The network was constructed using Cytoscape v3.4.0 software (Shannon et al., 2003 2003 ). In this network, we selected some genes, miRNAs, and TFs that showed a high node degree and closeness centrality (Su et al., 2017b 2017b ). Node degree is a measure of the number of interactions with other proteins. Closeness centrality is a measure of key node centrality in the network. The higher the degree and closeness centrality values, the more important is the protein in the network. The nodes with degree values above the average network degree value and with thresholds of closeness centrality of > 0.5 were identified. Functional enrichment analysis 2 To explore the functions and pathways of axon-associated target genes, we performed gene ontology (GO) enrichment analysis using the BiNGO tool of Cytoscape v3.4.0 software (Shannon et al., 2003 2003 ). P < 0.01 was set as the cutoff value. Results 1 Screening for DEGs 2 A total of 4025 genes were identified as being differentially expressed between whole hESC-Neuron and Axon hESC-Neuron samples. Of these regulated genes, 1227 were up-regulated and 2798 down-regulated in Axon hESC-Neuron compared with whole hESC-Neuron samples (Table S1 S1 ). Hierarchical cluster analysis showed that the three Axon hESC-Neuron samples were distributed within the Axon sample cluster and the three hESC-Neuron samples within the hESC-Neuron sample cluster (Fig. 1 1 ). This demonstrated that the data could be used directly for further analysis. The horizontal axis shows sample names. GSM2254871, GSM2254872, and GSM2254873 are Axon hESC-Neuron samples. GSM2254874, GSM2254875, and GSM2254876 are hESC-Neuron samples. The right vertical axis shows clusters of DEGs. Up-regulated genes are shown in red and down-regulated genes in green (Note: for interpretation of the references to color in this figure legend, the reader is referred to the web version of this article) Interaction of axon-related genes and miRNAs 2 To explore the co-regulatory network of miRNAs and genes in axon regeneration, we selected 59 genes (Genes Cluster 3) that overlapped between Genes Cluster 1 (axon-associated genes) and Genes Cluster 2 (miRNA-target genes) (Table 1 1 ). Interactions between the 59 axon-related genes and miRNAs are shown in Table 2 2 . Notably, one gene may be regulated by multiple miRNAs, and each miRNA can target more than one gene. For example, hsa-miR-16-5p, hsa-miR-15b-5p, and hsa-miR-137 all target cell division cycle 42 ( Cdc42 ), and hsa-miR-15b-5p also targets serine/threonine kinase 3 ( Akt3 ), Cdc42 , nuclear factor of activated T-cells 3 ( Nfatc3 ), nicotinamide nucleotide adenylyltransferase 2 ( Nmnat2 ), paired box 2 ( Pax2 ), RAB11 family interacting protein 2 ( Rab11fip2 ), ribosomal protein S6 kinase A3 ( Rps6ka3 ), and semaphorin 6D ( Sema6d ). Gene symbol Chromosome Fold change P -value FDR Gsk3b chr1 −31.75 0.000 498 0.013 168 Trim44 chr16 −14.62 0.001 713 0.021 677 Akt3 chr1 −12.17 0.000 105 0.006 809 Cdc42bpa chr19 −12.09 0.000 389 0.011 808 Wdr33 chr9 −10.42 0.002 029 0.023 013 Rprd1a chrX −10.09 0.000 351 0.011 204 Cdc42 chr12 −8.98 0.000 686 0.015 027 Nap1l1 chr9 −8.91 0.000 061 0.005 459 Wdr47 chr1 −8.52 0.000 212 0.009 014 Klf12 chr1 −8.37 0.000 200 0.008 827 Rap1a chr20 −8.01 0.002 805 0.026 899 Trim33 chr16 −7.01 0.000 052 0.005 198 Aktip chr17 −7.00 0.001 652 0.021 337 Ppp2cb chr4 −6.16 0.000 081 0.006 185 Robo1 chr15 −5.43 0.000 086 0.006 299 Arhgef12 chr3 −5.39 0.000 422 0.012 217 Creb1 chr7 −5.26 0.001 453 0.020 248 Rps6ka3 chr20 −5.20 0.004 840 0.034 406 Alcam chr2 −4.97 0.000 362 0.011 369 Sacs chr21 −4.90 0.005 379 0.035 973 Wdr82 chr8 −4.72 0.004 450 0.033 020 Nova1 chr19 −4.52 0.000 392 0.011 870 Ppp6c chr15 −4.51 0.000 466 0.012 752 Pten chr3 −4.47 0.004 651 0.033 739 Abce1 chr1 −4.32 0.003 454 0.029 535 Wdr26 chr2 −4.29 0.000 160 0.008 128 Clasp2 chr20 −4.26 0.001 469 0.020 333 Sacm1l chr5 −4.24 0.001 135 0.018 153 Ppp2r5e chr2 −4.07 0.001 179 0.018 369 Ppp4r1 chr16 −3.96 0.000 083 0.006 236 Arfgef1 chr15 −3.84 0.001 885 0.022 476 Ppp1r9a chr16 −3.84 0.000 033 0.004 571 Sema6a chr1 −3.75 0.001 977 0.022 756 Cdc7 chr9 −3.58 0.010 118 0.049 437 Rpp14 chr19 −3.55 0.009 971 0.049 070 Ppt1 chr10 −3.28 0.000 129 0.007 393 Rab11fip2 chr14 −3.22 0.005 658 0.036 905 Ccdc25 chr1 −3.19 0.001 153 0.018 233 Ndel1 chr2 −3.10 0.007 999 0.043 986 Rps6kb1 chr5 −3.03 0.006 484 0.039 601 Sema6d chr2 −2.95 0.008 081 0.044 224 Cdc23 chr1 −2.83 0.008 324 0.044 824 Cdc42se1 chr4 −2.77 0.001 443 0.020 183 Ppp3r1 chr3 −2.77 0.009 574 0.048 062 Tp53inp1 chr16 −2.46 0.000 496 0.013 132 Nfatc3 chr7 −2.38 0.005 497 0.036 417 Trim2 chr19 −2.32 0.001 761 0.021 893 Nmnat2 chr8 −2.14 0.005 352 0.035 897 Trim9 chr1 −2.12 0.000 355 0.011 273 Ppp1r15b chr9 −2.10 0.001 839 0.022 286 Nrp1 chrX −2.06 0.001 276 0.019 120 Sox8 chr6 2.12 0.003 220 0.028 611 Pax2 chr19 2.14 0.003 335 0.029 089 Creb3l1 chr1 2.18 0.007 181 0.041 678 Cspg4 chr19 2.18 0.008 308 0.044 805 Ldlrap1 chr2 2.19 0.008 209 0.044 562 Wnt3a chr2 2.25 0.004 912 0.034 633 Rab11fip5 chr2 2.47 0.005 251 0.035 664 Ppp1r1b chr5 2.57 0.001 708 0.021 640 miRNA Gene symbol miRNA Gene symbol miRNA Gene symbol hsa-miR-29c-3p Abce1 hsa-miR-15b-5p Nfatc3 hsa-miR-181b-5p Rps6ka3 hsa-miR-214-5p Akt3 hsa-miR-15b-5p Nmnat2 hsa-miR-15b-5p Rps6ka3 hsa-miR-15b-5p Akt3 hsa-miR-96-5p Nova1 hsa-miR-214-3p Rps6ka3 hsa-miR-181d-5p Akt3 hsa-miR-204-5p Nova1 hsa-miR-204-3p Rps6ka3 hsa-miR-541-3p Akt3 hsa-miR-138-5p Nova1 hsa-miR-124-3p Rps6kb1 hsa-miR-17-5p Aktip hsa-miR-181d-5p Nova1 hsa-miR-181b-5p Sacm1l hsa-miR-18a-5p Alcam hsa-miR-221-3p Nova1 hsa-miR-181d-5p Sacm1l hsa-miR-382-5p Arfgef1 hsa-miR-338-3p Nova1 hsa-miR-17-5p Sacs hsa-miR-26a-5p Arfgef12 hsa-miR-338-5p Nrp1 hsa-miR-26a-5p Sacs hsa-miR-96-5p Arfgef12 hsa-miR-338-3p Nrp1 hsa-miR-214-3p Sema6a hsa-miR-214-3p Arfgef12 hsa-miR-15b-5p Pax2 hsa-miR-124-3p Sema6a hsa-miR-34a-3p Arfgef12 hsa-miR-221-3p Ppp1r15b hsa-miR-16-5p Sema6d hsa-miR-17-5p Ccdc25 hsa-miR-221-5p Ppp1r1b hsa-miR-15b-5p Sema6d hsa-miR-16-5p Cdc23 hsa-miR-18a-3p Ppp1r1b hsa-miR-124-3p Sema6d hsa-miR-16-5p Cdc42 hsa-miR-19a-5p Ppp1r9a hsa-miR-214-3p Sox8 hsa-miR-15b-5p Cdc42 hsa-miR-214-3p Ppp2cb hsa-miR-17-5p Tp53inp1 hsa-miR-137 Cdc42 hsa-miR-132-3p Ppp2cb hsa-miR-125b-5p Tp53inp1 hsa-miR-96-5p Cdc42bpa hsa-miR-221-3p Ppp2r5e hsa-miR-29c-3p Tp53inp1 hsa-miR-221-5p Cdc42bpa hsa-miR-17-5p Ppp3r1 hsa-miR-18a-5p Trim2 hsa-miR-541-5p Cdc42bpa hsa-miR-96-5p Ppp3r1 hsa-miR-181b-5p Trim2 hsa-miR-19a-3p Cdc42bpa hsa-miR-204-5p Ppp3r1 hsa-miR-181d-5p Trim2 hsa-miR-125b-5p Cdc42se1 hsa-miR-221-3p Ppp3r1 hsa-miR-338-3p Trim33 hsa-miR-29c-3p Cdc42se1 hsa-miR-342-3p Ppp3r1 hsa-miR-329-3p Trim33 hsa-miR-204-5p Cdc7 hsa-miR-7-5p Ppp4r1 hsa-miR-17-5p Trim44 hsa-miR-29c-3p Cdc7 hsa-miR-17-5p Ppp6c hsa-miR-34a-3p Trim44 hsa-miR-96-5p Clasp2 hsa-miR-125b-5p Ppt1 hsa-miR-96-5p Trim9 hsa-miR-17-5p Creb1 hsa-miR-26a-5p Pten hsa-miR-29c-3p Wdr26 hsa-miR-34a-5p Creb3l1 hsa-miR-19a-3p Pten hsa-miR-214-3p Wdr33 hsa-miR-29c-3p Cspg4 hsa-miR-16-5p Rab11fip2 hsa-miR-29c-3p Wdr33 hsa-miR-199a-5p Gsk3b hsa-miR-18a-5p Rab11fip2 hsa-miR-16-5p Wdr47 hsa-miR-7-5p Klf12 hsa-miR-181b-5p Rab11fip2 hsa-miR-7-5p Wdr47 hsa-miR-221-5p Klf12 hsa-miR-15b-5p Rab11fip2 hsa-miR-199a-3p Wdr47 hsa-miR-338-5p Klf12 hsa-miR-181d-5p Rab11fip2 hsa-miR-17-5p Wdr82 hsa-miR-431-5p Klf12 hsa-miR-7-5p Rab11fip5 hsa-miR-7-5p Wdr82 hsa-miR-329-3p Klf12 hsa-miR-19a-3p Rap1a hsa-miR-204-5p Wdr82 hsa-miR-124-3p Ldlrap1 hsa-miR-29c-3p Robo1 hsa-miR-181d-5p Wdr82 hsa-miR-124-5p Nap1l1 hsa-miR-204-3p Rpp14 hsa-miR-214-3p Wdr82 hsa-miR-96-3p Ndel1 hsa-miR-142-3p Rprd1A hsa-miR-541-3p Wdr82 hsa-miR-16-5p Nfatc3 hsa-miR-16-5p Rps6ka3 hsa-miR-16-5p Wnt3a hsa-miR-204-5p Nfatc3 hsa-miR-17-5p Rps6ka3 hsa-miR-15b-5p Wnt3a Prediction of axon-related TFs and their interactions with axon-related genes and miRNAs 2 TFs are proteins that can regulate miRNAs and genes. We downloaded TFs from the TRED database and identified 27 TFs in the DEG cluster; 17 were up-regulated and 10 were down-regulated (Table 3 3 ). Gene symbol Chromosome Entrez ID Fold change P -value FDR   Hoxc13 chr12 3229 2.73 0.005 219 0.035 533   Relb chr19 5971 2.40 0.001 108 0.018 001   Wt1 chr11 7490 2.38 0.000 424 0.012 237   Sp6 chr17 80320 2.34 0.000 336 0.011 081   Hoxd12 chr2 3238 2.32 0.005 568 0.036 634   Fli1 chr11 2313 2.16 0.001 253 0.019 005   Hoxd10 chr2 3236 2.15 0.003 911 0.031 160   Sp7 chr12 121340 2.15 0.006 933 0.040 949   Pax2 chr10 5076 2.14 0.003 335 0.029 089   Hoxa7 chr7 3204 2.13 0.007 633 0.042 904   Pou5f1 chr6 5459 2.12 0.002 400 0.024 910   Tlx1 chr10 3195 2.11 0.008 127 0.044 313   Pax5 chr9 5079 2.10 0.003 859 0.030 993   Hoxb7 chr17 3217 2.05 0.000 383 0.011 731   Tfap2e chr1 339488 2.04 0.005 196 0.035 462   Tfap2b chr6 7021 2.01 0.000 176 0.008 461   Hoxd8 chr2 3234 2.01 0.001 178 0.018 369   Sp3 chr2 6670 −2.05 0.005 965 0.037 936   Atf4 chr22 468 −3.30 0.005 055 0.035 038   Smad5 chr5 4090 −3.67 0.007 672 0.043 027   Atf6 chr1 22926 −4.28 0.001 025 0.017 499   Smad4 chr18 4089 −4.50 0.000 632 0.014 514   Hoxa2 chr7 3199 −4.97 0.000 511 0.013 316   Creb1 chr2 1385 −5.26 0.001 453 0.020 248   Smad2 chr18 4087 −5.31 0.001 404 0.019 969   Atf2 chr2 1386 −5.90 0.000 344 0.011 143   Sp4 chr7 6671 −6.79 0.000 361 0.011 366 We predicted interactions between TFs and axon-related genes and miRNAs. We identified 18 overlapping TF interactions with 59 axon-related genes and 18 overlapping TF interactions with 51 axon-related miRNAs. The interactions between axon-related genes and TFs are shown in Table S5 S5 . One TF can regulate several genes, and one gene may be regulated by multiple TFs. The interactions between axon-related miRNAs and TFs are presented in Table S6 S6 . One miRNA can target several TFs, and TFs may in turn regulate miRNAs. Also, one TF may be regulated by multiple miRNAs. miRNA-TF-target gene network of axonal regeneration 2 To explore miRNA and TF co-regulatory networks involved in axon regeneration, we identified 359 FFL pairs among miRNAs, TFs, and axon-related genes (Table S7 S7 ). These FFLs are illustrated in Fig. 2 2 . Next, we constructed a network that showed the regulatory relationships among miRNAs, TFs, and target genes for axonal regeneration using Cytoscape software (Fig. 3 3 ). This network showed that one gene may be regulated by several TFs and that a TF may indirectly affect the expression of other genes by several miRNAs. This network included 109 nodes and 632 edges. Eighteen (30.5%) of the 59 axon-related genes, 39 (76.5%) of the 51 axon-related miRNAs, and 14 TFs were recruited. This miRNA-TF-target gene regulatory network may uncover new regulatory mechanisms involved in axonal regeneration. Rectangles indicate miRNA, target genes or transcription factors; trigonal arrows indicate activation; blue indicates down-regulation. TF: transcription factor; Tar: target gene (Note: for interpretation of the references to color in this figure legend, the reader is referred to the web version of this article) Squares: miRNAs; circles: target genes; triangles: transcription factors. Red indicates up-regulation; blue indicates down-regulation; large nodes indicate bigger degrees (Note: for interpretation of the references to color in this figure legend, the reader is referred to the web version of this article) Nodes with highly connected portions may play important biological roles in a network. Finally, we identified seven genes that were down-regulated, six miRNAs, and eight TFs including two that were up-regulated and six that were down-regulated (Table 4 4 ). Node Closeness centrality Degree Expression Average degree   Rps6ka3 0.524 271 84 12 Down 21 11.60   Rab11fip2 0.521 739 13 12 Down   Trim2 0.519 230 77 12 Down   Akt3 0.526 829 27 12 Down   Sema6d 0.521 739 13 12 Down   Nap1l1 0.529 411 76 12 Down   Arhgef12 0.521 739 13 13 Down  hsa-miR-204-5p 0.521 739 13 12  hsa-miR-124-3p 0.529 411 76 13  hsa-miR-26a-5p 0.524 271 84 13  hsa-miR-16-5p 0.529 411 76 15  hsa-miR-17-5p 0.516 746 41 15  hsa-miR-15b-5p 0.534 653 47 17  Creb1 0.548 223 35 21 Down  Smad5 0.500 000 00 27 Down  Smad4 0.606 741 57 54 Down  Sp3 0.635 294 12 56 Down  Sp6 0.635 294 12 57 Up  Wt1 0.646 706 59 57 Down  Fli1 0.720 000 00 67 Up  Smad2 0.843 750 00 92 Down Gene functional enrichment 2 GO enrichment analysis revealed a significant enrichment of axon-related genes in different biological processes, including the regulation of cellular component organization, axonogenesis, cell morphogenesis during neuronal differentiation, and morphogenesis of neuronal projections (Table S8 S8 ). Yellow nodes in Fig. 4 4 indicate significant biological processes. Yellow node: P < 0.01 (Note: for interpretation of the references to color in this figure legend, the reader is referred to the web version of this article) Discussion 1 Many studies have found that axon regeneration after injury can be induced by altering the internal environment. Comparing the regulation relationships of axonal transcriptomes has the potential to highlight axonal regeneration mechanisms. It is difficult to obtain nervous tissue or culture-isolated neurons, so direct identification of primary human neurons is severely restricted. Recent studies have shown that hESCs have the capacity to differentiate specific neuronal subtypes. So in this study, we used microarray data (GSE84975) related to neurons derived from hESCs and axons of hESCs differentiated into neurons from the Gene Expression Omnibus database (Bigler et al., 2016 2016 ). We sought to identify the axonal genes of hESCs that had differentiated into neurons and predict the regulatory network of miRNAs, TFs, and target genes in the axon. We investigated the regulation of genes associated with axon regeneration by TFs and miRNAs. We identified 359 FFL pairs and one regulatory network associated with axon regeneration using a bioinformatics database. In the co-regulatory network, we identified some nodes with high degree and closeness centrality, including seven genes ( Nap1l1 , Arhgef12 , Sema6d , Akt3 , Trim2 , Rab11fip2 , and Rps6ka3 ), six miRNAs (hsa-miR-204-5p, hsa-miR-124-3p, has-miR-26a-5p, hsa-miR-16-5p, hsa-miR-17-5p, and hsa-miR-15b-5p), and eight TFs (Smad2, Fli1, Wt1, Sp6, Sp3, Smad4, Smad5, and Creb1). The genes appear to be predominately down-regulated in the co-regulatory network, whereas selected TFs were up-regulated. These nodes may play important regulatory roles in axon regeneration and may further our understanding of gene regulatory mechanisms involved in axon regeneration. For example, AKT3 is the predominant AKT kinase isoform in the brain (Easton, 2005 2005 ). However, Miao et al. ( 2016 2016 ) have shown that axon regeneration is enhanced even when AKT3 is down-regulated, which suggests that high levels of AKT3 may be not required for axon regeneration. Our results were consistent with this finding. Regulatory FFL networks identified in our study suggest that AKT3 can be regulated by TFs such as POU class 5 homeobox 1 (POU5F1), Sma-and Mad-related (SMAD) family member 4 (SMAD4), SMAD family member 2 (SMAD2), Fli-1 proto-oncogene (FLI1), Wilms tumor 1 (WT1), Sp3 transcription factor (SP3), Sp6 transcription factor (SP6), and cyclic adenosine monophosphate (cAMP) responsive element binding protein 1 (CREB1), and several miRNAs including hsa-miR-15b-5p, hsa-miR-214-5p, hsa-miR-541-3p, and hsa-miR-181d-5p. The regulatory FFL networks also suggest that AKT3 acts with these TFs and miRNAs to form 17 pairs of FFLs in the related network (Table S7 S7 ). For example, SP6 regulates hsa-miR-15b-5p, which targets Akt3 . SMAD transcription modulators regulate multiple cellular processes (Dahle and Kuehn, 2016 2016 ). SMADs are divided into three classes: receptor-regulated Smads (SMAD1, SMAD2, SMAD3, SMAD5, and SMAD8/9), the common-mediator Smad (SMAD4), and inhibitory Smads (SMAD6 and SMAD7). Some members of this family have been proved to take part in the process of axon regeneration. The expression of Smad1 was induced by phosphatidylinositol 3 kinase (PI3K)-glycogen synthase kinase 3 (GSK3) signaling to prevent mammalian axon regeneration (Saijilafu et al., 2013 2013 ). SMAD4 is critical for dorsal neural patterning (Chesnutt et al., 2004 2004 ). Overexpression of SMAD6 blocked dorsal interneuron 1 axon outgrowth (Hazen et al., 2011 2011 ). Moreover, Smad6 can inhibit the activity of the receptor-regulated Smads to regulate axon regeneration (Hazen et al., 2011 2011 ). Smad2 knockdown by RNA interference (RNAi) stimulated axonal growth in neurons (Stegmuller et al., 2008 2008 ). Furthermore, SMAD2 co-operated with the ubiquitin ligase Cdh1-anaphase-promoting complex upstream of the transcriptional modulator SnoN to control axonal growth (Stegmuller et al., 2008 2008 ). In the present study, we found that SMAD2 (which had the highest degree in the co-regulatory network) was down-regulated during axon regeneration, indicating a regulatory role in axonal regeneration. We showed that SMAD2 regulated 53 genes and 34 miRNAs. However, only five miRNAs (hsa-miR-125b-5p, hsa-miR-132-3p, hsa-miR-18a-5p, hsa-miR-26a-5p, and hsa-miR-541-5p) targeted SMAD2 in the network. miR-15b is abundant in distal axons (Natera-Naranjo et al., 2010 2010 ). In the present study, we showed that hsa-miR-15b-5p targets nine genes and two TFs (HOXC13 and SMAD5). Together, these two TFs targeted hsa-miR-26a-5p and nucleosome assembly protein 1 like 1 (NAP1L1). CREB1 is a protein that binds the cAMP response element to stimulate transcription. The activation of CREB promotes the transcription of many genes associated with neuronal survival, cell differentiation and axonal growth (Freitas et al., 2013 2013 ). In our previous study, we showed that CREB1 plays pivotal roles in neuronal differentiation (Su et al., 2017a 2017a ). In our present study, the results showed that CREB1 plays pivotal roles in the negative regulation of axons. Functional enrichment analysis revealed that the biological processes are regulated by genes associated with axon regeneration. These processes and the target genes are listed in Table S8 S8 . Axon regeneration-associated genes are involved mainly in the regulation of cellular component organization, axonogenesis, and cell morphogenesis during neuronal differentiation. In conclusion, we have identified TFs, miRNAs, and target genes that are involved in axon regeneration. Furthermore, we have revealed that the regulatory relationships potentially regulate axon regeneration. However, these associations will need to be experimentally validated in future studies. Such studies will improve our understanding of the mechanisms of axon-associated miRNA-TF-target gene regulation. References KT Baldwin , RJ Giger . Insights into the physiological role of CNS regeneration inhibitors . Front Mol Neurosci, 8:23 , 2015 . DOI: 10.3389/fnmol.2015.00023 http://doi.org/10.3389/fnmol.2015.00023 . Insights into the physiological role of CNS regeneration inhibitors S Baronchelli , A La Spada , P Conforti , 等 . Investigating DNA methylation dynamics and safety of human embryonic stem cell differentiation towards striatal neurons . Stem Cells Dev , 2015 . 24 ( 20 ): 2366 - 2377 . 26132372 DOI: 10.1089/scd.2015.0057 http://doi.org/10.1089/scd.2015.0057 . Investigating DNA methylation dynamics and safety of human embryonic stem cell differentiation towards striatal neurons RL Bigler , JW Kamande , R Dumitru , 等 . Axonal mRNA in human embryonic stem cell derived neurons . BioRxiv, 066142 , 2016 . DOI: 10.1101/066142 http://doi.org/10.1101/066142 . Axonal mRNA in human embryonic stem cell derived neurons J Boch , U Bonas . Xanthomonas AvrBs3 family-type III effectors: discovery and function . Ann Rev Phytopathol , 2010 . 48 ( 1 ): 419 - 436 . DOI: 10.1146/annurev-phyto-080508-081936 http://doi.org/10.1146/annurev-phyto-080508-081936 . Xanthomonas AvrBs3 family-type III effectors: discovery and function C Chesnutt , L Niswander . Plasmid-based short-hairpin RNA interference in the chicken embryo . Genesis , 2004 . 39 ( 2 ): 73 - 78 . 15170691 DOI: 10.1002/gene.20028 http://doi.org/10.1002/gene.20028 . Plasmid-based short-hairpin RNA interference in the chicken embryo O Dahle , MR Kuehn . Inhibiting Smad2/3 signaling in pluripotent mouse embryonic stem cells enhances endoderm formation by increasing transcriptional priming of lineage-specifying target genes . Dev Dyn , 2016 . 245 ( 7 ): 807 - 815 . 27012147 DOI: 10.1002/dvdy.24407 http://doi.org/10.1002/dvdy.24407 . Inhibiting Smad2/3 signaling in pluripotent mouse embryonic stem cells enhances endoderm formation by increasing transcriptional priming of lineage-specifying target genes D Debanne , E Campanac , A Bialowas , 等 . Axon physiology . Physiol Rev , 2011 . 91 ( 2 ): 555 - 602 . 21527732 DOI: 10.1152/physrev.00048.2009 http://doi.org/10.1152/physrev.00048.2009 . Axon physiology DM Easton . Voltage-clamp predictions by gompertz kinetics model relating squid-axon Na + -gating and ionic currents . Int J Neurosci , 2005 . 115 ( 10 ): 1415 - 1441 . 16162448 DOI: 10.1080/00207450590956521 http://doi.org/10.1080/00207450590956521 . Voltage-clamp predictions by gompertz kinetics model relating squid-axon Na+-gating and ionic currents G Evan , E Harrington , A Fanidi , 等 . Integrated control of cell proliferation and cell death by the c- myc oncogene . Philos Trans R Soc Lond B Biol Sci , 1994 . 345 ( 1313 ): 269 - 275 . 7846125 DOI: 10.1098/rstb.1994.0105 http://doi.org/10.1098/rstb.1994.0105 . Integrated control of cell proliferation and cell death by the c-myc oncogene TA Ferguson , YJ Son . Extrinsic and intrinsic determinants of nerve regeneration . J Tissue Eng , 2011 . 2 ( 1 ): 1 - 12 . DOI: 10.1177/2041731411418392 http://doi.org/10.1177/2041731411418392 . Extrinsic and intrinsic determinants of nerve regeneration AE Freitas , DG Machado , J Budni , 等 . Fluoxetine modulates hippocampal cell signaling pathways implicated in neuroplasticity in olfactory bulbectomized mice . Behav Brain Res , 2013 . 237 176 - 184 . 23018126 DOI: 10.1016/j.bbr.2012.09.035 http://doi.org/10.1016/j.bbr.2012.09.035 . Fluoxetine modulates hippocampal cell signaling pathways implicated in neuroplasticity in olfactory bulbectomized mice CG Geoffroy , B Zheng . Myelin-associated inhibitors in axonal growth after CNS injury . Curr Opin Neurobiol , 2014 . 27 31 - 38 . 24608164 DOI: 10.1016/j.conb.2014.02.012 http://doi.org/10.1016/j.conb.2014.02.012 . Myelin-associated inhibitors in axonal growth after CNS injury VM Hazen , KD Phana , S Hudiburgha , 等 . Inhibitory Smads differentially regulate cell fate specification and axon dynamics in the dorsal spinal cord . Dev Biol , 2011 . 2 ( 356 ): 566 - 575 . . Inhibitory Smads differentially regulate cell fate specification and axon dynamics in the dorsal spinal cord W Hu , Y He , Y Xiong , 等 . Derivation, expansion, and motor neuron differentiation of human-induced pluripotent stem cells with non-integrating episomal vectors and a defined xenogeneic-free culture system . Mol Neurobiol , 2016 . 53 ( 3 ): 1589 - 1600 . 25663198 DOI: 10.1007/s12035-014-9084-z http://doi.org/10.1007/s12035-014-9084-z . Derivation, expansion, and motor neuron differentiation of human-induced pluripotent stem cells with non-integrating episomal vectors and a defined xenogeneic-free culture system ZR Huang , ZZ Hu , P Xie , 等 . Tyrosine-mutated AAV2-mediated shRNA silencing of PTEN promotes axon regeneration of adult optic nerve . PLoS ONE , 2017 . 12 ( 3 ): e0174096 28323869 DOI: 10.1371/journal.pone.0174096 http://doi.org/10.1371/journal.pone.0174096 . Tyrosine-mutated AAV2-mediated shRNA silencing of PTEN promotes axon regeneration of adult optic nerve RA Irizarry , B Hobbs , F Collin , 等 . Exploration, normalization, and summaries of high density oligonucleotide array probe level data . Biostatistics , 2003 . 4 ( 2 ): 249 - 264 . 12925520 DOI: 10.1093/biostatistics/4.2.249 http://doi.org/10.1093/biostatistics/4.2.249 . Exploration, normalization, and summaries of high density oligonucleotide array probe level data JJ Jiang , CM Liu , BY Zhang , 等 . MicroRNA-26a supports mammalian axon regeneration in vivo by suppressing GSK3β expression . Cell Death Dis, 6:e1865 , 2015 . DOI: 10.1038/cddis.2015.239 http://doi.org/10.1038/cddis.2015.239 . MicroRNA-26a supports mammalian axon regeneration in vivo by suppressing GSK3β expression D Kunik . Laser-based single-axon transection for high-content axon injury and regeneration studies . PLoS ONE , 2011 . 6 ( 11 ): e26832 22073205 DOI: 10.1371/journal.pone.0026832 http://doi.org/10.1371/journal.pone.0026832 . Laser-based single-axon transection for high-content axon injury and regeneration studies TI Lee , RA Young . Transcription of eukaryotic protein-coding genes . Ann Rev Genet , 2000 . 34 ( 1 ): 77 - 137 . 11092823 DOI: 10.1146/annurev.genet.34.1.77 http://doi.org/10.1146/annurev.genet.34.1.77 . Transcription of eukaryotic protein-coding genes Y Lin , Q Zhang , HM Zhang , 等 . Transcription factor and miRNA co-regulatory network reveals shared and specific regulators in the development of B cell and T cell . Sci Rep, 5:15215 , 2015 . DOI: 10.1038/srep15215 http://doi.org/10.1038/srep15215 . Transcription factor and miRNA co-regulatory network reveals shared and specific regulators in the development of B cell and T cell CG Lobe . Transcription factors and mammalian development . Curr Top Dev Biol , 1992 . 27 351 - 383 . 1424766 DOI: 10.1016/S0070-2153(08)60539-6 http://doi.org/10.1016/S0070-2153(08)60539-6 . Transcription factors and mammalian development L McKerracher , KM Rosen . MAG, myelin and overcoming growth inhibition in the CNS . Front Mol Neurosci, 8:51 , 2015 . DOI: 10.3389/fnmol.2015.00051 http://doi.org/10.3389/fnmol.2015.00051 . MAG, myelin and overcoming growth inhibition in the CNS L Miao , L Yang , H Huang , 等 . mTORC1 is necessary but mTORC2 and GSK3β are inhibitory for AKT3-induced axon regeneration in the central nervous system . eLife, 5:e14908 , 2016 . DOI: 10.7554/eLife.14908 http://doi.org/10.7554/eLife.14908 . mTORC1 is necessary but mTORC2 and GSK3β are inhibitory for AKT3-induced axon regeneration in the central nervous system O Natera-Naranjo , A Aschrafi , AE Gioio , 等 . Identification and quantitative analyses of microRNAs located inhe distal axons of sympathetic neurons . RNA , 2010 . 16 ( 8 ): 1516 - 1529 . 20584895 DOI: 10.1261/rna.1833310 http://doi.org/10.1261/rna.1833310 . Identification and quantitative analyses of microRNAs located inhe distal axons of sympathetic neurons CK Osborne , R Schiff , SA Fuqua , 等 . Estrogen receptor: current understanding of its activation and modulation . Clin Cancer Res , 2001 . 7 ( 12 Suppl ): 4338s - 4342s . 11916222 . Estrogen receptor: current understanding of its activation and modulation SS Pullamsetti , F Perros , P Chelladurai , 等 . Transcription factors, transcriptional coregulators, and epigenetic modulation in the control of pulmonary vascular cell phenotype: therapeutic implications for pulmonary hypertension . Pulm Circ , 2016 . 6 ( 4 ): 448 - 464 . 28090287 DOI: 10.1086/688908 http://doi.org/10.1086/688908 . Transcription factors, transcriptional coregulators, and epigenetic modulation in the control of pulmonary vascular cell phenotype: therapeutic implications for pulmonary hypertension Saijilafu , EM Hur , ZX Jiao , 等 . PI3K-GSK3 signaling regulates mammalian axon regeneration by inducing the expression of Smad1 . Nat Commun, 4:2690 , 2013 . DOI: 10.1038/ncomms3690 http://doi.org/10.1038/ncomms3690 . PI3K-GSK3 signaling regulates mammalian axon regeneration by inducing the expression of Smad1 P Shannon , A Markiel , O Ozier , 等 . Cytoscape: a software environment for integrated models of biomolecular interaction networks . Genome Res , 2003 . 13 ( 11 ): 2498 - 2504 . 14597658 DOI: 10.1101/gr.1239303 http://doi.org/10.1101/gr.1239303 . Cytoscape: a software environment for integrated models of biomolecular interaction networks J Stegmuller , MA Huynh , Z Yuan , 等 . TGFβ-Smad2 signaling regulates the Cdh1-APC/SnoN pathway of axonal morphogenesis . J Neurosci , 2008 . 28 ( 8 ): 1961 - 1969 . 18287512 DOI: 10.1523/JNEUROSCI.3061-07.2008 http://doi.org/10.1523/JNEUROSCI.3061-07.2008 . TGFβ-Smad2 signaling regulates the Cdh1-APC/SnoN pathway of axonal morphogenesis LN Su , XQ Song , HP Wei , 等 . Identification of neuron-related genes for cell therapy of neurological disorders by network analysis . J Zhejiang Univ-Sci B (Biomed & Biotechnol) , 2017 . 18 ( 2 ): 172 - 182 . DOI: 10.1631/jzus.B1600109 http://doi.org/10.1631/jzus.B1600109 . Identification of neuron-related genes for cell therapy of neurological disorders by network analysis LN Su , YB Wang , CG Wang , 等 . Network analysis identifies common genes associated with obesity in six obesity-related diseases . J Zhejiang Univ-Sci B (Biomed & Biotechnol) , 2017 . 18 ( 8 ): 727 - 732 . DOI: 10.1631/jzus.B1600454 http://doi.org/10.1631/jzus.B1600454 . Network analysis identifies common genes associated with obesity in six obesity-related diseases Y Tan , B Zhang , T Wu , 等 . Transcriptional inhibiton of Hoxd4 expression by miRNA-10a in human breast cancer cells . BMC Mol Biol, 10:12 , 2009 . DOI: 10.1186/1471-2199-10-12 http://doi.org/10.1186/1471-2199-10-12 . Transcriptional inhibiton of Hoxd4 expression by miRNA-10a in human breast cancer cells L Torres , U Juarez , L Garcia , 等 . External ear microRNA expression profiles during mouse development . Int J Dev Biol , 2015 . 59 ( 10-12 ): 497 - 503 . 26864490 DOI: 10.1387/ijdb.150124sf http://doi.org/10.1387/ijdb.150124sf . External ear microRNA expression profiles during mouse development I Venkatesh , MG Blackmore . Selecting optimal combinations of transcription factors to promote axon regeneration: why mechanisms matter . Neurosci Lett , 2017 . 652 ( 23 ): 64 - 73 . 28025113 DOI: 10.1016/j.neulet.2016.12.032 http://doi.org/10.1016/j.neulet.2016.12.032 . Selecting optimal combinations of transcription factors to promote axon regeneration: why mechanisms matter H Wang , Z Xu , M Ma , 等 . Network analysis of microRNAs, transcription factors, target genes and host genes in nasopharyngeal carcinoma . Oncol Lett , 2016 . 11 ( 6 ): 3821 - 3828 . 27313701 DOI: 10.3892/ol.2016.4476 http://doi.org/10.3892/ol.2016.4476 . Network analysis of microRNAs, transcription factors, target genes and host genes in nasopharyngeal carcinoma P Yu , YP Zhang , LB Shields , 等 . Inhibitor of DNA binding 2 promotes sensory axonal growth after SCI . Exp Neurol , 2011 . 231 ( 1 ): 38 - 44 . 21679705 DOI: 10.1016/j.expneurol.2011.05.013 http://doi.org/10.1016/j.expneurol.2011.05.013 . Inhibitor of DNA binding 2 promotes sensory axonal growth after SCI Compliance with ethics guidelines Li-ning SU, Xiao-qing SONG, Zhan-xia XUE, Chen-qing ZHENG, Hai-feng YIN, and Hui-ping WEI declare that they have no conflict of interest. This article does not contain any studies with human or animal subjects performed by any of the authors. 10.1631/jzus.B1700179 Project supported by the Key Project of Hebei North University (No. 120177) and the Science and Technology Bureau Research Development Plan of Zhangjiakou City in Hebei (No. 0911021D-4), China 2017-04-05 2017-08-16 2022-01-05 Copyright © Zhejiang University and Springer-Verlag GmbH Germany, part of Springer Nature 2018 2018 Journal of Zhejiang University-SCIENCE B (Biomedicine and Biotechnology) 04 19