Abstract ..
Previous studies showed both virus and host-derived MicroRNAs (miRNAs) played crucial roles in the pathology of virus infection. In this study, we use computational approaches to scan the SARS-CoV-2 genome for putative miRNAs and predict the virus miRNA targets on virus and human genome as well as the host miRNAs targets on virus genome. Furthermore, we explore miRNAs involved dysregulation caused by the virus infection. Our results implicated that the immune response and cytoskeleton organization are two of the most notable biological processes regulated by the infection-modulated miRNAs. Impressively, we found hsa-miR-4661-3p was predicted to target the S gene of SARS-CoV-2, and a virus-encoded miRNA MR147-3p could enhance the expression of TMPRSS2 with the function of strengthening SARS-CoV-2 infection in the gut. The study may provide important clues for the mechisms of pathogenesis of SARS-CoV-2.
...
Furthermore, given the importance of miRNA in regulating gene networks or pathways, viruses have evolved mechanisms to degrade, boost, or hijack cellular miRNAs to benefit the viral life cycle.
...
Moreover, recent studies indicate the virus-encoded miRNAs, host-encoded miRNAs, and miRNA targets together form a novel regulatory system between the virus and the host, which contributes to the outcome of infection25.
...
Computational identification of SARS-CoV-2- encoded miRNA
The complete genome sequence of SARS-CoV-2 (accession number: MN908947) was obtained from the National Center for Biotechnology Information (NCBI). The VMir software package (v2.2)26, an ab initio prediction program designed specifically to identify pre-miRNAs in viral genomes, was used to identify and visualize potential pre-miRNAs within the SARS-CoV-2 genome. The program examines all possible pre-miRNAs in the input sequence and consequently identifies a large number of potential candidates. The program slides a window of defined size across the entire sequence of interest, advancing each window by a given step size. The parameters of window and step sizes were set at 500 and 1 nt, respectively, as set in other literature. The Minimum Free Energy of the predicted hairpin structures were predicted with RNAFold algorithm from the Vienna package27. Hairpins with MEF > -20 kCal/mol were removed from the flowing analysis. Then HuntMi software was used to further distinguish miRNA hairpins from pseudo hairpins with a virus-based model28. Maturebayes web server29 was used to identify the mature miRNA within the pre-miRNAs based on sequence and secondary structure information of their miRNA precursors
..
Prediction of miRNA targets and functional annotation
Human 3′- and 5′-UTR, and the Enhancer sequence for five tissues were downloaded from the UTRdb database30 and EnhancerAtlas(35) database, respectively. Human miRNA sequence was downloaded from the TargetScan database. The targets of virus-encoded miRNA and human miRNAs were predicted with TargetFinder31 and miRanda32.
..
Data process for other coronaviruses
The genome sequence of SARS-C o V, M E R S-CoV-2, and HCoV-NL63 were obtained from NCBI with the accession of NC_004718.3, NC_005831.2, and NC_019843.3, respectively. The pipeline for miRNA prediction from the virus genome, miRNA target prediction, and functional annotation was the same as that applied to the SARS-CoV-2 as we described before.
Results and Discussion
Identification of SARS-CoV-2- encoded miRNAs
A total of 808 potential pre-miRNA in the SARS-CoV-2 genome were detected with VMir Analyzer26, an ab initio prediction program designed specifically to identify pre-miRNAs in viral genomes (Figure 1A). Candidate pre-miRNA sequences were further identified using HuntMi28 with the virus model and the sequences with minimum free energy (MFE) > -20 kCal/mol were removed. A total of 45 virus pre-miRNA candidates were finally obtained, with 30 were found in the direct orientation and 15 in the reverse orientation. The pre-miRNAs were evenly distributed across the SARS-CoV-2 genome (Figure 1B). The average length of the pre-miRNA sequence was 78 nt with an average MFE of -28.1 kCal/mol. Finally, the sequence and position of mature miRNA within the pre-miRNA candidate were identified with MatureBayes, resulting in 90 mature miRNAs.
Virus miRNAs suppress host gene expression
The classical mechanism of miRNAs to regulate their target gene is by binding to the 3′ UTR of mRNA to exert negative regulatory effects on gene expression.
Human 3′ UTR targets of the mature viral miRNAs were predicted
and 40 miRNA were predicted to be binding to the 3′ UTR of 73 genes.
Gene ontology analysis revealed that the top functions of the target genes were enriched in the function of
- Notch binding,
- single-stranded DNA endodeoxyribonuclease activity,
- deoxyribonuclease activity,
- cellular response to peptide hormone stimulus and
- regulation of fatty acid metabolic process (Figure2A).
However, consequently, viruses have evolved sophisticated molecular strategies to subvert host cell apoptotic defenses.
A couple of studies have investigated the mechanism by which viruses modulating apoptosis signaling16,38.
The targeted gene annotated to the related GO items included
CHAC1 and
RAD9A among others,
which were targeted by the virus-encoded miRNA, namely MD2-5p and MR147-3p, respectively (Figure 2B).
CHAC1 (cation transport regulator-like protein 1) is a proapoptotic enzyme and a regulator of Notch signaling.
RAD9A is a cell cycle checkpoint protein required for cell cycle arrest and DNA damage repair in response to DNA damage. RAD9A interacts with BCL-2/BCL-xL, the most important genes that regulate cell death, and promote apoptosis39
The suppressive role of the SARS-CoV-2-encode miRNAs on these genes suggested the possible role of them in reducing the host cell apoptotic to subvert host defense. On the other hand, the involvement of the targeted genes in the biological process of cellular response to peptide hormone stimulus and regulation of fatty acid metabolic process might imply the role of them in the pathological damage. Genes enriched in these biological processes included
FOXO3,
ADIPOQ, and
ADIPOR1, among others.
The FOXO family represents a group of transcription factors that are required for many stress-related transcriptional programs including antioxidant response, gluconeogenesis, cell cycle control, apoptosis, and autophagy40. FOXO3 is required for antioxidant responses and autophagy, and altered expression of FOXO3 was observed in Hepatitis C infection and fatty liver41. ADIPOQ is a gene encodes adiponectin, which regulates both glucose and lipid metabolism and exerts an insulin-sensitizing effect in the liver. The binding of adiponectin with its specific receptors, i.e. ADIPOR1/2, induces the activation of a proper signaling cascade that becomes altered in liver pathologies42. Zhang et. al reported that 2–11% of patients with COVID-19 had liver comorbidities and 14–53% cases reported abnormal levels of alanine aminotransferase and aspartate aminotransferase (AST) during disease progression7. However, it is unclear whether the liver damage of COVID-2019 patients is caused directly by the viral infection or by the drug toxicity. The analysis implied the liver dysfunction may be impaired by SARS-COV-2 virus infection through the release of virus miRNA.
Virus miRNA activate host gene expression
Besides the repression role of miRNA on the target genes, miRNAs were also found to target at promoter regions of genes to activate their expression43. Therefore, we also searched the targets of SARS-CoV-2 encoded miRNA on the 5′ UTR of human genes. The 11 virus miRNAs were identified to be bind to the 5'-UTR of 13 target genes, including the binding between MR385-3p and
- TGFBR3 (Transforming Growth Factor Beta Receptor 3)
Accumulating evidence has proved that miRNA can also activate gene expression by targeting enhancers region in the nucleus11. We then downloaded the tissue-specific enhancer sequences from the EnhancerAlta for five tissues which were reported to be affected by SARS-CoV-2 infection, i.e. lung, gut, spleen, liver, and heart. Target prediction revealed that the lung is with the largest number of genes targeted by miRNA on the enhancer, followed by spleen and gut (Figure 3A).
Functional annotation of these targeted genes revealed a notable enrichment in chemokine signaling pathways and cell skeleton related functions, such as actin filament-based process, actin cytoskeleton organization and cytoskeleton organization in the lung (Figure 3B). For example, MR147-5p targeted at the putative enhancer region of
- CXCL16 and
- ARRB2.
cytokine production in the lungs.
Another category of genes was related to the cell cytoskeleton organization, e.g.
- MYH9 (myosin-9) and
- ITGB5 (Integrin beta-5).
For example, MR66-3p was predicted to bind to the enhancer of TNF-α, one of the most important cytokines in the “cytokine storm”.
In the gut, 54 genes were predicted to be enhanced by 34 miRNAs. The most notable target of the virus miRNA is
- TMPRSS2
The most significant enriched gene function is “transport” (Figure 3D). Transporters involved in the movement of ions, small molecules, and macromolecules, such as another protein, across a biological membrane, and play an important role in the maintenance of intestine equilibrium52. In the liver, the virus miRNA mainly regulates genes involved in the function of actin filament severing and regulation of cellular protein metabolic process (Figure 3E). MR198-3p act on the enhancer of
- ADAR, an adenosine deaminases that act on RNA.
- The MYH9 (non-muscle myosin heavy chain 9)
- RXRA
- The putative upregulation of MYH9 and RARA
Virus genome hijacks host miRNA to regulate immune response
Several studies have reported that the virus genome can hijack host miRNA to modulate host biological processes.
For example, the virus transcript HSUR of Herpesvirus saimiri (HVS) base pairs with the host miRNA 27 (miR-27), leading to miRNA degradation in a sequence-specific and binding-dependent manner, and subsequent activation of infected T cells;
HCMV produces a bicistronic mRNA that mediates degradation of cellular miR-17/miR-20a family miRNAs similarly.
A total of 28 human miRNA were predicted to target at the SARS-Cov-2 genome. Since these miRNAs were hijacked by the virus genome, we could suppose that the genes that originally targeted by these miRNA could be affected. There were more than human 800 genes were predicted to be regulated by these miRNA (Figure 4A), and a notable enrichment at the immune system process was observed(Figure 4B). Among the hijacked miRNAs, miR-146 was reported as one of the key modulators of the immune response. Studies in mice have shown that miR-146b-/- mice displayed enlarged spleens, increased myeloid cell populations both in spleen and bone marrow and spontaneously develop intestinal inflammation56. Another miRNA, miR-939, though no target was predicted, was proposed to regulate proinflammatory genes by experiment. Increasing miR-939 levels should restore homeostasis by decreasing inflammatory protein synthesis57.
Our result posed a possibility that the hijack of human miRNA by the SARS-Cov-2 genome might contribute to the abnormal immunity activation in patients with COVID-19.
Virus genomic region targeted by virus and host miRNA
The SARS-Cov-2 genome consists of 11 open reading frames (ORF), including ORF1ab, which encoding 16 nonstructural proteins, followed by those encode structure proteins, i.e. S (spike protein), E (envelope protein), M (membrane protein), and N (nucleocapsid protein), and six accessory proteins(3a, 6, 7a, 7b, 8, and 10)58. Previous studies have shown that miRNAs derived from both virus genome and host can target viral transcript and regulate the virus infection14.
We then explore the targets of both virus and host cell-derived miRNAs on the SARS-CoV-2 genome. There were 27 SARS-CoV-2 encoded miRNA that can target the virus genome (Figure 5A). Forty-three targets were predicted for the 27 miRNAs and most of them were bind to the ORFlab region (Figure 5B), the longest ORF in the SARS-Cov-2 genome.
There were 2 miRNA targets at the 5′ UTR of the virus genome, and 2 at the S gene, which encodes a spike glycoprotein to bind its receptor ACE2 on human cells, and mediates membrane fusion and virus entry(5).
Twenty-eight(28) human miRNAs were predicted to have 30 targets on the SARS-CoV-2 genome (Figure 5B). Most of the human miRNAs were bind to the ORF1ab, where the enzymes for virus replication and translation were encoded. Followed by 7 human miRNA targeted at the N genes of the virus, and 2 at the 5′UTR. It has been reported that miR-323, miR-491, and miR-654 inhibit replication of the H1N1 influenza A virus through binding to the polymerase PB1 gene20. And miR-122 binding to the HCV genome alters the structure of the 5′ UTR in a manner that promotes viral RNA translation. A human miRNA, hsa-miR-4661-3p, was predicted to target at the genomic region 25,296-25,320 within S genes (21,563-25,384), which is the potential 3′ UTR of the S gene transcript (Figure 5C). This observation
suggested a possible repressor role of hsa-miR-4661-3p on the SARS-CoV-2 S gene expression. It might be an example of the antiviral mechanism except for immune response adopted by the host to defense virus.
These results implicated the possible role of both virus and host cell-derived miRNA in the life cycle and pathogenicity of SARS-CoV-2. Since the transcriptome landscape of the SARS-Cov-2 is largely unknown, further researches in characterizing the genomic and transcriptomic features of it may enlighten the role of miRNAs in regulating the virus gene expression and infection.
Mutations in virus genome affect miRNA function.
Mutations accumulated since the SARS-CoV-2 emerged, studies have implicated that the mutations on the virus genome may lead to an attenuated phenotype of SARS-CoV-2. We, therefore, explored whether miRNA encoded by either virus or human contributed to the fitness of SARS-CoV-2 to human.
A total of 810 mutations on the SARS-CoV-2 genome were collected from China National Center for Bioinformation, including 646 SNPs, 19 indels and 145 deletions (Figure 1B). There were 33 pre-miRNAs located within the mutant region, i.e. 30 pre-miRNAs overlapped with SNPs and 4 located in the deletion region, and the number of mutations within each pre-miRNA ranges from 1 to 5. We then examined the effect of genomic mutations on the generation of miRNA in SARS-CoV-2.
MFE is an indicator of the hairpin stability of pre-miRNA, the stability of pre-miRNA hairpin structure will be reduced with the increase of MFE. Compared with the wild-type pre-miRNA, the MFE of the mutant pre-miRNAs were significantly increased (paired Wilcox test p-value = 3.585e-05,Figure 6A), indicating an overall decreased stability of the mutant-genome encoded pre-miRNA. The most affected pre-miRNA is MR369 with MFE increased from -21.7 kCal/mol to -14.7 kCal/mol, which was considered to not be able to form stable hairpin structure after introducing the mutations (Figure 6B). We then focused on the miRNAs with increased MFE scores and those located in the deletion regions. In the lung, human genes targeted by these miRNAs at enhancers were enriched in the function of positive regulation of MHC class I biosynthetic process, which refers to two genes, i.e
- CIITA and
- NLRC5.
The protein encoded by CIITA is located in the nucleus and acts as a positive regulator of class II major histocompatibility complex gene transcription, and is referred to as the "master control factor" for the expression of these genes.
NLRC5 is a member of the NLR family that acts as a transcriptional activator of MHC class I genes. Silencing of NLRC5 resulted in increased IL-6, TNF, CXCL5, and IL-1β secretion，at the same time decreased secretion of the anti-inflammatory cytokine IL-1059.
The intact pre-miRNA is predicted to generate mature miRNA that enhancer the gene expression of both CIITA and NLRC5, resulting in an enhanced inflammation.
While in the mutant genome, the generation of these miRNAs might be disturbed, thus the enhanced effect was relieved.
Additionally, in the spleen, a couple of target genes were enriched in the cytokine-mediated signaling pathway that was missed by the mutant pre-miRNA, including
CCR6, OAS1, POMC, XAF1, TNFAIP3, CANX, NLRC5, SQSTM1, and CD27.
These results suggested that the mutation in the region overlapped with miRNA might contribute to the attenuated phenotype of SARS-CoV-2 by changing the miRNA repertories encode by the virus.
Next, we further explore the alteration of host miRNA target on the mutation virus genome.
Twenty-two human miRNA can bind to the mutant genome, with 16 overlapped with those binding to the wild type genome (Figure 6C). The genes targeted by the newly recruited miRNA were annotated to the function of
epithelial cell differentiation,
actin filament fragmentation and so like,
while the gene targeted by the miRNA specific binding to the wild type genome were annotated to the immuneresponse-related progress (Figure 6D).
The two miRNA as we described before, i.e. hsa-miR-939-5p and hsa-miR-146b-3p, which were involved in maintaining homeostasis by decreasing inflammatory, were predicted to be not capable of binding to the mutant wild type genome.
Taken together, we speculate that the mutant type of virus might relieve the high level of immunity response, even the cytokine storm, which leads to the severe symptom or dead of COVID-2019, but on the other hand, enhance the invasion of the virus by intensifying the regulation of the epithelial differential and actin filament fragmentation, through attracting different batches of host miRNAs.
miRNA view of the difference between SARS-CoV-2 with other human coronaviruses
Comparing with MERS-CoV and SARS-CoV, SARS-CoV-2 appears to have higher transmission rates, but lower mortality1,2. To explore whether and how miRNA implicated in this difference between them, we predicted the miRNAs encoded by MERS-CoV, SARS-CoV, as well as HCoV-NL63, a coronavirus that causes mild symptoms in human and uses the same acceptor with SARS-CoV and SARS-CoV-2 for cellular entry (ACE2), using the same pipeline as applied SARS-CoV-2.
The largest putative virus-encoded miRNA repository, i.e. 212 mature miRNAs, was detected from the MERS-CoV genome (Table 1), followed by SARS-CoV with 168 miRNAs and SARS-CoV-2 with 90 miRNAs. HCoV-NL63 encodes the least number of miRNAs (64 miRNAs). Furthermore, the number of targets on the virus genome binding by the host miRNA follows the same trend (Table 1).
Functional analysis revealed that high enrichment of genes targeted by the virus-encoded miRNAs in the metabolic process was observed in SARS-CoV and MERS-CoV, but not in SARS-CoV-2 and HCoV-NL630 (Figure S1A).
Enrichment in inflammatory and leukocyte was observed in SARS-CoV-2 and MERS-CoV-2, respectively (Figure S1A).
For the genes “hijacked” by virus genome, the functional keywords were enriched in chemotaxis, cytokine, chemokine and immune for the entire four viruses studied (Figure S1B), and HCoV-NL63.
But the miRNAs targeted at the virus genome were largely different from each other (Figure S2). These results implicated that differences in miRNA repositories and miRNA-target interaction might be involved in the pathological difference among the four coronaviruses.
Conclusions
In summary, this silico study revealed a systematic view of the possible role of the virus-encode miRNAs as well as the host-derived miRNAs influenced by the SARS-COV-2 infection, which provide insights from the perspective of miRNA into the biological process involved by the infection (Figure 7).
Impressively, the systematic analysis revealed that the immune response is the function most affected by the virus-infection introduced miRNA change in repository and targets, which might contribute to the immune evasion of the virus as well as the abnormal immunity activation in the host.
The SARS-CoV-2-encoded miRNA is also implicated in the cytoskeleton dynamics that facilitating the virus envision, trafficking within the cell and release.
The virus-encode miRNAs were also predicted to be able to repress the expression of genes involved apoptosis and regulating fatty acid metabolic process by binding at the 3′ UTR of host genes.
Additionally, the host miRNAs hijacked by the virus genome were mostly involved in the immune system process, such as hsa-miR-146b and hsa-miR-939.
Furthermore, the host miRNA hsa-miR-4661-3p was predicted to bind at the potential 3′ UTR of the S gene transcript with the possible role of a repressor on the expression of S gene.
Finally, the genomic mutation-introduced alterations on the coding ability of virus miRNA and the targets of both virus and host miRNA might contribute to the attenuated phenotype of SARS-CoV-2 during its evolution.
Finally, the comparison of miRNA repository and targets difference between SARS-CoV-2 with other human coronaviruses implicated the role of miRNAs in the shared and distinct clinical characteristics of them.
With further experiments to validate the role of candidate miRNAs in vivo, it is promissing to provide deep insight into the mechisms of SARS-CoV-2 pathogenesis from the pespective of miRNA.
Author contributionsXL conceived the project. XL and ZL design the project. ZL, XL, JW, YX, MG, and KM performed analysis. All authors did data analyses and interpretations; ZL and XL prepared and finished the manuscript.AcknowledgmentsThis work was supported by NSFC grant 81671983 and 81871628 to XL, NSFC grant 81703306 to ZL, NSFC grant 81902027 to JW and Young scientist funding from Jiangsu province to JW(BK20171045). Competing interests:Authors declare no competing interests
Inga kommentarer:
Skicka en kommentar