mRNA and miRNA profiling of Zika virus-infected human umbilical cord mesenchymal stem cells identifies miR-142-5p as an antiviral factor

ABSTRACT Zika virus (ZIKV) infection during pregnancy is associated with congenital brain abnormalities, a finding that highlights the urgent need to understand mother-to-fetus transmission mechanisms. Human umbilical cord mesenchymal stem cells (hUCMSCs) are susceptible to ZIKV infection but the underlying mechanisms of viral susceptibility remain largely unexplored. In this study, we have characterized and compared host mRNA and miRNA expression profiles in hUCMSCs after infection with two lineages of ZIKV, African (MR766) and Asian (PRVABC59). RNA sequencing analysis identified differentially expressed genes involved in anti-viral immunity and mitochondrial dynamics following ZIKV infection. In particular, ZIKV-infected hUCMSCs displayed mitochondrial elongation and the treatment of hUCMSCs with mitochondrial fission inhibitor led to a dose-dependent increase in ZIKV gene expression and decrease in anti-viral signalling pathways. Moreover, small RNA sequencing analysis identified several significantly up- or down-regulated microRNAs. Interestingly, miR-142-5p was significantly downregulated upon ZIKV infection, whereas cellular targets of miR-142-5p, IL6ST and ITGAV, were upregulated. Overexpression of miR-142-5p resulted in the suppression of ZIKV replication. Furthermore, blocking ITGAV expression resulted in a significant suppression of ZIKV binding to cells, suggesting a potential role of ITGAV in ZIKV entry. In conclusion, these results demonstrate both common and specific host responses to African and Asian ZIKV lineages and indicate miR-142-5p as a key regulator of ZIKV replication in the umbilical cords.


Introduction
Zika virus (ZIKV) is a mosquito-borne single-stranded RNA virus, which has been largely neglected for more than 60 years since its first isolation in 1947 [1]. The rapid expansion of ZIKV in South and Central America in 2015 has gained medical attention emphasizing the capacity of ZIKV to spread to non-endemic regions [2,3]. Although ZIKV infection is usually asymptomatic and self-limiting in non-pregnant individuals, ZIKV infection during pregnancy has been demonstrated to cause fetal developmental abnormalities and pregnancy loss [4][5][6]. ZIKV replication has been detected in the placenta, amniotic fluid, cerebrospinal fluid and fetal brain tissue [4,[7][8][9]. To date, effective antiviral agents or vaccines to treat or prevent ZIKVassociated disease are not available for human use. Although it is proposed that ZIKV pathogenesis is associated with its broad cell and tissue tropism during the first trimester pregnancy [10,11], the detailed mechanisms of how ZIKV crosses the physical barrier of the placenta are yet to be fully determined.
Placenta is essential in host defence against viral infections, and its protective functions range from a physical barrier of the multinuclear syncytium to the innate and adaptive immune responses [12,13]. The placenta, which is connected to the fetus by the umbilical cord, provides a rich source of mesenchymal stromal or stem cells (MSCs) that are multipotent, selfrenewable, and able to differentiate in vitro to bone, cartilage and adipose tissues [14]. A recent study by Tipnis et al suggested that umbilical cord-derived MSCs (hUCMSCs) possess essential immunomodulatory phenotypes mediated by the production of soluble factors such as tumour growth factor-β (TGF-β), interleukin-10 (IL-10) and other unidentified factors [15]. Additionally, El Costa et al provides evidence using ex vivo model that Asian lineage ZIKV is able to infect and damage the tissue architecture of maternal decidua basalis, fetal placenta and umbilical cord [10]. Although these studies show that the umbilical cord and hUCMSCs, in particular, were highly susceptible to ZIKV infection, the role and contribution of hUCMSCs in ZIKV infection remains to be determined.
MicroRNAs (miRNAs) are non-coding RNAs that are 20-25 nucleotides long. miRNAs have an important role in regulating stem cell self-renewal and differentiation by suppressing the translation of selected mRNAs [16]. During viral infections, miRNAs can be used either by viruses to promote virus replication, or by host cells to alter differential gene expression in response to viral infection [17]. Small RNA-sequencing represents a useful method to study genome-wide miRNA profiles [18], whereas high-throughput RNAseq technology provides a fast, integrative and effective way to reveal the dynamic changes in host gene expression during virus infection [19]. Little is currently known about the role of miRNAs and miRNAregulated networks in ZIKV-induced pathogenesis.
In this study, we used next generation sequencing technology to characterize and compare host mRNA and miRNA expression profiles in hUCMSCs infected with two lineages of ZIKV. ZIKV-induced mitochondrial dysfunction was observed by transmission electron microscopy, confocal microscopy, and mitochondrial flux analysis. Furthermore, we analysed the role of selected mRNA and miRNAs in ZIKV replication.

Cells and viruses
Umbilical cords were collected from healthy donors and mesenchymal stem cells were isolated from the umbilical cords using an in-house established protocol. All methods were carried out in accordance with relevant Institutional Review Board (IRB) guidelines and regulations (Korea University Guro Hospital, Seoul, Korea). In brief, the cord was cleaned with phosphate-buffered saline (Invitrogen, Carlsbad, CA, USA), blood clots were removed, and the cord was dissected into smaller explants and placed on tissue culture dishes in Dulbecco's modified Eagle's medium supplemented with 10% fetal bovine serum (FBS) (Corning Mediatech) and 1% non-essential amino acids (Sigma, St Louis, MO, USA). The cells were allowed to grow out from the explants, expand into monolayers and supplemented with fresh media every other day. hUCMSCs from passage numbers 3-10 were used for all the experiments.
A549 adenocarcinoma cells, JAR choriocarcinoma cells, Hela cells and Vero African green monkey kidney epithelial cells were obtained from the American Type Culture Collection (Manassas, VA, USA). A549 and JAR cells were cultured in RPMI supplemented with 10% FBS and antibiotics, whereas Hela and Vero cells were cultured in DMEM supplemented with 5% FBS and antibiotics. Human neural progenitor cells (hNPCs) were generated as previously described [20][21][22].
ZIKV MR766 (African lineage) and PRVABC59 (Asian lineage) strains were purchased from ATCC and propagated in Vero cells. Viral titers were determined using a standard plaque assay as described previously [20,21]. For the measurement of 50% tissue culture infective dose (TCID50), supernatants from ZIKV-infected hUCMSCs were collected at 24, 48, and 72 hpi. TCID50 was determined using the Spearman-Kärber method and expressed as TCID50 units/mL.

RNA sequencing and bioinformatics analysis
Total RNA was isolated using Direct-zol RNA mini Prep Kit (Zymo Research, Germany) following manufacturer's instructions. RNA sequencing was performed using an Illumina HiSeq 2500 sequencer (Illumina, San Diego, CA) as previously described [18,23,24]. Sequence reads were aligned to the reference genome (Gene bank accession no. hg19) using Tophat (v 1.4.1, Baltimore, MD, USA). The raw RNA-Seq data files were deposited in NCBI's Gene Expression Omnibus (GEO) and are accessible via GEO Series accession number GSE122193.
The gene expression level was calculated based on fragments per kilobase of exon per million mapped reads using Cufflinks v2.1.13 from Ensembl release 72. Functional classification of the genes was performed for gene ontology (GO) in Database for Annotation, Visualization and Integrated Discovery (DAVID). The GO database classifies genes according to the three categories of biological process, cellular component, and molecular function, and predicts the function of the selected genes. Pathway analysis was performed using Kyoto Encyclopaedia of Genes and Genomes (KEGG), allowing the visualization of metabolic pathways and molecular interaction networks captured in the de novo transcriptome [25]. Use of KEGG analysis was approved with the official permission from Kanehisa laboratory.

Small RNA sequencing analysis
Small RNA sequencing libraries were constructed using the NEXTflex Small RNA sample preparation protocol with an initial input of 70 ng of total RNA as previously described [18]. Reads were then mapped onto the human reference genome [26], using bowtie [27]. The reads for the microRNAs were counted using HTSeq [28] with the mod of 'intersection-nonempty' based on miRBase 20 release [29]. Given the read counts of each microRNAs, EdgeR [30] was applied to analyse the differential expression patterns. Differentially expressed microRNAs were identified with a significance of q-value <0.05 and log 2 fold change > ±2. The raw small RNA-seq data files are accessible via GEO Series accession number GSE149983.
miRNA-mRNA integration analysis miRTarVis+ was used as an interactive visual analysis tool for miRNA-mRNA expression profile data [31]. miRTarVis+ uses prediction algorithms that are based on both sequence and expression profile data. miRTarVis+ supports two sequence-based prediction algorithms, TargetScan and microRNA.org, which are two of the most cited miRNA target prediction algorithms.

Quantitative real-time RT polymerase chain reaction (qRT-PCR)
First-strand cDNA was synthesized from 0.5 μg of total RNA using ImProm-II Reverse Transcription System (Promega, Madison, WI, USA) according to manufacturer's instructions. QuantStudio 6 Flex Real-time PCR system (Thermo Fisher Scientific, Waltham, MA, USA) was utilized for cDNA amplification with Power SYBR® Green Master Mix (Invitrogen) under the following conditions: 95°C for 10 min, followed by 40 cycles of 95°C for 30 s and 60°C for 1 min. Relative mRNA levels were determined based on the ΔΔC t method, and normalized against β-actin mRNA.
To validate miRNA expression changes induced by ZIKV infection, miRNAs from hUCMSCs were isolated using miRNeasy mini kit (Qiagen, Germantown, MD, USA). cDNA from miRNA was synthesized using miScript II RT kit (Qiagen) according to manufacturer's protocol. miRNAs were detected using primers specific to miRNA (Qiagen). SNORD68 small nuclear RNA was used as internal control for normalization. qPCR was performed for detection of levels of specific miRNAs using miScript SYBR Green PCR kit (Qiagen). Quantification was carried out on QuantStudio 6 Flex Real-time PCR system.

Transmission electron microscopy (TEM)
Mock (PBS) or ZIKV-infected hUCMSCs were pelleted and washed twice with PBS. Fixation was performed with PBS containing 2.5% glutaraldehyde and 2% paraformaldehyde for 30 min at 4°C. The pellets were rinsed twice with cold PBS, post-fixed in 1% osmium tetroxide, dehydrated in ethanol series, incubated twice with propylene oxide for 20 min and embedded in Epon mixture. Ultrathin sections of 70 nm were obtained using a Reichert-Jung Ultracut E ultramicrotome (Leica, Wetzlar), mounted on copper grids, and counterstained with uranyl acetate and lead citrate. The specimens were observed with a Hitachi H-7600 electron microscope (Hitachi, Japan) at 80 kV acceleration voltage.
Mitochondrial metabolic flux analysis hUCMSCs were plated on XF-96 plates (Seahorse BioSciences, Billerica, MA, USA) at a density of 8,000 cells/well. On the day of mitochondrial flux analysis, media was changed to Krebs-Henseleit buffer at pH 7.4 and incubated at 37°C in a non-CO 2 incubator for 1 h. All media and injection reagents used in this assay were adjusted to pH 7.4. Using the XFp Extracellular Flux analyser, oxygen consumption rate (OCR) was measured as previously described [32].

Transfection of plasmids, miRNA mimics and siRNAs
Tom20-GFP plasmid was kindly gifted by Dr. Woong Sun (Korea University School of Medicine) and ZIKV NS-FLAG tagged plasmids were gifted by Dr. Vaithi Arumugaswami (Addgene plasmid #79639, 79636, and 79640). For transfection of miR-142-5p mimic, cells were seeded in 12-well plates and incubated for 24 h to reach 80% confluency. Cells were transiently transfected with miR-142-5p mimic or a mimic negative control (Bioneer, Daejeon, Korea) using Lipofectamine 2000 transfection reagent (Invitrogen). After 24 h, RNA was extracted to determine transfection efficiency and knockdown by qRT-PCR. For siRNA transfection, cells were seeded in 6-well plates and cultured until 60∼80% confluency on the day of transfection. Transient transfections with control scrambled, ITGAV-, IL6ST or AXL-specific siR-NAs (Bioneer) were performed with RNAiMAX transfection reagent (Invitrogen) according to manufacturer's protocol.

Virus entry assay
Virus entry assay was performed as previously described [33]. Cells were absorbed with ZIKV MR766 for 1 h at 4°C. After washing with cold DMEM, cells were moved to 37°C incubator for 1 h to allow for virus entry. Cells were then treated with acid glycine to inactivate non-internalized virus for 5 min and washed with cold HBSS. Without permeabilization, cells were stained with anti-pan-flavivirus envelope (E) monoclonal antibody (1:200 dilution; Abcam, Cambridge, UK), followed by anti-mouse Alexa 594 conjugated antibody (1:200 dilution; Invitrogen).

Luciferase reporter assay
Wildtype or mutant of 3 ′ UTR sequences of ITGAV or IL6ST were ligated into the pRL-tk reporter plasmid (a kind gift of Dr. Hoseok Song, Korea University School of Medicine). To construct the mutant ITGAV or IL6ST 3 ′ UTR, the sequences that interact with the seed sequence of miR-142-5p were mutated (ITGAV:  from ATACTAAAGACTTTATAACTGCAT-GAACTT to ATACTAAAGACTT CG TAACTG-CATGAACTT; IL6ST: from GGTAGATGAAC ACTT TATACAGTATATATCT to GGTAGATGAA-CACTTTACGCAGTATATATCT). HEK-293 T cells were seeded in 96-well plates and co-transfected with 100 ng ITGAV or IL6ST reporter plasmids and miR-142-5p mimic or negative control using Lipofectamine 2000 transfection reagent. At 24 h post-transfection, the cells were harvested according to the manufacturer's instructions (Promega, Madison, WI, United States), and the luciferase activities were measured using a Dual-Luciferase Reporter Assay System (Promega).

Statistical analysis
Statistical analysis was performed using Graphpad Prism (Graphpad Software, La Jolla, CA, USA). Data are represented by means ± standard deviation (SD) and analysed using Student's t-test to determine the significance of differences between the means of two groups. A two-tailed p value of less than 0.05 was considered significant.

ZIKV infection in hUCMSCs
Human umbilical cord derived MSCs (hUCMSCs) were prepared as described previously. Similar to previously reported studies [10], hUCMSCs used in our study expressed a major MSC marker, CD73 ( Figure  1(A)). In addition, hUCMSCs also expressed another MSC marker, CD90, as analysed by flow cytometry. (Supplementary Figure 1(A)). Cells were infected with ZIKV (MR766 or PRVABC59) at a MOI of 1. As depicted in Figure 1(A), hUCMSCs were permissive to ZIKV infection, as shown by the detection of ZIKV protein production at 24 h post-infection (hpi). We next evaluated ZIKV replication kinetics in hUCMSCs by qRT-PCR at different time points after infection. The African lineage MR766 strain and Asian lineage PRVABC59 strain were used to infect hUCMSCs and placenta-derived choriocarcinoma cell lines (JAR) (Figure 1(B)). Expression levels of viral polymerase NS5 were measured by qRT-PCR at 4, 24, 48, and 72 hpi. ZIKV NS5 transcript levels gradually increased in both MR766 and PRVABC59-infected hUCMSCs in a time dependent manner. Consistent with previous reports of strain-specific differences in viral replication kinetics of ZIKV [34], approximately 2∼3-fold higher levels of ZIKV NS5 expression were detected in response to infection of human choriocarcinoma JAR cells with MR766 than with PRVABC59. As depicted in Figure 1(C), immunostaining of ZIKV E protein further confirmed increased rates of ZIKV infection at 72 hpi in hUCMSCs. We also measured infectivity titers of ZIKV by TCID50 assay using supernatants from hUCMSCs. Similar to qRT-PCR data, our data showed high virus infectivity titers in both MR766-and PRVABC59-infected cells in a time dependent manner (Figure 1(D)). Viral titers were moderately higher in MR766 compared with that of PRVABC59. Next, we measured and compared ZIKV-induced cytotoxicity levels in hUCMSCs and human neuronal progenitor cells (hNPCs). A slight reduction in cell viability was observed at 48 hpi in ZIKV-infected hUCMSCs in comparison to hNPCs (Supplementary Figure S1(B)). Collectively, our data indicate that hUCMSCs are highly susceptible to ZIKV infection.

Comprehensive transcriptome analysis of ZIKVinfected hUCMSCs
Next we used high throughput RNA-sequencing and small RNA-sequencing to analyse differential expressions of mRNAs and miRNAs in ZIKV-infected hUCMSCs (Figure 2(A)). Each sample had relatively high sequencing coverage range. The volcano plot indicates that fold change differences correlate with significance based on p-values (i.e. genes with a large fold change difference had a low p-value in the group-wise comparison) (Supplementary Figure S2).
We compared the host transcriptional signatures in samples collected at early (4 hpi) and late (48 hpi) time points following MR766 or PRVABC59 virus infection to identify differentially expressed genes (DEGs) with a fold change of ± 2 and a q-value < 0.05. Figure 2   whereas top 20 down-regulated genes are listed in Supplementary Table 2. We also analysed top 20 up-or down-regulated DEGs specific to each strain of ZIKV at 48 hpi (Supplementary Table 3), but most of the DEGs were without any previously known functions.
To determine the gene categories broadly affected by ZIKV, Gene ontology (GO)-based bioinformatics analysis was performed. ZIKV-infected cells were differentially enriched with genes associated with the cell cycle, developmental process, and immune response (Supplementary Figure S3). DEGs were also mapped to the Kyoto Encyclopaedia of Genes and Genomes (KEGG) database, which allows mapping of known biological pathways and explicitly labelling of each proven molecular-level participant in these maps based on sequencing results. KEGG identified functional groups that were significantly enriched in response to ZIKV infection. Top 10 upregulated KEGG pathways activated in ZIKV-infected hUCMSCs include activation of "Cell cycle," "DNA replication" and "Oxidative phosphorylation" (Supplementary Table 4). In addition, KEGG analysis of DEGs involved in the phagosome pathway showed that ZIKV infection induced expression of integrins, such as αvβ3 and αvβ5, at early timepoints (Supplementary Figure S4).
In order to confirm the RNA-seq data, qRT-PCR was used to analyse the expression levels of selected DEGs. Similar to data obtained from RNA-sequencing, qRT-PCR validated that ZIKV infection induces up-regulation of interferon gamma induced protein 10 (IP-10), interleukin-6 (IL-6), 2'-5'-oligoadenylate synthetase 1 (OAS1) and interferon-β (IFN-β) genes (Figure 2(D)). Genes that were upregulated during infection with both viruses were mainly associated with IFN signalling. Thus, we next analysed genes responsible for sensing RNA viruses in hUCMSCs. Expression levels of viral RNA sensors, such as retinoic acid inducible gene-I (RIG-I), melanoma differentiation associated gene (MDA5), and toll-like receptor 3 (TLR3) were significantly up regulated after ZIKV infection (Figure 2(D)). Following ZIKV infection, phospho-AKT expression was rapidly upregulated starting at 0.5 hpi and remained increased until 8 hpi. Additionally, in correlation with our RNA-seq and qRT-PCR data, RIG-I protein expression gradually increased over time, suggesting a potential role of RIG-I in sensing ZIKV RNA in hUCMSCs (Figure 2(E)). Given that RIG-I is a well-established viral RNA sensor, we investigated whether the knockdown of RIG-I leads to modulation of ZIKV infection and host gene expression. As shown in supplementary Figure S5, inhibition of RIG-I expression led to enhanced ZIKV NS5 transcript levels. In addition, expression of downstream genes of RIG-I signalling was significantly reduced.

ZIKV infection alters mitochondrial dynamics in hUCMSCs
Our transcriptome analysis demonstrate the differential expression of several genes involved in mitochondrial dynamics following ZIKV infection (Figure 3 (A)). In particular, upregulation of mitochondrial fusion-related genes such as mitofusin-1(MFN1), mitofusin-2 (MFN2), and OPA1 encoding dynamin-related protein, and downregulation of mitochondrial fissionrelated genes like mitochondrial fission protein 1 (FIS1) were observed during RNA-seq analysis. In correlation with the sequencing results, qRT-PCR data also show increased expression levels of MFN1 and OPA1 following ZIKV infection (Figure 3(B)). To investigate the structure of mitochondria in cells infected with ZIKV, we stained hUCMSCs with mito-Tracker Green, a fluorescent dye that accumulates in active mitochondria. Confocal analysis revealed that ZIKV-infected cells exhibited larger and fused mitochondria at 48 hpi, whereas mock-infected control cells showed both fused and fragmented mitochondria (Supplementary Figure S6).
Transmission electron microscopy (TEM) also indicated that ZIKV-infected cells mostly showed fused and elongated mitochondria at 48 hpi (Figure 3(C)), similar to DENV infection [35,36]. Quantification of TEM analysis demonstrated that ZIKV infection induces an increase in mitochondrial length per cell at 48 hpi (Figure 3(D)). To determine whether the alterations in mitochondrial morphology during ZIKV infection affect mitochondrial function, we measured the oxygen consumption rate (OCR), an indicator of cellular respiration, in mock-and ZIKVinfected hUCMSCs. Compared to mock-infected cells, ZIKV-infected cells displayed increased basal levels of OCR, which suggests that ZIKV induces increased respiration in hUCMSCs (Supplementary Figure S7). Oxygen consumption rate (OCR) in mock-infected and ZIKV-infected hUCMSCs was measured in real-time under basal conditions by XF Seahorse analyser. OCR was performed every 7 min following ZIKV infection in real time. We measured basal OCR level within 2 h after ZIKV infection. There was no cytotoxic effect by ZIKV infection at 2 h, causing no change in cell density. To address the relevance of the alterations in mitochondrial dynamics in ZIKV infection, we treated hUCMSCs with various concentrations of Mdivi-1, a mitochondrial fission inhibitor. Mdivi-1 treatment led to a dose-dependent increase in ZIKV NS5 transcript levels, suggesting that mitochondrial fusion may facilitate viral replication (Figure 3(E)). To characterize the detailed mechanisms by which mitochondrial fusion may alter viral replication, expression levels of RIG-I, mitochondrial antiviral signalling protein (MAVS), and interferon regulatory factor 3 (IRF3) in mdivi-1-treated cells were measured by Western blot analysis. Figure 3(F) shows that high dose of mdivi-1 leads to the suppression of RIG-I, MAVS and IRF3 expression levels, suggesting that impairment of mitochondrial fission in ZIKV-infected hUCMSCs can be a part of the immune evasion strategies utilized by ZIKV in order to facilitate viral pathogenesis (Figure 3(F)).

Distinct cellular miRNA expression profiles in response to ZIKV infection
MicroRNAs (miRNAs) are non-coding RNAs known to regulate various cellular process [37]. The Illumina genome analyser platform was used for high throughput sequencing of short cDNA libraries corresponding to 8-to 35-nucleotide RNAs. Only miRNAs that experienced at least q-value <0.05 and log 2 fold change > ±2 were considered as significant. Pie charts were drawn to depict the frequency of different classes of RNA species present in the small RNA datasets. As expected, the majority of small RNA-seq reads were mapped to human encoded miRNAs, as well as long non-coding RNAs and small non-coding RNAs (Figure 4). We then compared the expression patterns of miRNAs in MR766 and PRVABC59 infected hUCMSCs at 4 and 48 hpi. More than 13 million filtered high quality reads, representing more than 500 known miRNAs, were obtained for both MR766  Table  5). Top 10 differentially up-regulated and down-regulated miRNAs in cells infected with MR766 and PRVABC59 are listed in Tables 1 and 2. At 48 hpi, 13 and 18 miRNAs were up-or downregulated by more than 2 fold in MR766-infected cells and 7 and 21 miRNAs were up-or downregulated by more than 2 folds in PRVABC59-infected cells, respectively.

and PRVABC59-infected cells (Supplementary
To further elucidate the function of differentially expressed miRNAs in ZIKV infection, biological analysis of mRNA targets of the significantly regulated miR-NAs was carried out using TargetScan and microRNA.org. Among the top ten most differentially expressed miRNAs, miR142-5p was selected for further investigation following a comprehensive analysis of the predicted targets. The map of miRNA-mRNA network analysis show downregulated miRNAs (square) targeting the upregulated putative mRNA targets (circles) based on miRTarVis+ analysis. Among many predicted targets of miR-142-5p, IL6ST and ITGAV expression and function were further analysed ( Figure 5(A)).
IL6ST is an IL-6 signal transducer, which plays a critical role in regulating myocyte apoptosis, whereas ITGAV is an integrin subunit alpha V, which functions in cell surface adhesion and signalling. To validate the small RNA-seq data of differentially regulated miRNAs, a subset of miRNAs that were commonly downregulated by both MR766 and PRVABC59 strains was selected for further studies using qRT-PCR. In accordance with the small RNA-seq profiles, both miR-142-5p and miR-219-5p decreased significantly upon ZIKV infection ( Figure 5(B)). Consistent with the data from miRNA-mRNA network analysis, qRT-PCR data showed that the expression of miR-142-5p negatively correlated with the mRNA levels of IL6ST and ITGAV in both ZIKV-infected hUCMSCs ( Figure 5(C)).  miR-142-5p regulates ZIKV replication through the target genes IL6ST and ITGAV To confirm if IL6ST and ITGAV are directly targeted and regulated by miR-142-5p, luciferase reporter genes with IL6ST or ITGAV of 3' UTR and the mutant counterparts at the miR-142-5p binding regions were co-transfected with miR-142-5p mimic.
To examine the potential function of miR142-5p during ZIKV infection, we further investigated the effect of miR-142-5p mimic on ZIKV replication in A549 cells. Overexpression of miR-142-5p in ZIKVinfected A549 cells significantly attenuated ZIKV NS5 gene expression, while it markedly reduced IL6ST and ITGAV expressions at the mRNA level ( Figure 6 (C)), while no effects were observed on off-target genes such as integrin β1 (ITGB1), IL6, and TNF-α ( Figure 6(D)). To validate qRT-PCR data, we also examined the effect of miR-142-5p mimic on ZIKV NS1 protein expression by Western blot analysis. We performed quantitative densitometric analysis of Western blot analysis to quantify the protein expression levels. As shown in Figure 6(E), there was more than 50% decrease in ZIKV NS1 expression levels following miR-142-5p mimic treatment compared with that of control mimic treated cells. Therefore, the Western blot data is positively correlated with qRT-PCR data, suggesting significant change in viral replication by miR-142-5p mimic treatment. The release of live viral particles was also quantified by TCID50 assay, and For further evaluation of the effects of IL6ST and ITGAV on viral replication efficiency, IL6ST or ITGAV siRNA was transfected into A549 cells. Given that AXL is a well-known entry receptor for ZIKV, AXL siRNA was used as a positive control. First, qRT-PCR results showed downregulation of both IL6ST and ITGAV expression following specific siRNA treatment (Figure 7(A and C)). Blocking ITGAV by siRNA transfection resulted in significantly suppressed gene expression of ZIKV RNA (Figure 7 (B)). Furthermore, we performed viral binding assay to determine whether miR-142 enhances viral binding/entry through ITGAV. Supplementary Figure S9 indicates that ITGAV-specific siRNA treatment led to diminished ZIKV E binding to cells. This finding suggests a potential role of ITGAV in ZIKV entry. We also investigated the effect of IL6ST on ZIKV replication. Interestingly, knockdown of IL6ST significantly upregulated ZIKV replication, highlighting the possible anti-viral function of IL6ST (Figure 7(C)). The effect of IL6ST siRNA treatment downstream of IL6ST signalling was also investigated by measuring the level of STAT3 expression using western blot. Reduced expression of STAT3 was observed following IL6ST siRNA transfection, suggesting that IL6ST also affects JAK/STAT signalling pathways (Figure 7(D)).

Discussion
In light of the rapid spread of ZIKV epidemic since 2015, a global collaborative effort has been made to understand the mechanisms of ZIKV maternal-fetal transmission. However, to this date, earliest events of ZIKV transmission in the maternal uterine environment have yet to be fully determined. In this study, we identified significantly modulated mRNAs and miRNAs during the course of ZIKV infection in the hUCMSCs using next-generation sequencing. Our data revealed a similar pattern in the regulation of common mRNA and miRNAs between MR766 and PRVABVC59-infected cells. Using bioinformatics analysis coupled with gene expression profiles, we identified various host pathways involved in antiviral immunity and mitochondrial dynamics potentially involved in ZIKV pathogenesis. Further, our data suggest that miR-142-5p along with its target, IL6ST and ITGAV, plays an important role in the regulation of ZIKV replication.
ZIKV can be mainly classified into the African and Asian lineages, the latter of which has been responsible for the recent epidemics. Several groups have noted phenotypic differences between Asian and African lineage ZIKV strains depending on the cell type. In particular, infection with African lineage strains of human neuronal stem cells resulted in a higher rate of infection and virus production, as well as elevated degree of cell death and antiviral response [38][39][40][41][42]. Similar findings were observed in African-lineage ZIKV-infected trophoblasts [43]. Consistent with previous studies, our results also suggest that MR766-infected cells showed relatively higher rates of infection than PRVABC59infected cells (Figure 1). In terms of the host transcriptome, top 20 upregulated DEGs were enriched in functions related to anti-viral immunity. The similarities, in terms of both quantity and magnitude, between the expression profiles of these top 20 upregulated DEGs are likely due to the similar rates of viral replication observed within hUCMSCs. Interestingly, both RNAsequencing and qRT-PCR assay validated that IP-10, IL-6, OAS1, and IFN-β gene expression levels were elevated to a greater extent (1.5∼2 fold) in PRVABC59-infected cells than in MR766-infected cells at 48 hpi (Figure 2(D)). Thus, it is highly possible that enhanced anti-viral response following PRVABC59 infection contributed to the low viral replication efficiency in MR766-infected cells.
Our transcriptome analysis highlighted an upregulation of several genes involved in mitochondrial dynamics such as MFN1, MFN2, and OPA1 (Figure 3 (A)). Mitochondrial dynamics can be a major determinant of a viral infection outcome as they influence cellular homeostasis, metabolism and innate immune signalling [44]. Other flaviviruses such as DENV and Japanese encephalitis virus have been shown to trigger structural and functional alterations in mitochondria. In particular, DENV infection induces mitochondrial elongation by blocking dynamin-related protein 1 triggered mitochondrial fission [35,36]. Specifically, NS4b of DENV induces mitochondrial elongation by blocking Drp1 activation [45]. Similarly, NS2b3 of DENV partially cleaved mitofusin 1 (Mfn1) and 2 (Mfn2), leading to the reduction of IFN activation [36]. In correlation with this, we also observed that NS4b of ZIKV co-localizes with a mitochondrial marker (Supplementary Figure S8) and mitochondria of ZIKV-infected cells had an elongated morphology compared with mitochondria of mock-infected cells (Figure 3). Moreover, real-time measurement of mitochondrial respiration profiling exhibited an increase in oxygen consumption rate in ZIKV-infected cells compared with mock-infected controls. Given that mdivi-1 (mitochondrial fission inhibitor) treatment enhanced viral gene expression and reduced RIG-I-mediated activation of IFN signalling pathways, mitochondrial elongation driven by ZIKV infection may dampen the early IFN response in favour of viral replication. Further studies on which ZIKV proteins alter virusmodulated morphodynamics of the mitochondria are necessary to understand the strategies utilized by ZIKV to subvert the innate immune response.
Few studies have shown that ZIKV infection induces changes in the expression levels of miRNAs, including miR-30e-3p, miR-30e-5p, and miR-17-5p in astrocytes [46], and miR-124-3p in human neural stem cells [47]. Here, our study identified miR-142-5p as a potential contributing factor that controls ZIKV pathogenesis. Transfection with hsa-miR-142-5p mimic in A549 cells effectively reduced ZIKV replication, as evident from the reduced ZIKV NS5 transcript and NS1 protein levels in the presence of the mimic ( Figure 6). In support of our findings, a recent report by Berrien-Elliott et al. suggested miR-142 as an important anti-viral defence player that maintains homeostasis and function of type I innate lymphoid cells [48]. Given that miR-142 is enriched in hematopoietic tissue [49], this miRNA may be involved in the lineage differentiation of hematopoietic cells. Interestingly, miR-142 has gained considerable attention for its essential role in regulating the immune response. For example, Wong et al. has reported that a fever can lead to the induction of miR-142-5p, which targets many immune-related genes, such as IL6ST, TNF, toll-like receptor 2, and prostaglandin E receptor 2 [50]. Considering that miR-142-5p may be a critical modulator of viral pathogenesis, further studies are necessary to understand the molecular mechanisms through which miR-142 mediate cellular response to ZIKV.
We also analysed the potential role of putative miR-142-5p target genes, IL6ST and ITGAV, during ZIKV infection. IL6ST (gp130) is a signal transduction molecule that mediates downstream signalling of IL-6 family of cytokines. Harker et al. emphasizes the importance of gp130 signalling in T cells, in terms of efficient T cell response against viral infection, indicating a potential antiviral role for IL6ST [51]. In accordance, our result also suggest that IL6ST may activate JAK/STAT pathway-mediated antiviral immunity to suppress ZIKV replication. Another target of miR-142 is integrin subunit alpha V (ITGAV). In supplementary Figure 4, KEGG analysis of upregulated DEGs revealed an enrichment of the phagosome pathway and upregulation of integrin αvβ3 and αvβ5 during ZIKV infection at 4 hpi. As highly conserved molecules expressed by almost every cell type, integrins serve as entry receptors of several viruses, including Adenovirus [52], Foot-and-mouth disease virus [53], hantavirus [54], human coxsackievirus A9 [55], human echovirus 9 [56], and WNV [57]. A very recent study by Wang et al. identified integrin avb5 as an entry receptor for ZIKV to establish an infection in neural stem cells [58]. Similarly, our data from ZIKV entry assay confirm that ITGAV is essential for the initiation of ZIKV internalization. It will be interesting to examine whether other integrins play a role in viral binding and entry, and if integrins work together with AXL receptor during ZIKV entry.
Overall, this study provides a comprehensive view of the changes in host mRNA and miRNAs induced by ZIKV infection and highlights the importance of miR-NAs in the discovery and characterization of cellular factors involved in the modulation of virus replication.