Eleutheroside B1 mediates its anti-influenza activity through POLR2A and N-glycosylation.

Influenza viruses represent a serious threat to human health. Although our research group has previously demonstrated the antiviral and anti‑inflammatory activities of eleutheroside B1, a detailed explanation of the mechanism by which it is effective against the influenza virus remains to be elucidated. In the present study, the transcriptomic responses of influenza A virus‑infected lung epithelial cells (A549) treated with eleutheroside B1 were investigated using high‑throughput RNA sequencing, and potential targets were identified using a molecular docking technique, reverse transcription‑quantitative polymerase chain reaction (RT‑qPCR) assay, and DNA methylation analysis. The transcriptomic data revealed that there are 1,871 differentially expressed genes (DEGs) between the cells infected with the influenza virus strain variant PR8, and the cells infected with PR8 and treated with eleutheroside B1. Among the DEGs, RNA polymerase II subunit A (POLR2A; encoding the largest subunit of RNA polymerase II) and mannosidase α class II member 1 (MAN2A1) were selected from the molecular docking analysis with eleutheroside B1. The docking score of Drosophila melanogaster MAN2A1 (3BVT) was 11.3029, whereas that of POLR2A was 9.0133. The RT‑qPCR results demonstrated that the expression levels of host genes (MAN2A2, POLR2A) and viral genes (PA, PB1, PB2, HA) were downregulated following eleutheroside B1 treatment. Bisulfite‑sequencing PCR was performed to investigate whether eleutheroside B1 was able to modify the DNA methylation of POLR2A, and the results suggested that the average proportion of methylated CpGs (‑222‑72 bp) increased significantly following treatment with eleutheroside B1. Taken together, these findings suggested that eleutheroside B1 may affect N‑glycan biosynthesis, the chemokine signaling pathway, cytokine‑cytokine receptor interaction and, in particular, may target the POLR2A to inhibit the production of influenza virus genes.

Introduction

Influenza A viruses cause worldwide outbreaks of influenza and seasonal pandemics, and pose serious risks to public health (1). Influenza A viruses belong to the Orthomyxoviridae family, with a negative single-stranded, segmented RNA genome, and have been categorized into different subtypes on the basis of their hemagglutinin and neuraminidase antigens (2). Their negative-sense RNA genomes have no proof-reading mechanism during replication, so they are extremely error-prone, giving rise to a high mutation rate (3). Therefore, it is a lengthy process to update the matched vaccine for these quickly modifying viruses, and it is necessary to continually develop novel effective antiviral drugs against influenza viruses in order to prepare for the continual seasonal outbreaks.

Currently available anti-influenza virus drugs target the viral life cycle, including amantadine, rimantadine, oseltamivir, zanamivir and peramivir. However, prolonged treatment and the resulting immuno-compromised status of patients lead to increases in drug-resistant mutations among influenza viruses worldwide (4). Previous studies have indicated that the influenza strains H3N2 and pdmH1N1 are resistant to adamantanes (5), and the latest outbreak of the H7N9 virus is also resistant to oseltamivir (6). On the other hand, novel antiviral drugs derived from traditional Chinese medicine (TCM) do not tend to lead to the development of drug resistance among the viruses so easily. Therefore, numerous studies have focused on developing anti-influenza drugs using natural resources, such as traditional medicines, which offer hopeful new prospects for influenza management (1,7).

Coumarin is a fragrant organic chemical compound of the benzopyrone chemical class that is a natural substance found in many plant species, which exhibits a variety of potent pharmacological activities, including antioxidant, antibacterial, anti-inflammatory, antitumor and antiviral activities (8-13). Previous studies by our research group have demonstrated that eleutheroside B1, a coumarin compound, has a wide spectrum of anti-human influenza virus efficacy, with an IC50 value (i.e., the concentration which leads to half-maximal inhibition) of 64-125 µg/ml in vitro, and it inhibited the mRNA expression of several chemokine genes and the influenza nucleoprotein (NP) gene, while exhibiting low cytotoxicity (14).

Over the past few decades, computational chemistry, molecular biology, pharmacognosy and biotechnology have become major scientific areas for the research of natural products. Some modern technologies, including RNA sequencing and molecular docking, have also been used to identify novel molecules for the effective treatment of diseases, and to investigate the underlying mechanisms of action and the specific targets, as well as DNA, RNA, protein and enzyme interactions, associated with natural products (15). RNA sequencing is a genome-wide analytical technology that has been used to analyze the transcriptome of the host response to human or avian influenza virus infection (16,17). A previous study has also used this technology to provide a comprehensive analysis of the pharmacological effects of natural products (18). The molecular docking approach has become an increasingly important tool in pharmaceutical research, and can be used to model the interaction between a small molecule and a protein at the atomic level, which enables the identification of potential drug targets, as well as the ability to characterize the behavior of small molecules in the binding site of target proteins (19). In the present study, the pathway profiles of influenza A virus-infected lung epithelial (A549) cells following eleutheroside B1 treatment were assessed, which has enabled us to gain a comprehensive understanding of the mechanism of eleutheroside B1 activity against influenza A viral infection. From the RNA sequencing results, it was determined that eleutheroside B1 may exert its pharmacological effects on multiple targets, including the immune system, the glycan biosynthesis and metabolism pathways, signaling molecules and their interactions, transcriptional regulation of the chemokine signaling pathway, various types of N-glycan biosynthesis, and cytokine-cytokine receptor interaction. A molecular docking approach was subsequently used to investigate the potential target receptor proteins for eleutheroside B1. Finally, a variety of molecular biology techniques, including reverse transcription-quantitative polymerase chain reaction (RT-qPCR) assay, bisulfite treatment and DNA methylation analyses, were used to evaluate the results from the RNA sequencing and molecular docking experiments.

Materials and methods

Compounds, cells and viruses, and software

Eleutheroside B1 was purified from Sarcandra glabra extract (Si Chuan, China) and characterized by high-resolution mass spectrometry and 1H and 13C nuclear magnetic resonance spectroscopy, as described previously (13). The purity of the compound exceeded 98%, according to analysis by ultra-performance liquid chromatography/time-of-flight mass spectrometry. Eleutheroside B1 was dissolved in dimethyl sulfoxide (Sigma-Aldrich; Merck KGaA, Darmstadt, Germany) as a stock solution of 50 mg/ml, and stored at −20°C until use. A549 cells were purchased from the American Tissue Culture Collection (ATCC; Manassas, VA, USA). The cells were grown in Dulbecco's modified Eagle's medium with 10% fetal bovine serum under standard conditions at 37°C in 5% CO2 humidified air. The influenza virus strain A/PR/8/34 (H1N1) was also purchased from ATCC. The influenza viruses were propagated in the allantoic cavities of chicken eggs. ORCID iD of ChemDraw Ultra 8.0 and SYBYL-X2.1.1 software was used in this study (no. 0000-0003-1628-7416), kindly provided by Dr Xin-An Huang (Tropical Medicine Institute, Guangzhou University of Chinese Medicine, Guangzhou, China).

Cell culture, virus infection and sample preparation

The A549 cells were grown in a monolayer up to 80% confluency and detached from the flask using 10 mM EDTA (pH 7.4) and 0.25% trypsin. The cells were harvested, and 6×105 A549 cells were seeded in 6-well tissue culture plates. On the following day, the cells were washed twice with PBS and infected with A/PR/8/34 [H1N1; 0.1 multiplicity of infection (MOI)] using serum-free medium for 2 h at 37°C. The inoculum was removed, and the cells were treated with or without eleutheroside B1 at a concentration of 100 µg/ml (14). At 24 h post-infection, the cells were lysed in TRIzol reagent (Thermo Fisher Scientific, Inc., Waltham, MA, USA) and stored at −80°C.

RNA isolation, cDNA library construction and sequencing

Total RNA extracts from each sample were obtained with TRIzol reagent, according to the manufacturer's instructions (Thermo Fisher Scientific, Inc.). The total RNA quality was analyzed using agarose electrophoresis (1% gels). The A260/A280 ratio was determined using a NanoDrop spectrophotometer (NanoDrop Technologies; Thermo Fisher Scientific, Inc.). RNA integrity was assessed by Agilent 2100 TapeStation analysis (Agilent Technologies, Santa Clara, CA, USA). An A260/A280 ratio between 1.8 and 2.0 and an RNA integrity number >7 were considered acceptable parameters for RNA integrity. RNA sequencing was performed on an Illumina X-ten RNA-Seq sequence production system (Illumina, Inc., San Diego, CA, USA).

Pathway analysis of differentially expressed genes (DEGs)

In order to obtain a list of DEGs, Gene Ontology (GO) and pathway enrichment analyses were performed. In addition, GO terms, Interpro (protein sequence analysis and classification) terms, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways that were significantly enriched on our list of genes with altered expression (P<0.05) were identified. Additionally, the integrated pathways with statistical values were computed from our list of DEGs using a reactome pathway analysis (http://www.reactome.org/). For over-representation analysis, the number of altered genes associated with a pathway was compared between the treatment and control groups. To determine pathway perturbation, the significance of a particular gene to a pathway was considered in determining the overall impact on the pathway by examining all annotated functions/interactions of the gene.

Targeting and molecular optimization

Two-dimensional [2D(a)] and three-dimensional [3D(b)] structural information of the compound were obtained from the National Centre for Biotechnology Information (NCBI; (https://pubchem.ncbi.nlm.nih.gov/compound/12302276#section=Top), and these were saved as mol.2 files. The 3D-formatted molecule was imported into the Pharmaceutical Target Seeker (PTS; http://www.rcdd.org.cn/PTS/result), and a search was conducted. Molecular modeling was conducted with the SYBYL program package (20). Eleutheroside B1 was used as the template. Its conformation was searched and identified by energy minimization using the MMFF94 force field with the Powell conjugate gradient minimization algorithm, and a convergence criterion of 0.005 kcal/(mol·Å). The MMFF94 charge was used to calculate partial atomic charges. The parameters were set as defaults.

Protein preparation and docking

Through the PTS, select proteins were identified that were potentially matched with eleutheroside B1. Given our goal of attempting to find a potential target of eleutheroside B1 that would be effective against the influenza virus, and considering our previous results, Golgi mannosidase α class II member 1 (MAN2A1) and RNA polymerase II subunit A (POLR2A) were selected as candidate targets. The four-letter PDB code (3BVT) of Drosophila melanogaster MAN2A1 (with resolution at 1.3 Å) and the PDB code (3D9N) of POLR2A (with resolution at 1.60 Å) were retrieved from the Protein Data Bank (PDB; http://www.pdb.org/pdb/home/home.do). Structural optimization of the protein was conducted using SYBYL software. The crystal structure of D. melaogaster MAN2A1 was subjected to geometry optimization using the Swiss PDB viewer v4.1.0 program. The compound structure of eleutheroside B1 was prepared for docking using ChemDraw Ultra 8.0 and SYBYL-X2.1.1 software. The execution of molecular docking analysis and visualization were conducted using Surflex-Dock software. The removal of solvent molecules, hydrogen addition, and the AMBER7 FF99 charges calculation were performed. Furthermore, all default parameters were used for docking simulations by SYBYL Vina analysis.

RT-qPCR assay

A549 cells were seeded in 6-well plates at 37°C in an atmosphere of 5% CO2, then infected with influenza virus (strain A/PR/8/34; 0.1 MOI) and subsequently treated with eleutheroside B1 at a diluted concentration. At 24 h post-infection, the cells were collected for mRNA expression analysis of the host genes (MAN2A2, POLR2A) and the viral genes [polymerase acid (PA), polymerase basic (PB) 1 and 2 and hemagglutinin (HA)] by RT-qPCR. Total RNA was extracted using TRIzol reagent (Invitrogen; Thermo Fisher Scientific, Inc.). Total RNA (1 µg) was reverse transcribed into cDNA using the Prime-Script RT-PCR kit (Takara Biotechnology Co., Ltd., Dalian, China) at 37°C for 15 min, 85°C for 5 sec. The samples were amplified in SYBR Green qPCR buffer (ABP Biosciences, Jiemei long biological Co., Ltd., Wuhan, China). qPCR was performed using an ABI7500 Real-time PCR System (Applied Biosystems; Thermo Fisher Scientific, Inc.) with the following thermocycling conditions: 95°C for 3 min, followed by 40 cycles of 95°C for 5 sec and 60°C for 30 sec (16). Relative gene (PA, PB1, PB2, HA, POLR2A, MAN2A2 and GAPDH) expression levels were calculated using the 2−ΔΔCq method (21). Primer sequences are listed in Table I.

Table I

Primer sequences.

Gene	Primer	Sequence (5′-3′)
PA	Forward	ACACTACAGGGGCTGAGAAA
	Reverse	TGAACGAGAAAATGTGGATG
PB1	Forward	AGTTTTGGTGTGTCTGGGA
	Reverse	TTCGGGTTTGTATTTGTGTG
PB2	Forward	ACCCAGATGAAGGCACAG
	Reverse	TAGAGTCCCGTTTTCGTTTC
POLR2A	Forward	GATGAACTGAAGCGAATGTCT
	Reverse	GTCGTCTCTGGGTATTTGATG
HA	Forward	TGAACAGGGAAAAGGTAGATG
	Reverse C	AGGGAGACCAAAAGCAC
MAN2A2	Forward	GCCCTCATTTTCTGTTTATTG
	Reverse C	TGCCCTATTTACCCATCAC
GAPDH	Forward	GCTGAGTATGTTGTGGAGTC
	Reverse	GCAGAAGGAGCAGAGATGA

POLR2A, RNA polymerase II subunit A; HA, hemagglutinin; MAN2A1, mannosidase α class II member 1.

CpG island (CGI) prediction

A Hidden Markov Model (HMM)-based CGI prediction for the MAN2A1 and POLR2A gene loci was obtained from the browser (http://www.urogene.org/cgi-bin/methprimer/methprimer.cgi) (22).

Bisulfite treatment and DNA methylation analysis

Genomic DNA was extracted from eleutheroside B1-treated and non-treated A549 cells following an influenza virus infection using a DNA Extraction of Cell and Tissue kit (cat. no. GK0122; Generay Biotech Co., Ltd., Shanghai, China). The bisulfite conversion was performed with 1 µg genomic DNA using an EpiTect Fast DNA Bisulfite kit (cat. no. 59824; Qiagen China Co., Ltd., Shanghai, China), according to the manufacturer's instructions. The bisulfite-modified DNA was amplified by PCR with primers that were designed to cover two regions One primer extended from −545 bp upstream and −213 bp downstream to the transcription start site of MAN2A1, whereas the other extended from -199 bp upstream and 139 bp downstream. For the POLR2A gene, 66 CpG sites ranging from −606 to +476 bp relative to exon 1 in a CpG island were selected for subsequent analysis, and for the design of primers. Subsequently, the first-stage PCR product was amplified under the following conditions: 95°C for 4 min; 40 cycles of 94°C for 30 sec, 55°C for 30 sec, and 72°C for 40 sec, followed by a final extension at 72°C for 5 min. PCR products were purified using a PCR purification kit (cat. no. GK2043; Generay Biotech Co., Ltd.). The PCR products were cloned into the TA cloning vector (cat. no. GV6021; Generay Biotech Co., Ltd.). Ten clones from each donor were sequenced using an M13 reverse primer.

Statistical analysis

All data are expressed as the mean ± standard deviation of at least three separate experiments. Differences between two groups were analyzed using Student's t-test, whereas differences between multiple groups were analyzed using one-way analysis of variance followed by Fisher's Least Significant Difference post-hoc test. Statistical analyses were performed using SPSS 18.0 (SPSS, Inc., Chicago, IL, USA). P<0.05 was considered to indicate a statistically significant difference.

Results

Sequencing and alignment

To investigate the molecular mechanism of eleutheroside B1 action against influenza A virus infection, RNA sequencing of influenza A virus-infected human lung epithelial A549 cells, in the presence or absence of eleutheroside B1 treatment, was performed using an Illumina X-ten RNA-seq sequence production system. The RNA sequences were subsequently aligned against the human gene sequences. Three different experimental groups were established: Lung epithelial A549 cells without infection (A549 group), A549 cells infected with influenza strain A/PR8/34/(H1N1) alone (PR8 group), and A549 cells both infected with A/PR8/34/(H1N1) and treated with 100 µg/ml eleutheroside B1 (PR8 + eleu group). At 24 h after the influenza virus infection, samples were harvested and subjected to library construction and RNA sequencing. An excess of 13 million raw reads for each sample were generated from the constructed RNA sequencing libraries using the Illumina X-ten RNA-seq sequence production system. After stringent data filtering, 34.75 M (97.11%; A549 group), 33.97 M (95.93%; PR8 group), or 34.79 M (95.45%; PR8 + eleu group) sequencing reads could be unambiguously mapped against the reference genome (Table II), which confirmed the high quality of the sequencing of these samples, and a satisfactory variable for subsequent analyses.

Table II

RNA-seq overview: Reads mapping quality summary.

Sample name	Total reads	Total bases	Mapped reads	Mapped rate (%)	Proper paired mapped	Singletons	MAPQ≥5 rate (%)	Discordantly mapped
A549	34753376	5.21E+09	33748010	97.11	33748010	0	95.26	0
PR8	33971006	5.1E+09	32589740	95.93	32589740	0	94.10	0
PR8+eleu	34797038	5.22E+09	33215018	95.45	33215018	0	93.40	0

Effect of eleutheroside B1 treatment on pathway regulation in influenza A virus-infected lung epithelial A549 cells

The gene lists generated for influenza A virus-infected human lung epithelial A549 cells in the absence or presence of eleutheroside B1 treatment revealed large differences in the number and the type of genes that were transcriptionally activated. To gain an improved understanding of the host cell response to influenza virus infection at 24 h, the 1,871 DEGs (958 upregulated and 913 downregulated genes) between the PR8 and the PR8 + eleu groups were characterized by performing an enrichment analysis according to biological process, cellular composition and molecular function. GO enrichment analysis indicated that numerous DEGs were involved in molecular function and biological processes, and these results are presented in Table III. The GO enrichment results revealed that eleutheroside B1 recruited molecular function and biological process genes in the host to inhibit the infection of influenza virus. KEGG pathway enrichment analysis suggested that these DEGs were associated with the immune system, glycan biosynthesis and metabolism, signaling molecules and their interaction, transcription (e.g., of the chemokine signaling pathway), various types of N-glycan biosynthesis, cytokine-cytokine receptor interaction, and RNA polymerase function, which are host-response pathways that are activated in response to influenza virus infection (Fig. 1). Twelve DEGs were enriched in cytokine-cytokine receptor interactions. Of these DEGs, bone morphogenetic protein receptor type 1B (BMPR1B), interleukin-6 signal transducer (IL6ST), interferon γ receptor 2 (IFNGR2), interleukin-1 receptor-associated protein (IL1RAP) and erythropoietin receptor (EPOR) have been reported to be involved in the host immune response to virus infections (Fig. 2). BMPR1B is a host gene, the levels of which change during the process of influenza virus infection (23), and this was expressed at a 548-fold higher level compared with the PR8 group. IFNGR2 is the upstream molecule in the interferon γ (IFN-γ) signaling pathway, which is affected by the influenza virus (24). Compared with the PR8 group, there was a 1,369-fold upregulation of the IFNGR2 gene in the PR8 + eleu group. EPOR is crucial for the phosphorylation of signal transducer and activator of transcription (STAT) 5, which acts as an important molecule in the Janus kinase (JAK)-STAT signaling pathway (22). Influenza A virus directly interferes with IFN-γ-stimulated signal transduction via this pathway, which enables the virus to inhibit IFN-γ-inducible gene expression (24). In the present study, levels of the EPOR gene were found to be upregulated by 255-fold. IL1RAP is a co-receptor of the IL-1 receptor (IL1R1) and the IL-33 receptor (ST2), and is required for signaling through both receptor complexes. The inhibition of IL1RAP is also capable of blocking IL-1 signaling, which disrupts the activation of nuclear factor (NF)-κB mediated by IL-1 and IL-33, and may inhibit influenza virus-induced pro-inflammatory reaction (25).

Table III

Significantly enriched GO terms in response to eleutheroside B1.

GO description	P-value	Number of genes
Oxidoreductase activity, acting on the CH-NH group of donors, NAD or NADP as acceptor	0.001378	6
Heterocycle biosynthetic process	0.001468	119
Aromatic compound biosynthetic process	0.001584	118
Organic cyclic compound biosynthetic process	0.001615	119
Oxidoreductase activity, acting on the CH-NH group of donors	0.002303	6
DNA ligase activity	0.002306	4
DNA ligation	0.002306	4
DNA ligase (ATP) activity	0.002306	4
Binding	0.002597	754
Molecular function	0.002715	952
Intra-Golgi vesicle-mediated transport	0.003632	6
Establishment of protein localization to Golgi	0.003993	5
Protein targeting to Golgi	0.003993	5
Retrograde transport, vesicle recycling within Golgi	0.003993	5
Protein localization to Golgi apparatus	0.003993	5
Double-stranded RNA-specific ribonuclease activity	0.004917	4
Ribonuclease III activity	0.004917	4
Methylenetetrahydrofolate dehydrogenase (NADP+) activity	0.004917	4
Negative regulation of transcription, DNA-templated	0.005146	5
Negative regulation of gene expression	0.005146	5
Negative regulation of nucleic acid-templated transcription	0.006515	5
Negative regulation of RNA metabolic process	0.006515	5
Negative regulation of RNA biosynthetic process	0.006515	5
Ligase activity, forming phosphoric ester bonds	0.006757	4
Nucleobase-containing compound biosynthetic process	0.007534	110
RNA phosphodiester bond hydrolysis, endonucleolytic	0.00998	5
Endoribonuclease activity	0.00998	5
NADP biosynthetic process	0.010598	3
NAD+ kinase activity	0.010598	3
Formate-tetrahydrofolate ligase activity	0.010598	3
Regulation of neurotransmitter levels	0.011107	2
Argininosuccinate synthase activity	0.011107	2
Methylenetetrahydrofolate reductase (NAD(P)H) activity	0.011107	2
Gamma-tubulin binding	0.011107	2
Pteridine-containing compound biosynthetic process	0.011677	4
Folic acid-containing compound biosynthetic process	0.011677	4
Folic acid-containing compound metabolic process	0.011677	4
DNA topoisomerase type I activity	0.015169	3
Ligase activity, forming carbon-nitrogen bonds	0.015411	9
Regulation of Ras protein signal transduction	0.016281	16
Ras protein signal transduction	0.016281	16
regulation of cellular process	0.016648	193
Ion binding	0.016903	192
Regulation of small GTPase mediated signal transduction	0.017061	18
Regative regulation of cellular macromolecule biosynthetic process	0.017259	5
Regulation of developmental process	0.017259	5
Metal ion binding	0.017719	188
Endoribonuclease activity, producing 5′-phosphomonoesters	0.018441	4
NADP metabolic process	0.018441	4
Cation binding	0.019467	188
Regulation of biological process	0.019795	193
Negative regulation of biosynthetic process	0.020303	5
Negative regulation of cellular biosynthetic process	0.020303	5
Negative regulation of nitrogen compound metabolic process	0.020303	5
Negative regulation of nucleobase-containing compound metabolic process	0.020303	5
Negative regulation of macromolecule biosynthetic process	0.020303	5
DNA-directed RNA polymerase II, core complex	0.021299	2
Pteridine-containing compound metabolic process	0.022568	4
Exonuclease activity	0.023676	5
Biological regulation	0.02412	201
Heterocycle metabolic process	0.024922	168
Cellular aromatic compound metabolic process	0.025948	168
Quanyl-nucleotide exchange factor activity	0.02599	18
Sulfur compound transmembrane transporter activity	0.027145	3
Organic cyclic compound metabolic process	0.028103	168
Cellular nitrogen compound biosynthetic process	0.028363	129
ARF guanyl-nucleotide exchange factor activity	0.031456	5
Regulation of ARF protein signal transduction	0.031456	5
ARF protein signal transduction	0.031456	5
GDP binding	0.031599	18
Hydrolase activity, acting on ester bonds	0.031838	40
Regulation of intracellular signal transduction	0.031966	19
Transmembrane receptor protein serine/threonine kinase activity	0.032389	4
Regulation of multicellular organismal development	0.032389	4
Taurine transmembrane transporter activity	0.034046	2
Negative regulation of blood vessel morphogenesis	0.034046	2
Calcium-dependent protein binding	0.034046	2
Taurine:sodium symporter activity	0.034046	2
Negative regulation of angiogenesis	0.034046	2
Histamine receptor activity	0.034046	2
Negative regulation of vasculature development	0.034046	2
Neurotransmitter transporter activity	0.034554	3
Neurotransmitter:sodium symporter activity	0.034554	3
Mannose metabolic process	0.034554	3
Organic acid:sodium symporter activity	0.034554	3
Hormone receptor binding	0.0381	4
Wnt signaling pathway	0.0381	4
Nuclear hormone receptor binding	0.0381	4
Molecular function regulator	0.038432	32
DNA biosynthetic process	0.039752	6
Phosphatase regulator activity	0.040667	5
Intrinsic component of plasma membrane	0.040951	7
Integral component of plasma membrane	0.040951	7
DNA topological change	0.042894	3
Aspartate family amino acid metabolic process	0.042894	3
Signal transduction	0.042957	97
Coenzyme biosynthetic process	0.044016	8
Pyridine nucleotide biosynthetic process	0.044349	4
Pyridine-containing compound biosynthetic process	0.044349	4
Nicotinamide nucleotide biosynthetic process	0.044349	4
Anatomical structure morphogenesis	0.044349	4
Single organism signaling	0.044355	97
Cellular response to stimulus	0.045164	111
Signaling	0.046516	97
Transcription, DNA-templated	0.048815	90
IMP dehydrogenase activity	0.048989	2
Negative regulation of developmental process	0.048989	2
Nucleoside phosphate biosynthetic process	0.049767	13
Nucleotide biosynthetic process	0.049767	13

GO, Gene Ontology.

Figure 1

KEGG pathways enriched in response to eleutheroside B1 treatment. Enrichment analysis results of DEGs between the eleutheroside B1 treatment group (PR8 + eleu) and the virus-infected alone group (PR8). DEG, differentially expressed gene; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Figure 2

Enrichment analysis results of DEGs regulated by eleutheroside B1 treatment in cytokine-cytokine receptor interactions pathway. The heat-map shows the 12 DEGs that were regulated by eleutheroside B1 treatment among three samples in the cytokine-cytokine receptor interaction pathway (high levels of expression, red; low expression, blue). DEG, differentially expressed gene.

Influenza A virus infection in the lung is associated with a robust host inflammatory response and sustained elevated levels of inflammatory immune mediators (26), whereas certain natural product compounds have been demonstrated to exhibit anti-inflammatory effects (27,28). The DEG data associated with the chemokine signaling pathway (Fig. 3) confirmed our previous finding that eleutheroside B1 exhibits anti-inflammatory activity against the influenza virus (14). The DEGs identified in the N-glycan biosynthesis pathway suggest other ways in which eleutheroside B1 may work against the influenza virus. In the N-glycan biosynthesis pathway, there are only a small number of glycosyltransferases catalyzing different glycosidic linkages on different N-glycans that exert a crucial role in cellular activity (29). The synthesis of influenza virus hemagglutinins has been reported to be involved in this process (30). In the present study, 12 enriched DEGs were identified in the major N-glycan biosynthesis pathway, 9 enriched DEGs were found in various other types of N-glycan biosynthesis pathways, and 3 genes were found in other glycan degradation pathways (Fig. 4), all of which are involved in glycan biosynthesis and metabolism. These included: UDP-N-acetylglucosaminyltransferase subunit (ALG13), mannosidase α class 2A member 2 (MAN2A2), mannosidase α class 2C member 1 (MAN2C1), chitobiosyldiphosphodolichol β-mannosyltransferase (ALG1), and α-1,2-mannosyltransferase (ALG9). MAN2C1 has dual functions: One is associated with glycan catabolism, whereas the other is in apoptotic signaling, and its downregulation has been reported to delay cancer growth by inducing either mitotic arrest or apoptosis (31). MAN2C1 overexpression also results in synthesis of the non-glucosylated Man9GlcNAc2 oligosaccharide precursor, which is transferred to the protein substrates more inefficiently, thereby resulting in protein underglycosylation (32). In the present study, MAN2C1 was upregulated 9-fold following eleutheroside B1 treatment.

Figure 3

Enrichment analysis results of DEGs regulated by eleutheroside B1 treatment in the chemokine signaling pathway. The heat-map shows the 24 DEGs that were regulated by eleutheroside B1 treatment among three samples in the chemokine signaling pathway (high levels of expression, red; low expression, blue). DEG, differentially expressed gene.

Figure 4

Enrichment analysis results of DEGs regulated by eleutheroside B1 treatment in various N-glycan biosynthesis pathways. The heat-map shows the 9 DEGs that were regulated by eleutheroside B1 treatment among three samples in various types of N-glycan biosynthesis pathway (high levels of expression, red; low expression, blue). DEG, differentially expressed gene.

RNA synthesis is a fundamental process in gene expression, and POLR2A, a housekeeping gene coding for the large subunit of RNA polymerase II (RNAPII) in RNA biosynthesis (33), is a host factor involved in RNAPII-dependent transcription and mRNA processing of the influenza virus (34). The RNA sequencing results of the present study revealed that POLR2A was downregulated by 16,384-fold following eleutheroside B1 treatment.

POLR2A may interact with eleutheroside B1, as determined from the molecular docking assay

Over the course of the last three decades, specific inhibitors of certain enzymes associated with the N-glycan biosynthesis pathway have demonstrated their potential for antiviral and tumor therapy (35). Golgi MAN2A1 is a key enzyme, which acts late in the N-glycan biosynthesis pathway and is a target in the development of anticancer and antiviral therapies (36). Mannostatin A effectively blocked processing of the influenza viral hemagglutinin in cultured MDCK cells, causing the accumulation of hybrid-type protein-linked oligosaccharides, a result which is in agreement with blocking Golgi mannosidase II (36). MAN2A1 inhibitors are also extracted from plant sources and, therefore, the present study aimed to investigate whether MAN2A1 may be a potential target of eleutheroside B1. However, we did not identify the structure of human MAN2A1 and is not available in databases; therefore, the structure of D. melanogaster MAN2A1 was used in the present study. In the initial results, D. melanogaster MAN2A1 was docked with eleutheroside B1. In addition, the crystal structure of human MAN2A1 was established through homology modeling, and this predicted structure was used for docking assays with eleutheroside B1. Furthermore, a previous study by our research group had shown that eleutheroside B1 inhibited the polymerase activity of influenza virus (14). POLR2A is an important host factor that is involved in the polymerase activity of influenza virus, and is downregulated by eleutheroside B1. Therefore, POLR2A was also selected for further investigation.

Small molecules were docked to these selected proteins, and the docking scores were analyzed. In general, a compound that has a high docking score with the protein (Total Score and CScore) indicates that the compound may have potential activities against the targeted receptor (37). From the results of the docking study, it was possible to determine that the selected molecules exhibited good binding affinity towards the homologous proteins in the model. The specific scores are listed in Table IV.

Table IV

Score of molecular docking.

Protein	Total_Score	CSCORE	UNIFIED-CECORER
MAN2A1	11.3029	5	3
POLR2A	9.0133	5	2

MAN2A1, mannosidase α class II member 1; POLR2A, RNA polymerase II subunit A.

To determine its potency, the best-docked conformation of eleutheroside B1 was analyzed, and this molecular compound was found to bind to the inside of D. melanogaster MAN2A1(3BVT) and human POLR2A(3D9N), based upon the hydrogen and hydrophobic interactions. If such features are identified, it demonstrates that the compound in question may be responsible for mediating biological activities. The docked conformations of eleutheroside B1 in the D. melanogaster MAN2A1(3BVT) and human POLR2A(3D9N) structures are illustrated in Figs. 5 and 6, respectively. The interaction analysis revealed that there are 13 hydrogen bonds formed by eleutheroside B1 in the binding site of D. melanogaster MAN2A1(3BVT), and four hydrogen bonds formed by eleutheroside B1 in the binding site of human POLR2A(3D9N), together with two carbon-hydrogen bonds and two conventional hydrogen bonds. Other interactions were also investigated, including van der Waals' forces, which contact the surrounding residues in the binding site of bovine serum antigen, including Arg-893, Leu-925, Arg-565, Leu-568, Arg-777, Glu-753, Trp-531, Arg-427 and Asp-930. Inspection of the molecular interactions revealed that hydrogen bonds and van der Waals' forces were the major factors influencing the stabilities of the eleutheroside B1-MAN2A1 and the eleutheroside B1-POLR2A complexes, which are illustrated in Figs. 5 and 6, respectively.

Figure 5

Docking results of eleutheroside B1 and MAN2A1. (A) The predicted three-dimensional structure of eleutheroside B1 (rendered by sticks) binding to MAN2A1. (B) The binding interaction of MAN2A1 and the surrounding residues, in which hydrogen donors are shown in pink and hydrogen acceptors are shown in green. The rest of the surface is white. (C) Details of the interaction between eleutheroside B1 and MAN2A1, as well as the hydrophobicity, which is shown in a different color, from the highest lipophilic area (brown) to the highest hydrophilic area (blue). MAN2A1, mannosidase α class II member 1.

Figure 6

Docking results of eleutheroside B1 and POLR2A. (A) The predicted three-dimensional structure of eleutheroside B1 (rendered in sticks) binding to POLR2A. (B) The interaction between eleutheroside B1 and POLR2A, as well as the hydrophobicity, which is shown in a different color, from the highest lipophilic area (brown) to the highest hydrophilic area (blue). (C) The binding interaction of POLR2A and the surrounding residues, in which hydrogen donors are shown in pink and hydrogen acceptors are shown in green. The rest of the surface is white. POLR2A, RNA polymerase II subunit A.

MAN2A1(3BVT) in this docking assay belongs to D. melanogaster, which has high sequence identity with human MAN2A1 (41% identity; 61% similarity) (38). A predicted crystal structure for the human MAN2A1 was also generated through homology modeling with docked eleu-theroside B1, however this had a lower docking score (data not shown). Since the protein POLR2A is associated with a higher molecular docking score and is also of human origin, this protein was identified as a target protein for eleutheroside B1, and further experimental verification was performed at the cellular level.

Eleutheroside B1 modifies the CpG methylation levels of the POLR2A and MAN2A1 genes

Even in the molecular docking assay, human MAN2A1 did not match precisely with eleutheroside B1; however, the possibility of their biological interaction could not be excluded. Therefore, MAN2A1 and POLR2A were selected to examine the CpG methylation levels modified by eleutheroside B1. CpG islands were identified on the POLR2A and MAN2A1 promoters by sequence analysis (http://www.urogene.org/cgi-bin/methprimer/methprimer.cgi; Fig. 7). To investigate whether eleutheroside B1 modifies DNA methylation of POLR2A and MAN2A1, bisulfite-sequencing PCR (BSP) was performed to quantitatively analyze the methylation levels of CpG sites in these genes in A549 cells incubated with or without eleutheroside B1 following infection with A/PR/8/34 (H1N1). Ten clones were selected for the analysis of each CpG site. The BSP analysis revealed that the average proportion of methylated CpGs of the POLR2A gene increased significantly (P<0.05) in A549 cells treated with eleutheroside B1 compared with the PR8 group (Fig. 7). However, no significant differences were identified for the proportion of methylated CpGs of the MAN2A1 gene between A549 cells with or without eleutheroside B1 treatment following infection with A/PR/8/34 (H1N1) (data not shown).

Figure 7

Prediction of CpG islands and effects of eleutheroside B1 on DNA methylation of POLR2A in virus-infected A549 cells. POLR2A, RNA polymerase II subunit A. The white circles are unmethylated sites and the black circles are methylated sites. *P<0.05 compared with the PR8 group.

mRNA expression of host genes and virus genes are influenced by eleutheroside B1 treatment in A549 cells

To identify genes involved in novel pathways, such as various types of N-glycan biosynthesis and RNA polymerase function, the expression levels of host genes (MAN2A2, POLR2A) were detected. In glycan biosynthesis and metabolism, MAN2A1, MAN2A2, and MAN2C1 were identified as the important targets. MAN2A1 did not display a potential interaction with eleutheroside B1; MAN2A2, however, was involved in the various types of N-glycan biosynthesis pathways, whereas MAN2C1 was involved in other glycan degradation pathways. POLR2A demonstrated a marked potential interaction with eleutheroside B1 in the present study, and significant levels of expression were observed following influenza infection in A549 cells treated with eleutheroside B1, and influenza infection of the A549 cells in the RNA sequencing studies. Therefore, MAN2A2 and POLR2A were selected for determination of their expression levels. In addition, virus genes (PA, PB1, PB2, HA) in A/PR/8 (0.1 MOI)-infected A549 cells were also subjected to mRNA expression studies. The results demonstrated that the expression of host genes (MAN2A2, POLR2A) and virus genes (PA, PB1, PB2, HA) were decreased following eleutheroside B1 treatment at a concentration of 100 µg/ml (Fig. 8). The expression levels of the influenza virus ribonucleoprotein (RNP) genes (PA, PB1, PB2) were down-regulated by eleutheroside B1, a finding which was consistent with the results of a previous study, showing that eleutheroside B1 may target the influenza virus RNP and RNA polymerases in the host (14).

Figure 8

Effect of eleutheroside B1 on host and viral gene expression. (A and B) mRNA expression levels of the host genes (MAN2A2, POLR2A) and (C-F) the viral genes (HA, PB2, PB1, PA) was examined in the virus-infected A549 cells, with or without eleutheroside B1 treatment. *P<0.05 compared with the PR8 group. MAN2A1, mannosidase α class II member 1; POLR2A, RNA polymerase II subunit A.

Discussion

Since a comprehensive analysis of the mechanisms underlying specific agents of TCM acting against influenza virus at the molecular level is difficult, RNA sequencing approaches may compensate for such draw-backs, and provide a systematic analysis for the study of TCM pharmacology in order to identify the key molecular events that are linked with their efficacy (18).

In the RNA sequencing experiments of the present study, the 1,871 DEGs identified between the PR8 and PR8 + eleu groups were characterized by performing an enrichment analysis according to biological processes, cellular composition, and molecular function. These DEGs were involved in various types of N-glycan biosynthesis, the chemokine signaling pathway, cytokine-cytokine receptor interaction, and RNA polymerase function, which are host responses to influenza viral infection. The biosynthesis of N-glycans is very complex in mammals, as is the glycosylation-mediated quality control of protein folding by N-glycans (39). Glycosylation has become an area of intensive study in recent years, due to its ability to impact viral biology. A previous study demonstrated that N-glycosylation of the influenza hemagglutinin serves an important role in the life cycle of the influenza virus (35). The present RNA sequencing results demonstrated that certain DEGs, including ALG13, MAN2A2, MAN2C1, ALG1, and ALG9, were enriched in glycan biosynthesis and metabolism. MAN2C1 overexpression also results in the synthesis of non-glucosylated Man9GlcNAc2 oligosaccharide precursor, which is transferred to the protein substrates rather inefficiently, thereby resulting in suboptimal protein glycosylation. In the present study, MAN2C1 was upregulated by 9-fold following eleutheroside B1 treatment. This finding indicated that eleutheroside B1 may regulate the N-glycosylation of the influenza hemagglutinin via this pathway, although further studies are required to confirm this phenomenon.

Cytokines are soluble extracellular proteins or glycoproteins that are crucial intercellular regulators and are involved in innate, as well as adaptive, inflammation of host defenses, cell growth, differentiation, cell death, angiogenesis, development and repair processes (40). Cytokines and cytokine receptor interaction networks are regarded as essential components of inflammation (41). In the current RNA sequencing studies, 12 DEGs were enriched in cytokine-cytokine receptor interactions. These DEGs included BMPR1B, IL6ST, IFNGR2, IL1RAP, and EPOR, which have been reported to be involved in host immune response to viral infections. IL-6ST encodes membrane glycoprotein 130 (gp130). It is the common signal transducing component of the IL6 cytokine receptor family, which forms complexes with the ligand-binding receptor subunits of multiple IL6 cytokine family members, serving as an important activator of JAK-STAT signaling. IL-6ST has also been reported to be involved in the host response for the infection of influenza virus (42,43). These results suggested that eleutheroside B1 may inhibit influenza virus through cytokine and cytokine receptor interaction networks, although more evidence is required to conclusively demonstrate this finding.

Chemokines are small soluble molecules that regulate cellular homing through molecular gradients (44). A previously published study reported that the chemokines, C-C motif chemokine ligand 2 (CCL2)/monocyte chemoattractant protein 1 (MCP-1), C-X-C motif chemokine ligand 8 (CXCL8)/interleukin 8 (IL8), CXCL9/monokine induced by interferon-γ (MIG), and CXCL10/interferon γ-induced protein 10 (IP-10), are also associated with mortality caused by human and avian influenza viral infections (45). In the RNA sequencing results of the present study, certain DEGs were also identified in the chemokine signaling pathway, a finding which was consistent with our previous study showing that eleutheroside B1 exhibited anti-inflammatory activity against the influenza virus (14).

The influenza virus polymerase performs numerous functions during the virus' life cycle, suggesting that many cellular factors interact with this complex and are required for the viral genome's transcription and replication. POLR2A is one of the influenza virus polymerase-interacting proteins, and is required for viral replication and transcriptional activity of the viral polymerase (34). In the RNA sequencing results, the levels of POLR2A were downregulated by 21,007-fold [logFC (PR8+eleu/PR8) = −14.3585910971] following eleutheroside B1 treatment. The influenza virus transcribes and replicates its genome inside the nucleus of infected cells, a process that requires viral RNA polymerase, including the three subunits PA, PB1, and PB2, as well as host cell factors to transcribe and replicate the viral genome. The cellular transcription factor POLR2A is required for virus replication and the transcriptional activity of the viral polymerase. A previous study demonstrated that the depletion of POLR2A could reduce by 2-fold the viral polymerase transcription of the influenza virus (34). In the current transcriptome data, eleutheroside B1 treatment also reduced the expression of POLR2A, which demonstrated that POLR2A may be a potential drug target for the inhibition of influenza virus. This result was consistent with our previous finding that eleutheroside B1 inhibited RNP mRNA expression and transcriptional activity of the viral polymerase (14).

Molecular docking methodology explores the behavior of small molecules in the binding pocket of a target protein, and calculates the ability of the compounds to act against a particular target (and their specificity). In modern drug discovery, molecular docking methodology serves an important role in predicting the orientation of the ligand, and provides potential leads for researchers to identify potential drugs and drug targets (46). In the molecular docking results in the present study, eleutheroside B1 was found to bind to the inside of the D. melanogaster MAN2A1(3BVT) and human POLR2A, but this was not the case for the predicted crystal structure of human MAN2A1. These findings indicated that POLR2A may be the potential target for the eleutheroside B1 against influenza virus infections in humans, which was in accordance with the RNA sequencing results.

MAN2A2 is an important enzyme in the N-glycan biosynthesis pathway, and its proper functioning is required for the glycosyltransferases of influenza virus hemagglutinins (36). POLR2A is a housekeeping gene in RNA polymerase II-dependent transcription in hosts, and it is involved in the mRNA processing of influenza virus (33). PA, PB1, PB2, HA are four RNA segments of the influenza virus genome that are crucial for the genome-packaging process of influenza virus (47). Similar to human genes, influenza virus must also pass its genome accurately on to the next generation to ensure the survival of the species (48). The results of the RT-qPCR assay revealed that the mRNA expression levels of MAN2A2, POLR2A, PA, PB1, PB2, and HA, which are principally involved in viral RNP function and virus replication, were all downregulated following eleutheroside B1 treatment. Even in the cells receiving treatment with eleutheroside B1 alone (i.e., no virus), effects on the expression levels of MAN2A2 and POLR2A were also observed. These results indicated that eleutheroside B1 may affect certain factors associated with the N-glycan biosynthesis pathway and mRNA processing in the host, which leads to the inhibition of the production of influenza virus proteins.

Through recruiting proteins involved in gene repression or inhibiting the binding of transcription factors to DNA, DNA methylation is involved in the regulation of gene expression, mRNA splicing, and genomic stability (49). In general, DNA methylation may block the expression of a gene (50). DNA methylation analysis revealed that eleutheroside B1 increased the average proportion of methylated CpGs of the POLR2A gene significantly (P<0.05) in A549 cells with influenza virus infection, compared with the virus control cells. On the basis of our previous study (14) and the results of the present study, it can be hypothesized that eleutheroside B1 may affect the POLR2A gene DNA methylation pattern to inhibit the production of the influenza virus genome. For the D. melanogaster MAN2A1(3BVT), the molecular docking results demonstrated that it was a good match with eleutheroside B1, but this was not true for the human MAN2A1. Even in the DNA methylation analysis, the DNA methylation levels of human MAN2A1 gene were not influenced by eleutheroside B1 treatment. Therefore, the POLR2A gene is more likely to be the target of eleutheroside B1 compared with the MAN2A1 gene. However, additional studies are required in order to confirm the detailed mechanism of action.

In conclusion, using RNA sequencing technology, the present study demonstrated that eleutheroside B1 may inhibit influenza virus via the chemokine signaling pathway, cytokine-cytokine receptor interactions, N-glycan biosynthesis, and RNA polymerase, and that several DEGs involved in these pathways are interconnected as a part of a network, as illustrated in Fig. 9. The changes in the chemokine signaling pathway and in cytokine-cytokine receptor interactions revealed in the present study are consistent with our previous study, which demonstrated that the expression levels of several cytokines were all downregulated by eleutheroside B1 (14). Previous studies have also reported that compounds extracted from herbs exert a role in the inflammatory process, and influence the immune system for the prevention of infectious disease (51,52). In the present study, it was observed that eleutheroside B1 targeted the POLR2A gene in order to affect the expression of viral polymerase. Therefore, eleutheroside B1 may be a promising agent for the prevention and treatment of influenza A virus.

Figure 9

Protein interaction map for differentially expressed genes in target pathways. The map was constructed using the String website (https://string-db.org/).