Identification of differentially expressed genes
The microarray expression dataset GSE57647 was obtained from GEO and was analyzed in R using bioconductor packages. The dataset comprised with two pairs of each uninfected and JE infected for 6 hours, 24 hours and 48 hours samples of miRNA and mRNA in human microglial cells (CHME3) infected with Japanese encephalitis virus (JEV), P20778 strain. The principal components analysis revealed the relatedness between samples of each category as shown in the Fig. 1A and 1B. The dataset was normalized by using RMA normalization approach and DEGs were identified by using the threshold of logFC 1 and a Benjamini & Hochberg adjusted p-value cut- off of 0.05 as shown in the Fig. 1C and 1D. The Venn diagram plots clearly classify the distinct expression of genes in uninfected and JEV infected samples at 6hrs and 48hrs (Fig. 1C and D). Each cluster represents the sample with same feature within it and differentiates with the samples in other cluster. Next, DEGs were identified that were further used for network construction.
GGIN, enrichment analysis and identification of target genes information
With the help of GeneMANIA tool, the DEGs were found to have well-established role in the pathogenesis of JEV based on co-expression, co-localization, and pathways. The GGIN obtained from STRING consisted of 1509 nodes and 21993 edges, with average local clustering coefficient of 0.381 and average node degree of 29. The CytoHubba plugin of Cytoscape 3.7.1 software identified central elements of the GGIN. Hubs are the nodes with maximum interactions and occupy central position in the interaction networks. Maximum Clique Centrality measure identified the hub genes in upregulated dataset includes BCLAF1, VAMP3, BIRC5, MEST, EHD1, WNT5A, SYNCRIP, CDH13, SS18, HDAC8, GNG4, MGAT2, ATL3, TAF9B genes and downregulated dataset includes SPAG9, SLCO4A1, VEGFA, TM4SF1, NABP1, ASPH genes (Fig. 2A and 2B). Further these genes showed the least genetic alterations therefore selected as target genes for the further study. The enrichment analysis of target genes revealed that there is an association between biological process (biological regulation), cellular component (membrane), & molecular function (protein binding) as shown in Fig. 2C.
Identification of validated miRNAs and Transcriptions Factors (TF) for target genes
RegNetworkas and miRtarBase tools were used to identify the transcription factors (TF) and validated miRNAs for the regulatory elements of target genes (Table 1). The experimentally validated miRNAs and TF were identified the target genes i.e. VEGFA and WNT5A; for other target genes, no validated regulators were found. Our data reveled that EGR1, EPAS1, ETS1, FOS, HIF1A were the main transcription factors that regulate the expression of VEGFA (Table 1). Similarly, transcription factors GLI3, GLI2, GLI1 was found to regulate the expression of WNT5A gene (Table 1). miRtarBase tool identified several miRNAs for both VEGFA and WNT5A that are listed in Table 1 (few top ranked miRNAs).
Table 1
The table shows the list of transcription factors and miRNAs that are known to regulate the expression of VEGFA and WNT5A genes.
Gene | Regulatory Elements |
VEGFA | Transcription Factor | miRNA |
EGR1, EPAS1, ETS1, FOS, HIF1A | hsa-miR-126-3p |
hsa-miR-29b-3p |
hsa-miR-15a-5p |
mmu-miR-20b-5p |
hsa-miR-125a-5p |
hsa-miR-361-5p |
hsa-miR-20b-5p |
hsa-miR-106a-5p |
hsa-miR-106b-5p |
hsa-miR-205-5p |
WNT5A | GLI3, GLI2, GLI1 | hsa-miR-30a-5p |
hsa-let-7c-3p |
hsa-miR-190a-3p |
hsa-miR-487b-3p |
hsa-miR-4635 |
hsa-miR-5187-3p |
hsa-miR-4680-3p |
hsa-miR-4427 |
hsa-miR-33a-3p |
Seed Pairing and Construction of Gene Regulatory Network reveals miRNA-Mediated gene regulation
Next, we further characterized the interaction of gene and miRNAs observed in this study. Our data revealed that the seed pairing of VEGFA-hsa-miR-205-hsa-miR-20b forms a triplex with a binding free energy of -43.76kcal/mol (Table 2) to regulate the expression of VEGFA. Furthermore, WNT5A-hsa-miR-330-5p- hsa-miR-455-5p also forms a canonical triplex with a binding free energy of -36.86kcal/mol (Table 2). The table shows the top ten best gene-miRNAs pairing for both VEGFA and WNT5A. Our data show that the negative binding free energy of VEGFA with its miRNAs were more than that of WNT5A-miRNAs (Table 2) indicating that the interaction of VEGFA-miRNAs is more stable than WNT5A-miRNAs.
Table 2
The table shows the interaction of RefSeq (VEGFA and WNT5A), Gene-miRNA (Seed pairing) along with the seed distance and free energy of gene-miRNAs interactions.
Gene | RefSeq ID | miRNA1ID | miRNA2ID | Seed distance (nt) | Binding Free energy (kcal/mol) |
VEGFA | NM_001033756 | hsa-miR-205 | hsa-miR-20b | 27 | -43.76 |
AF323587 | hsa-miR-150 | hsa-miR-339-5p | 30 | -42.86 |
NM_001033756 | hsa-miR-17 | hsa-miR-205 | 27 | -42.06 |
NM_001033756 | hsa-miR-106a | hsa-miR-205 | 27 | -41.26 |
NM_001033756 | hsa-miR-20a | hsa-miR-205 | 27 | -41.16 |
NM_001033756 | hsa-miR-93 | hsa-miR-205 | 27 | -40.16 |
NM_001033756 | hsa-miR-205 | hsa-miR-106b | 27 | -39.46 |
NM_001033756 | hsa-miR-205 | hsa-miR-302e | 26 | -37.86 |
NM_001033756 | hsa-miR-199a-5p | hsa-miR-150 | 33 | -37.86 |
NM_001033756 | hsa-miR-410 | hsa-miR-590-3p | 25 | -19.26 |
WNT5A | NM_003392 | hsa-miR-330-5p | hsa-miR-455-5p | 18 | -36.86 |
NM_003392 | hsa-miR-129-5p | hsa-miR-181b | 25 | -35.06 |
NM_003392 | hsa-miR-129-5p | hsa-miR-181d | 25 | -33.26 |
NM_003392 | hsa-miR-136 | hsa-miR-429 | 23 | -32.56 |
NM_003392 | hsa-miR-200b | hsa-miR-539 | 32 | -31.16 |
NM_003392 | hsa-miR-326 | hsa-miR-455-5p | 18 | -30.56 |
NM_003392 | hsa-miR-136 | hsa-miR-200c | 23 | -29.36 |
NM_003392 | hsa-miR-200b | hsa-miR-136 | 23 | -28.56 |
NM_003392 | hsa-miR-205 | hsa-miR-142-3p | 15 | -24.36 |
NM_003392 | hsa-miR-129-5p | hsa-miR-200a | 14 | -23.66 |
Analytical modeling of gene regulation to understand the regulatory pathogenesis of JEV
Since, the seed paring of VEGFA with miRNAs were more stable; therefore, for subsequent study VEGFA was taken. In the mathematical modeling of the miRNA mediated transcriptional regulatory network which was associated with the pathogenesis of JEV, the involved factors are VEGFA gene, activated by the TFs: EGR1, EPAS1, ETS1, FOS, HIF1A and repressed by the miRNA hsa-miR-205 with the cooperation of hsa-miR-20b. Since the cooperative miRNA- hsa-miR-20b act as a catalyst in driving the pathogenesis, the core miRNA, hsa-miR-205 is preferred in the modeling and the ordinary differential equation of the mathematical model that initiates the regulatory pathogenesis in JEV is illustrated as:[d (VEGFA)/ dt = k(Synthesis of EGR1, EPAS1, ETS1, FOS, HIF1A) + VEGFA + k(Up-regulation of EGR1, EPAS1, ETS1, FOS, HIF1A) (hsa-miR-205) - k(Degradation of EGR1, EPAS1, ETS1, FOS, HIF1A) (EGR1, EPAS1, ETS1, FOS, HIF1A) where k represents the rate of synthesis/degradation of TFs. Similar analysis can be done for WNT5A gene-miRNA pairing.
mRNA-miRNA-TF regulatory analysis through feed forward loops (FFLs)
Comparative studies have identified genes regulated by miRNAs that may have considerable contribution in cellular processes (Friedman, R.C., et al. 2008). The available experimental evidences indicate that miRNAs mediate gene regulation via translational repression with or without mRNA decay. However, the variation in these contributions over time remains undefined and also remained unnoticed in case of TFs (Duk MA, et al. 2014). Therefore, we design a transcription network utilizing one of the most significant motifs i.e. Feed Forward Loop (FFL) motif (Dubitzky W., et al.2013). In our study, the FFL is composed of a transcription factor (EGR1, EPAS1, ETS1, FOS, HIF1A) which regulates the expression of hsa-miR-205, then TFs and miRNA both bind at the regulatory region of the target gene (VEGFA) and jointly modulate its transcription rate. This FFL collectively has three transcription interactions, which could be either activation or repression. In this study, the FFL has two possible structure configurations in both coherent and incoherent FFL specified as type 1 or 2 coherent FFLs and type 1 or 2 incoherent FFLs respectively as shown in Fig. 3(a, b, c, d) where hsa-miR-205 is miRNA and VEGFA is the target protein. In 3(a) represents Type 1 incoherent FFL where TF activates both target mRNA and miRNA synthesis. 3(b) Type 1 coherent FFL, miRNA and TF represses target mRNA but activates miRNA synthesis. 3(c) Type 2 incoherent FFL, TF represses both target mRNA and miRNA synthesis. 3(d) Type 2 coherent FFL, TF activates target mRNA and represses miRNA synthesis.