Abstract

   Myasthenia gravis (MG) is an autoimmune disease and the most common
   type of neuromuscular disease. Genes and miRNAs associated with MG have
   been widely studied; however, the molecular mechanisms of transcription
   factors (TFs) and the relationship among them remain unclear. A
   TF–miRNA–gene network (TMGN) of MG was constructed by extracting six
   regulatory pairs (TF–miRNA, miRNA–gene, TF–gene, miRNA–TF, gene–gene
   and miRNA–miRNA). Then, 3/4/5-node regulatory motifs were detected in
   the TMGN. Then, the motifs with the highest Z-score, occurring as
   3/4/5-node composite feed-forward loops (FFLs), were selected as
   statistically significant motifs. By merging these motifs together, we
   constructed a 3/4/5-node composite FFL motif-specific subnetwork
   (CFMSN). Then, pathway and GO enrichment analyses were performed to
   further elucidate the mechanism of MG. In addition, the genes, TFs and
   miRNAs in the CFMSN were also utilized to identify potential drugs.
   Five related genes, 3 TFs and 13 miRNAs, were extracted from the CFMSN.
   As the most important TF in the CFMSN, MYC was inferred to play a
   critical role in MG. Pathway enrichment analysis showed that the genes
   and miRNAs in the CFMSN were mainly enriched in pathways related to
   cancer and infections. Furthermore, 21 drugs were identified through
   the CFMSN, of which estradiol, estramustine, raloxifene and tamoxifen
   have the potential to be novel drugs to treat MG. The present study
   provides MG-related TFs by constructing the CFMSN for further
   experimental studies and provides a novel perspective for new
   biomarkers and potential drugs for MG.

   Subject terms: Data mining, Neurological disorders

Introduction

   Myasthenia gravis (MG) is an autoimmune disease in which antibodies
   destroy acetylcholine receptors on the postsynaptic membrane of the
   neuromuscular junction^[52]1. MG and its various subtypes are the major
   causes that affect the neuromuscular junction^[53]1. With the
   development of high-throughput data, research on MG from
   protein-encoding genes to noncoding RNAs, especially miRNAs, has been
   growing exponentially. However, the ultimate molecular mechanisms
   underlying the pathogenesis of MG remain to be further explored.

   Among many genetic factors, miRNAs and transcription factors (TFs) are
   two types of key gene regulators, and they both participate in many
   important cellular processes that share a common regulatory logic in
   the coregulation of target genes. MiRNAs mainly regulate gene
   expression at the posttranscriptional level, while TFs modulate gene
   transcription at the transcriptional level^[54]2. Substantial evidence
   has demonstrated that miRNAs participate in explaining the mechanism of
   MG. For example, the abnormal expression of miR-15a facilitates
   proinflammatory cytokine production by regulating the expression of
   C-X-C motif chemokine 10 (CXCL10), thereby contributing to the immune
   response in MG^[55]3. miR-139-5p and miR-452-5p were approved to
   negatively regulate regulator of G protein signalling 13 (RGS13)
   expression by assessing the miRNA and mRNA profiles of the MG
   thymus^[56]4. Moreover, it has been discovered that dysregulation of
   TFs, which leads to significant alterations in gene expression, is a
   novel phenomenon in MG. For instance, Yong et al. detected that altered
   expression of the transcription factors IRF4 and IRF8 is crucial for
   the counterbalancing mechanisms controlling the differentiation of
   plasma cells in patients with MG^[57]5. Li et al. verified that the
   overexpression of FRA1, also known as FOSL1, a FOS member of the
   activator protein 1 (AP-1) transcription factor, disrupted inflammatory
   cytokine secretion by medulla thymic epithelial cells (mTECs) in the MG
   thymus^[58]6. Nevertheless, as gene regulators, how miRNAs and TFs
   cooperate to regulate gene expression to cause the pathogenesis of MG
   has not yet been investigated.

   Recently, comprehensive analyses have been proposed to elucidate these
   gene regulators (miRNAs and TFs) as motifs, such as the feed-forward
   loop (FFL), which consists of two regulators, one of which regulates
   the other, and both collectively regulate a target gene^[59]7. It has
   been reported that FFLs can form recurrent network motifs to enhance
   the robustness of gene regulation in mammalian genomes^[60]8. The
   primary type was a three-node FFL consisting of a miRNA, a TF, and a
   common gene target^[61]7,[62]9. A three-node FFL could be extended to
   generate a four-node FFL, and a coexpressed gene was added on the basis
   of a three-node FFL as a common target of both miRNA and TF^[63]9. In
   the same way, a five-node FFL was created by introducing an additional
   miRNA–miRNA interaction to the existing four-node FFL^[64]10. According
   to the regulatory module between miRNA and TF, all FFLs are classified
   into the following three main types: miRNA FFL, TF FFL, and composite
   FFL^[65]9.

   A multitude of studies have been published to elucidate the underlying
   molecular mechanism by analysing FFL in many human diseases. For
   example, Engel et al. revealed neuronal activity-dependent P2X7R
   expression, which is induced by the transcription factor Sp1 and
   repressed in a calcium-dependent manner by microRNA-22^[66]11. Shi et
   al. proposed a novel method for identifying dysregulated miRNA–TF FFLs
   and excavating potential biomarkers for hypertrophic cardiomyopathy
   (HCM)^[67]12. A five-node FFL network was constructed to delineate
   miRNA, TF, and gene interactions in ischaemic stroke, and NFKB and STAT
   were identified as the chief regulators to explain the underlying
   mechanism of ischaemic stroke^[68]10. However, TF–miRNA–gene FFLs of MG
   have not been explored. Therefore, FFLs can be used to decipher the
   mechanisms of MG to provide new clues for understanding the
   pathogenesis and improving the treatment of MG.

   In this study, we designed a network-based analysis pipeline to
   decipher the potential mechanism and underlying drugs for MG, which is
   summarized in Supplementary Fig. [69]S1 online. First, a TF–miRNA–gene
   network of MG was constructed by extracting six regulatory pairs
   (TF–miRNA, miRNA–gene, TF–gene, miRNA–TF, gene–gene, and miRNA–miRNA)
   from several public databases. Then, 3-node, 4-node and 5-node
   regulatory motif types were detected in the network. Following the
   criteria we set to define 3-node, 4-node and 5-node motifs, the motif
   with the highest Z-score (3-node composite FFL, 4-node composite FFL
   and 5-node composite FFL) was selected as the statistically significant
   motif. By merging these motifs together, we constructed a 3/4/5-node
   composite FFL motif-specific subnetwork (CFMSN) and extracted related
   genes, TFs and miRNAs for further enrichment analysis. We found that
   the genes, TFs and miRNAs in the CFMSN were mainly enriched in
   cancerous and infectious pathways. In addition, the genes, TFs and
   miRNAs in the CFMSN were also applied to identify potential drugs to
   greatly improve the treatment of MG. Therefore, this study provides a
   novel perspective on the pathogenesis and treatment of MG.

Results

Construction of a global view of the TF–miRNA–gene network for MG

   By analysing the six regulatory relationships (miRNA–TF, miRNA–gene,
   TF–miRNA, TF–gene, gene–gene and miRNA–miRNA), we discovered several
   peculiarities. Among 276 miRNA–gene pairs, 82 MG risk genes were found
   to be validated targets for 73 MG risk miRNAs. Among 73 miRNAs,
   hsa-miR-125b-5p targeted the most genes, while VEGFA was the top target
   gene of miRNAs. Other top miRNAs were hsa-miR-145-5p and
   hsa-miR-29b-3p, while the top target genes were BCL2, IGF1R and MYC. Of
   39 miRNA–TF pairs, hsa-miR-29b-3p was found to target the greatest
   number of TFs, and MYC was targeted by the largest number of miRNAs. In
   addition, among 69 regulatory pairs of TF–miRNAs, the transcription
   factor MYC was found to regulate the highest number of miRNAs, while
   hsa-miR-17-5p was the top regulated miRNA in the list. In addition, the
   transcription factor MYC was identified to target the highest number of
   MG risk genes, and MYC and VEGFA were the most targeted genes among 51
   TF–gene pairs.

   We constructed a global view of the TF–miRNA–gene network (TMGN) for MG
   by merging the six regulatory relationships (miRNA–TF, miRNA–gene,
   TF–miRNA, TF–gene, gene–gene and miRNA–miRNA) identified above, and the
   network is visualized in Fig. [70]1A. The network included 248 nodes
   (16 TFs, 76 miRNAs and 156 genes) and 799 edges. Among these edges, 276
   belonged to miRNA–gene pairs, 39 belonged to miRNA–TF pairs, 47
   belonged to TF–gene pairs, 60 belonged to TF–miRNA pairs, 222 belonged
   to gene–gene pairs, and 194 belonged to miRNA–miRNA pairs.

Figure 1.

   [71]Figure 1
   [72]Open in a new tab

   The basic characteristics of TMGN with all six types of regulatory
   pairs (miRNA–gene, miRNA–TF, TF–miRNA, TF–gene, gene–gene and
   miRNA–miRNA). (A) A global of the network of TMGN. TFs, miRNAs and
   genes are colored green, orange and blue, respectively. (B–E) The basic
   features of the network include degrees, clustering coefficients,
   topological coefficients and neighborhood connectivity of the network.
   (F) A sub-network with MYC as the central node. (G) A sub-network with
   ESR1 as the central node. TMGN: TF–miRNA–gene network; TF:
   transcription factor.

   To examine the global view of the TMGN, we first analysed the
   topological features of this network, including degrees, clustering
   coefficients, topological coefficients and neighbourhood connectivity
   (Fig. [73]1B–E). We observed that a majority of nodes had a low degree,
   and only a few nodes had a high degree. Like many other biological
   networks, the degree distribution of this TMGN displayed a power law
   distribution f(x) = 78.338x^(− 1.177) with an R^2 of 0.873, indicating
   that TMGN followed a scale-free distribution and presented a
   small-world phenomenon^[74]13. The nodes with a larger degree are often
   network hubs and are considered to play important roles in maintaining
   the overall connectivity of the network^[75]14. In the TMGN, MYC and
   ESR1 had the top two largest numbers of transcription factors, so we
   generated subnetworks with MYC and ESR1 as the central node and their
   first neighbours (Fig. [76]1F,G); MYC was found to be connected to 5
   TF&genes, 16 genes and 39 miRNAs, while ESR1 was connected to 2
   TF&genes, 7 genes and 15 miRNAs. These results indicate that a global
   view of the TMGN could be a useful background for studies of MG.

Detection of 3/4/5-node motifs in the TMGN and discovery of feed-forward
loops

   By combining the miRNA–TF, miRNA–gene, TF–miRNA, TF–gene, gene–gene and
   miRNA–miRNA interactions, a TF–miRNA–gene network (TMGN) was
   constructed. Then, 3-node regulatory motifs were detected in the
   resulting network (TMGN) with four regulatory pairs. Here, we opted to
   investigate only motif types having a Z-score higher than 2 and p value
   less than 0.05. As a result, 15 different types of 3-node network
   motifs were identified using the FANMOD tool, and the results are
   visualized in Fig. [77]2A. Following the definition of 3-node motifs
   and compared with other motif types, the three-node composite FFL
   (miRNA–TF, miRNA–gene, TF–miRNA, and TF–gene) motif (surrounded by a
   red square in Fig. [78]2A) was selected as the statistically
   significant motif (Z-score: 2.4327, p value: 0.009).

Figure 2.

   [79]Figure 2
   [80]Open in a new tab

   3-node, 4-node and 5 node motifs. (A) 15 different types of 3-node
   regulatory motifs. (B) 29 different types of 4-node regulatory motifs.
   (C) 39 screened different types of 5-node regulatory motifs. The motifs
   are composed of miRNAs, TFs, and target genes and their Z-scores and p
   value s are presented. Orange triangles represent miRNAs, green rounds
   represent TFs and purple squares represent genes. An arrow ending in a
   circle represents the regulatory relationship, while, an arrow ending
   in a “T” represents the repression relationship. Six types of
   relationships involved in these motifs: miRNA–gene (miRNA represses
   gene expression); miRNA–TF (miRNA represses TF gene expression);
   TF–miRNA (TF regulates miRNA expression); TF–gene (TF regulates gene
   expression); gene–gene (gene and gene interaction); and miRNA–miRNA
   (miRNA and miRNA interaction). Motifs surrounded by red squares
   (3-node, 4-node and 5-node composite feed-forward loop) as significant
   motifs were merged to form 3-node, 4-node and 5-node regulatory
   sub-networks, respectively. TF: transcription factor; FFL: feed-forward
   loop.

   Similarly, 4-node regulatory motifs were also detected with five
   regulatory pairs (TF–miRNA, miRNA–gene, TF–gene, miRNA–TF, and
   gene–gene) using the FANMOD tool, and cut-off values of a Z-score
   higher than 2 and a p value less than 0.05 were set. As a result, a
   total of 29 different types of 4-node motifs were detected (shown in
   Fig. [81]2B). Following the definition of 4-node motifs, two four-node
   composite feed-forward loop (FFL) (TF–miRNA, miRNA–gene, TF–gene,
   miRNA–TF, and gene–gene) motifs (surrounded by red squares in
   Fig. [82]2B) met the criteria, and the motif with the highest Z-score
   was selected as the statistically significant motif (Z-score: 7.0735, p
   value: 0).

   Analogously, 5-node regulatory motifs were also detected with six
   regulatory pairs (TF–miRNA, miRNA–gene, TF–gene, miRNA–TF, gene–gene
   and miRNA–miRNA) using the FANMOD tool, and the same cut-off values of
   a Z-score higher than 2 and a p value less than 0.05 were set. As a
   consequence, 1775 different types of 5-node motifs were detected.
   Following the definition that we set for a 5-node motif, 39 different
   types of 5-node motifs were selected (visualized in Fig. [83]2C).
   Finally, a five-node composite feed-forward loop (FFL) (TF–miRNA,
   miRNA–gene, TF–gene, miRNA–TF, gene–gene and miRNA–miRNA) motif with
   the highest Z-score (Z-score: 50.147, p value: 0) was selected as the
   most statistically significant motif.

   Motifs with the highest Z-score (3-node, 4-node and 5-node composite
   feed-forward loop (FFL)) were selected as significant motifs and merged
   to form 3-node, 4-node and 5-node regulatory subnetworks, respectively.

Generation of a 3/4/5-node composite feed-forward loop (FFL) motif-specific
subnetwork

   The regulatory subnetwork, as well as the 3/4/5-node composite FFL
   motif-specific subnetwork (CFMSN), was visualized by Cytoscape 3.6.1.,
   which is presented in Fig. [84]3A. We obtained 3/4/5-node FFLs with 13
   miRNAs, 3 TFs, and 5 genes. The CFMSN comprised miR-20b-5p, miR-451a,
   miR-17-5p, miR-145-5p, miR-155-5p, miR-34a-5p, miR-20a-5p, miR-29b-5p,
   miR-221-5p, miR-29a-5p, let-7a-5p, let-7c-5p and let-7 g-5p as the
   principal miRNAs with three transcription factors (TFs), MYC, ESR1 and
   BCL6, and five genes, BCL2, VEGFA, KRAS, IL6 and MAPK1. Then, we
   extracted examples of the 3-node, 4-node and 5-node motifs (as well as
   3/4/5-node composite FFLs) from the CFMSN separately, as shown in
   Fig. [85]3B–D. These examples are consistent with the motifs we
   discovered in the above results. We also analysed the topological
   features of the CFMSN, including degrees, clustering coefficients,
   topological coefficients and neighbourhood connectivity
   (Fig. [86]3E–H). MYC has the largest degree in the CFMSN, indicating
   that as a transcription factor, MYC plays a central role in the
   mechanism of MG.

Figure 3.

   [87]Figure 3
   [88]Open in a new tab

   3/4/5-node composite FFL motif-specific sub-network (CFMSN). (A) A
   global view of 3/4/5-node composite FFL motif-specific sub-network
   (CFMSN). Orange triangle represents miRNA, green round represents TF
   and blue square represents gene. A line ending in a circle represents
   the regulatory relationship, while, a line ending in a “T” represents
   the repression relationship. (B–D) The examples of 3/4/5 node composite
   FFL. (E–H) The basic features of the network include degrees,
   clustering coefficients, topological coefficients and neighborhood
   connectivity of the network. TF: transcription factor; FFL:
   feed-forward loop.

Validation of the expression of miRNAs and regulatory pairs in CFMSN

   To verify the expression of the 13 miRNAs we identified through the
   CFMSN, differential expression analysis was performed using GEO2R with
   the dataset [89]GSE103812^[90]4. The result of the analysis is shown in
   Table [91]S1. We can conclude that all 13 miRNAs were expressed in MG
   thymus tissue. hsa-miR-145-5p was considered to be the differentially
   expressed miRNA. Since this dataset is a study of ectopic germinal
   centres in the thymus of MG, the expression of miRNAs differs to some
   extent and is for reference only. Moreover, the correlation of
   expressions between these relationship pairs in the CFMSN was also
   assessed to further evaluate the accuracy of the CFMSN. The results of
   the correlation analysis demonstrated that 41 out of 96 regulatory
   pairs (42.7%) in the CFMSN showed significant correlations, which
   indirectly showed that the extracted network (CFMSN) is meaningful.
   Detail information is shown in Table [92]S2.

Gene Ontology and pathway analyses of genes and miRNAs in the CFMSN

   With the aid of the DAVID database, we performed Gene Ontology analysis
   and KEGG pathway analysis using genes and targets of miRNAs in the
   CFMSN. Genes that underwent enrichment analysis included MAPK1, IL6,
   KRAS, BCL2 and VEGFA. As a result, 41 significant KEGG pathways
   (p < 0.05) and 40 significant GO terms (p < 0.05) were enriched, and
   the top 10 pathways and GO terms are shown on the left of
   Fig. [93]4A,B. Similarly, the target genes of 13 miRNAs in the CFMSN
   were also enriched, and 107 KEGG pathways (p < 0.05) and 1127 GO terms
   (p < 0.05) were significantly enriched. The top 10 pathways and GO
   terms are shown on the right of Fig. [94]4A,B. In comparison with the
   left and right of Fig. [95]4A,B, seven pathways (PI3K-Akt signalling
   pathway and pathways in cancer, hepatitis B, bladder cancer, colorectal
   cancer, pancreatic cancer and prostate cancer) and two gene functions
   (positive regulation of cell proliferation and cytoplasm) were
   collectively enriched in both gene and miRNA analyses. Among these
   coenriched pathways, five were associated with cancer, and one was
   associated with infectious diseases, which is in accordance with our
   previous findings^[96]15–[97]17.

Figure 4.

   [98]Figure 4
   [99]Open in a new tab

   Gene Ontology and KEGG pathway analysis using genes and targets of
   miRNAs in CFMSN. (A) The top 10 pathways enriched by genes (left) and
   target genes of miRNAs (right). Pathways colored in red are common ones
   in both left and right. (B) The top 10 GO terms enriched by genes
   (left) and target genes of miRNAs (right). GO terms colored in red are
   common ones in both left and right. (C) A part of PI3K-Akt signaling
   pathway. (D) A part of Hepatitis B pathway. Genes with a yellow
   background is the risk genes from CFMSN. GO: Gene Ontology; KEGG: Kyoto
   Encyclopedia of Genes and Genomes; CFMSN: 3/4/5-node composite FFL
   motif-specific sub-network.

   Next, three significant pathways (hsa04151:PI3K-Akt signalling pathway,
   hsa05161: hepatitis B pathway and hsa05200:pathways in cancer) were
   selected for further analysis. BCL2 alone was found to be enriched in
   the PI3K-Akt signalling pathway (part of the pathway is shown in
   Fig. [100]4C) and hepatitis B pathway (part of the pathway is shown in
   Fig. [101]4D), taking part in cell survival and antiapoptosis. BCL2,
   VEGFA and IL-6 were found to be enriched in pathways in cancer.
   Furthermore, the PI3K-Akt signalling pathway was found to be
   simultaneously involved in apoptotic regulation of the hepatitis B
   pathway and pathways in cancer. Our results indicate that the PI3K-Akt
   signalling pathway may be involved in regulation of apoptosis in the
   pathogenesis of MG.

Excavation of potential drugs for MG with the CFMSN

   Based on the construction of the CFMSN, drug analysis was performed.
   Thus, 3 TFs, 5 genes and 13 miRNAs identified from the CFMSN were
   applied to discover potential drugs. The three TFs (MYC, ESR1 and BCL6)
   screened from the CFMSN were applied to identify potential drugs
   through the DrugBank database; as a result, only ESR1 was found to be
   targeted by 38 approved drugs. We discarded drugs that have not been
   approved by the FDA, are not available on the market, can lead to
   carcinogenesis, are illegal, have only mixed products with other drugs,
   are vaginally administered, or lack APRD. After screening for the above
   criteria, 12 drugs targeting ESR1 remained (shown in Fig. [102]5A).

Figure 5.

   [103]Figure 5
   [104]Open in a new tab

   CFMSN and associated drugs. (A) A sub-network of TF in CFMSN and
   selected drugs. (B) A sub-network of genes in CFMSN and selected drugs.
   (C) A sub-network of miRNAs in CFMSN and selected drugs. (D) A global
   view of CFMSN and selected drugs through TFs, genes and miRNAs in
   CFMSN. (E) A sub-network of selected drugs targeted with screened TFs,
   genes and miRNAs. Orange triangle represents miRNA, green round
   represents TF, blue square represents gene and red “V” represents drug.
   (F) Part of pathway 07,226: Progesterone, androgen and estrogen
   receptor agonists/antagonists. Drugs colored in red are from 21 drugs
   we selected. CFMSN: 3/4/5-node composite FFL motif-specific
   sub-network; TF: transcription factor.

   We also identified five genes (BCL2, VEGFA, KRAS, IL6 and MAPK1) by
   constructing the CFMSN. By searching the DrugBank database, except
   KRAS, the other four out of these five genes were approved as drug
   targets. The drugs and their target genes are shown in Fig. [105]5B.
   These drugs are consistent with the criteria above.

   The CFMSN contained miR-20b-5p, miR-451a, miR-17-5p, miR-145-5p,
   miR-155-5p, miR-34a-5p, miR-20a-5p, miR-29b-5p, miR-221-5p, miR-29a-5p,
   let-7a-5p, let-7c-5p and let-7 g-5p as the principal miRNAs. To
   determine how these screened miRNAs can affect the therapy and the drug
   targets of MG, miRNA–drug pairs were obtained by merging miRNA–gene
   pairs and drug-gene pairs from DrugBank. Then, we discarded the drugs
   with fewer than 2 target genes. A cumulative hypergeometric
   distribution was performed to identify significant miRNA–drug pairs. We
   considered results with a p value less than 0.05 to be statistically
   significant. As a result, 4 miRNAs, 13 drugs and 27 miRNA–drug pairs
   were identified. Only 4 drugs (estradiol, estramustine, raloxifene and
   tamoxifen) met the criteria we set above; thus, 2 miRNAs and 8
   miRNA–drug pairs were selected (shown in Fig. [106]5C).

   A global view of the CFMSN and selected drugs through TFs, genes and
   miRNAs in the CFMSN is visualized in Fig. [107]5D, and a subnetwork of
   selected drugs targeted with screened TFs, genes and miRNAs is
   visualized in Fig. [108]5E. A total of 1 TF, 2 miRNAs, 4 genes and 21
   drugs were identified. Among the 21 drugs, estradiol, estramustine,
   raloxifene and tamoxifen had the largest degrees. Detailed information
   on these four drugs is summarized in Table [109]1.

Table 1.

   Detail information about Estradiol, Estramustine, Raloxifene and
   Tamoxifen.
   Drug Description References