extfieldsofstudy: Biology Medicine

text:
Plasminogen kringle 5 suppresses gastric cancer via regulating HIF-1α and GRP78

Inhibition of tumour angiogenesis has an important role in antitumour therapy. However, a recent study indicates that antiangiogenesis therapy may lead to glucose-related protein 78 (GRP78) associated antiapoptotic resistance. The present study aims to elucidate the dual effects of plasminogen kringle 5 (K5) on tumour angiogenesis and apoptosis induction by targeting hypoxia-inducible factor 1α (HIF-1α) and GRP78. Co-immunoprecipitation and western blotting were used for examining the ubiquitination of HIF-1α and analysing angiogenesis and apoptosis-associated proteins. K5 promoted the sumo/ubiquitin-mediated proteasomal degradation of HIF-1α by upregulating von Hippel-Lindau protein under hypoxia, resulting in the reduction of vascular endothelial growth factor and thus suppressing tumour angiogenesis. Furthermore, K5 decreased GRP78 expression via downregulation of phosphorylated extracellular-regulated protein kinase, leading to caspase-7 cleavage and tumour cell apoptosis. Blocking voltage-dependent anion channel abrogated the effects of K5 on both HIF-1α and GRP78. K5 significantly inhibited the growth of gastric carcinoma xenografts by inhibiting both angiogenesis and apoptosis. The dual effects suggest that K5 might be a promising bio-therapeutic agent in the treatment of gastric cancer, particularly in patients who exhibit the induction of GRP78.

Gastric cancer is an aggressive malignancy that is frequently diagnosed at an advanced stage with poor prognosis. Although surgery and/or a combination of chemotherapy improve the survival rates, the 5-year relative survival rates of the patients receiving these treatments remains low at 30% and that of the patients with advanced disease is o1 year. 1,2 Therefore, it is necessary to develop more effective therapeutic strategies.
As the growth and survival of most tumour cells is dependent on oxygen and nutrient supply, angiogenesis has been a long-standing and attractive therapeutic target for treating malignant tumours. 3 Angiogenesis is induced by hypoxic conditions and regulated by hypoxia-inducible factor 1α (HIF-1α) and its major target gene, vascular endothelial growth factor (VEGF). 4 The expression of HIF-1α and VEGF is increased in many human cancers. 5,6 Close association between high-level expression of HIF-1α and VEGF and patient mortality have been demonstrated in cancers of brain, breast, cervix, oropharynx, ovary, and uterus. 7 Preclinical studies have shown that inhibition of VEGF pathway impedes tumour growth. VEGF-neutralizing antibody Avastin (bevacizumab) and VEGF receptor tyrosine kinase inhibitors (sorafenib and sunitinib) inhibit primary tumour growth in clinical applications and have been used as anticancer treatments in several tumour types. 8 However, clinical observations indicate that these therapies may have limited efficacy because the response has been mostly transient, in addition to the development of drug resistance. 9 Glucoserelated protein 78 (GRP78) could be induced by severe glucose and oxygen deprivation resulting from antivascular and antiangiogenesis therapies, which could lead to drug resistance in an HIF-1α-independent manner. 10 In conditions of stress or hypoxia associated with oncogenesis, GRP78, a major endoplasmic reticulum chaperone, has an essential role in counteracting the apoptosis, which promotes cancer cell proliferation, survival, metastasis, and resistance to a wide variety of therapies. 11,12 K5, the fifth kringle domain in human plasminogen, has been shown to display the most potent inhibitory effect on endothelial cells and angiogenesis. 13 Some studies have demonstrated the antiangiogenesis and antitumour effects of K5 both in vivo and in vitro. 14,15 Our previous study also showed that downregulation of HIF-1α by K5 is responsible for the decreased VEGF expression in endothelial cells and lung cancer cells. 16 Furthermore, K5 has been shown to directly induce apoptosis and cause cell cycle arrest in proliferating endothelial cells and stressed tumour cells through cell surface voltage-dependent anion channel (VDAC) 17,18 or GRP78. 19,20 However, the receptors and the downstream signals by which K5 downregulates HIF-1α and GRP78 in cancer cells have remained unknown. The present study was designed to investigate the effects and the mechanism underlying these effects of K5 on both HIF-1α-induced angiogenesis and GRP78-dependent apoptosis resistance in gastric cancer cells.

Results
High expression levels of CD34, HIF-1α, VEGF, and GRP78 in human gastric cancer. To evaluate angiogenesis-related molecules and antiapoptosis molecule, GRP78, in gastric cancer, we first examined microvessel density (MVD; CD34), HIF-1α, VEGF, and GRP78 in 42 cases of gastric carcinoma. The clinical-pathological characteristics observed in immunohistochemical analyses are shown in Figure 1a. The levels of CD34 and HIF-1α in the tumour and para-tumour groups and in the tumour and normal groups are significantly different (Table 1). The difference in VEGF expression between tumour and normal groups is statistically significant, whereas there was no significant difference between tumour and para-tumour groups and between para-tumour and normal groups. The differences in the expression of GRP78 between tumour and normal groups and between tumour and para-tumour groups are statistically significant, whereas there was no significant difference between normal and para-tumour groups (Table 1). These results suggested that tumour angiogenesis was stimulated by HIF-1α and VEGF, and GRP78 may facilitate human gastric cancer growth.
K5 inhibits angiogenesis in heterotopic tumours in athymic mice, resulting from human gastric cancer transplants. To evaluate the effect of K5 on tumour growth, heterotopic tumour xenografts of human gastric cancer were established in athymic mice. The animals received intraperitoneal injection of K5 or PBS and the tumour growth was monitored ( Figure 2a). Statistical analysis indicated that tumour weight was significantly lower in the K5-injected mice compared with that in the PBS-treated mice. The mean suppression of primary tumour growth in the K5-treated group was 67% (Po0.01, n = 6, Figure 2b). To understand the mechanism underlying the suppression of tumour growth, tumour sections were analysed for MVD by staining with anti-CD34 antibody. In contrast to the PBS group, K5-treated group showed a dramatic reduction in MVD (Po0.01; Figure 2c).
K5 downregulates the expression of HIF-1α and VEGF in vivo and in vitro. A large body of data indicates that HIF-1α and VEGF contribute to tumour angiogenesis. To determine whether HIF-1α-VEGF signalling pathway has a role in tumour angiogenesis, tumour sections were analysed for VEGF and HIF-1α. The results of immunohistochemical staining and western blotting analysis (Figures 2d and e) showed that the expression of VEGF and HIF-1α were downregulated by K5. To further explore the effect of K5 on VEGF and HIF-1α, in vitro experiments were performed using SGC-7901. The results indicated that the expression of VEGF and HIF-1α in SGC-7901 cells was upregulated under hypoxic conditions. K5 markedly suppressed the expression of VEGF and HIF-1α in a dose-dependent manner (Figures 2f  and g). The downregulation of VEGF does not contribute to the pro-apoptotic effects of K5, because there are no VEGF receptors (Flt-1 and KDR) on SGC-7901 cells (Supplementary Figure S1). These data suggested that K5 inhibited the neovascularization of human gastric cancer by downregulating the VEGF and HIF-1α expression.
K5 promotes the SUMO/ubiquitin-mediated proteasomal degradation of HIF-1α. To investigate whether VEGF is regulated by HIF-1α in SGC-7901 cells, specific siRNA targeting HIF-1α was used. Knockdown of HIF-1α inhibited the accumulation of VEGF protein dramatically (Figure 3a). We found that HIF-1α accumulated significantly after MG132 treatment ( Figure 3b). Moreover, the downregulation of HIF-1α by K5 under hypoxia was blocked by MG132 ( Figure 3b). These observations suggested that the Figure 1 Expression of CD34, HIF-1α, VEGF, and GRP78 in normal human gastric tissues, para-tumour tissues, and tumour tissues. Immunohistochemical staining of paraffin-embedded serial tissue sections from 30 normal gastric tissue specimens, 17 gastric para-tumour specimens, and 42 gastric carcinoma specimens, with antibodies against CD34, HIF-1α, VEGF, and GRP78. Immunoreactivity is indicated by brown staining. The red arrow represent positive, and the black arrow represent negative. Scale bars represent 20 μm downregulation of HIF-1α by K5 may be through ubiquitinproteasomal degradation pathway. Recent data showed that SUMOylation of HIF-1α facilitated polyubiquitination, leading to the proteasomal degradation of HIF-1α. 21 To examine whether K5 has the same role in this process, HIF-1α was precipitated by a specific anti-HIF-1α antibody. We observed that K5 promoted the SUMOylation and ubiquitination of HIF-1α, which eventually resulted in its degradation ( Figure 3c).

K5 increases VHL expression and promotes HIF-1α/VHL binding via VDAC.
Our results also showed that K5 upregulated von Hippel-Lindau protein (VHL) (Figure 3d), which could bind to HIF-1α and function as an E3 ubiquitin ligase in the proteasomal degradation of HIF-1α under hypoxia. 22,23 To further confirm the specific contribution of VHL to K5-induced proteasomal degradation of HIF-1α, coimmunoprecipitation (co-IP) assay and VHL-specific knockdown were performed. The results showed that K5 promoted the binding of VHL to HIF-1α ( Figure 3e). Furthermore, the effect of K5 on HIF-1α was abolished when VHL was specifically knocked down ( Figure 3f). However, we found no changes in prolyl hydroxylase domain proteins 1 (PHD1) and PHD3. PHD2 did exhibit an increase but the increase was not statistically significant (Supplementary Figure S2), suggesting VHL, rather than PHDs, could be the target of K5. Earlier studies have proposed VDAC as the receptor of K5 in K5-induced endothelial cell apoptosis. 17,18 We therefore tried to demonstrate the role of VDAC in VHL regulation. Our results showed that specific knockdown of VDAC using siRNA abrogated the induction of VHL by K5 (Figure 3g). Taken together, these data strongly suggest that K5 promoted the SUMO/ubiquitin pathway-mediated proteasomal degradation of HIF-1α by inducing VHL and binding to it, which is mediated by VDAC. K5 downregulates GRP78 in vivo and in vitro. GRP78 is induced in a wide variety of cancer cells and has an important role in inhibiting apoptosis. 12 Our results showed that GRP78 was downregulated by K5 in the xenografts (Figure 4d). The effect of K5 on GRP78 expression was further confirmed by RT-PCR. The result showed that GRP78 mRNA was decreased in a time-dependent manner ( Figure 4e). Consistent with the changes in mRNA levels, both immunocytochemistry and western blotting analysis indicated that the protein level of GRP78 was upregulated under hypoxic condition in SGC-7901 cells and reduced significantly by K5 in a dose-dependent manner (Figures 4g and g).
Pro-apoptotic effect of K5 is mediated via downregulation of GRP78 in SGC-7901 cells under hypoxia. To confirm the relevance of the effects of K5 to apoptosis and GRP78, GRP78 was knocked down using specific siRNA. We observed that GRP78 knockdown significantly increased the apoptotic rate in SGC-7901 cells (Figure 5a), accompanied by the cleavage of pro-caspase-7 ( Figure 5b). K5 induced the cleavage of pro-caspase-7 ( Figure 5c). These data confirm that the pro-apoptotic effect of K5 is mediated via downregulation of GRP78 in SGC-7901 cells under hypoxia.
K5 downregulates GRP78 via VDAC and ERK phosphorylation. To identify the possible additional signalling pathways involved in the apoptosis induced by K5, the total protein levels and the phosphorylation status of proteins involved in GRP78 synthesis were examined. The results showed that K5 downregulated p-ERK at 1 h ( Figure 6a). We then used U0126 to evaluate the roles of p-ERK and GRP78. The results indicated that U0126 suppressed p-ERK protein expression, along with the downregulation of GRP78 (Figure 6b). In addition, neither Jun N-terminal kinase (JNK) nor p38 was necessary for the K5-induced apoptosis (Supplementary Figure S4). VDAC was specifically knocked down to validate the role of VDAC in the process of K5-  We found that specific knockdown of VDAC prevented the downregulation of GRP78 by K5 (Figure 6c). Furthermore, the downregulation of p-ERK was also abolished by siVDAC (Figure 6d). The results of all these experiments together demonstrate that VDAC has a vital role in K5-induced gastric cell apoptosis via the p-ERK-GRP78 pathway.

Discussion
In the present study, we confirmed the antitumour activity of K5, which significantly inhibits gastric cancer growth via the dual effects of antiangiogenesis and pro-apoptosis. Importantly, we found that K5 promoted the SUMO/ubiquitinmediated proteasomal degradation of HIF-1α by inducing Tumour angiogenesis is the core process in tumourigenesis. One of the major drivers of tumour angiogenesis is hypoxia, which characterizes the tumour microenvironment. At the molecular level, the master switch orchestrating the cellular response to low O 2 environment is generally considered to be the transcription factor HIF, which acts as the angiogenic switch and increases the formation of pro-angiogenic factors, such as VEGF. 24,25 Some studies have indicated that gastric carcinoma patients with VEGF-positive tumours had worse prognosis than those with VEGF-negative tumours. 5, 6 We examined in this study gastric carcinoma tissues from 42 patients and found that CD34 (MVD), VEGF, and HIF-1α were significantly higher in these tissues than those in para-tumour or normal tissues (Figure 1). Our results also demonstrated that K5 significantly inhibited the expression of HIF-1α, VEGF, and CD34 in xenograft gastric tumours in mice and in SGC-7901 tumour cells (Figure 2). Knockdown of HIF-1α led to the corresponding reduction of VEGF (Figure 3a). These data strongly supported the antiangiogenic effect of K5 in gastric tumours. This result is consistent with that of our previous studies, which showed that K5 suppressed the growth of hepatocellular carcinoma by targeting tumour angiogenesis. 14 The hypoxic microenvironment contributes to tumour progression through the stabilization of the potent transcriptional factor HIF-1α. Under normoxia, HIF-1α is hydroxylated by the oxygen-sensitive PHDs, which leads to the binding of E3 ubiquitin ligase, VHL, to HIF-1α in a complex with elongin B, elongin C, and Cul2. 22,23,26 When cells are exposed to a hypoxic microenvironment, this hydroxylation-mediated degradation pathway is blocked and results in HIF-1α translocation and accumulation in the nucleus, where it binds to the constitutively expressed HIF-1β to form a heterodimer and activates hypoxia-responsive genes, such as VEGF and Glu-1. 27,28 Consistent with previous results, HIF-1α accumulated under hypoxia ( Figure 3d). As PHDs and VHL have an important role in the proteolytic degradation of HIF-1α, we next assessed which of these contributed to the K5-induced proteasomal degradation of HIF-1α. We found that K5 treatment under conditions of hypoxia upregulated VHL and downregulated HIF-1α (Figure 3d). The regulation of VHL by K5 in normoxia condition was measured; VHL was upregulated a little but much less than in hypoxia. The reason may come from the high constitutive expression of VHL under normoxia; the regulation under normoxia was not as obvious as under hypoxia. Our previous study showed similar results that K5 downregulated VEGF much significantly in hypoxia than in normoxia in the oxygen-induced retinopathy model. 29 The common features of diabetic retinopathy and tumour microenvironment are hypoxia, therefore we speculate that the effects of K5 on VHL, HIF-1α, and VEGF are more significant and suitable for the treatment of the diseases with hypoxia microenvironment (Supplementary Figure S5). However, there was no change in the levels of PHD1 and PHD3. PHD2 demonstrated an increasing trend, but the increase was not statistically significant (Supplementary Figure S2), suggesting VHL, rather than PHDs, could be the target of K5. We further demonstrated that K5 facilitated the binding of HIF-1α to VHL and subsequent degradation of HIF-1α under hypoxia (Figure 3e). A recent study showed that SUMOylation of HIF-1α, which promotes its binding to VHL, leads to its ubiquitination and degradation via PHD-independent mechanism. 21 Similarly, we observed that K5 treatment led to the decrease in HIF-1α by promoting HIF-1α SUMOylation and ubiquitination (Figure 3c). These results suggested that K5 stimulated the SUMOylation of HIF-1α and upregulation of VHL, which promoted HIF-1α/VHL binding and consequent ubiquitination and degradation of HIF-1α.  The receptor that mediates K5-induced VHL upregulation remains unknown. Several groups have independently demonstrated that VDAC serves as the receptor for K5 in human endothelial cells, 1-LN human prostate tumour cells, and retinal Müller cells. 18,30 Recently, our group showed that VDAC was upregulated by K5 in endothelial cells and served as the plasma membrane receptor to amplify the apoptotic effect of K5, acting as a 'positive feed-back loop'. 17 To investigate whether VDAC mediated the regulation of VHL by K5 in gastric cancer cells, we blocked VDAC expression by specific siRNA and the results showed that the effect of K5 on VHL was abolished (Figure 3g), which indicated that VDAC acts as the receptor in gastric cancer cells and mediates the antiangiogenesis effect of K5.
Antivascularization and antiangiogenesis therapies have been used to treat cancer. However, recent clinical observations indicated that antiangiogenesis agents such as Avastin and Sorafenib may have limited efficacy because of the induction of apoptosis resistance. 10,25 Interestingly, we found that K5 could induce apoptosis in SGC-7901 gastric cancer cells and gastric tumour xenografts in mice, in addition to exerting an antiangiogenesis effect ( Figure 4). Downregulation of VEGF does not contribute to the pro-apoptotic effects of K5 because of the absence of VEGF receptors (Flt-1 and KDR) on SGC-7901 cells (Supplementary Figure S1). However, antiangiogenesis may lead to severe glucose and oxygen deprivation as well as GRP78 expression, which are HIF independent and closely related to apoptosis resistance. 10 It is well known that GRP78 promotes tumour proliferation, survival, and metastasis. 11 GRP78 overexpression confers resistance to a wide variety of chemotherapeutic agents in multiple tumour types, such as lung, bladder, stomach, breast, gastric, and epidermoid carcinoma. Knocking down GRP78 sensitizes the tumour cells to drug treatment. 12,31,32 Our results showed that the level of GRP78 was elevated in all the 42 human gastric tumour specimens included in this study (Figure 1). Silencing GRP78 for 48 h resulted in extensive apoptosis of SGC-7901 cells and increased the levels of cleaved caspase-7 (Figures 5a and  b), which indicates that GRP78 had a key role in the apoptotic process in gastric cancer. K5 could significantly downregulate GRP78 expression in vivo and in vitro (Figure 4), increasing the cleavage of caspase-7 (Figure 5c). These results are consistent with the recent reports that GRP78 protects cells from apoptosis during stress by binding to procaspase-7, blocking its activation, and decreasing stress-induced cellular apoptosis. 20,30 Zhang et al. 33 and Song et al. 34 indicated that activation of ERK by ER stress such as chronic hypoxia is necessary for the induction of GRP78, which protects gastric cancer cells against apoptosis. 33,34 As shown in Figures 6a and b, phosphorylation of ERK1/2 and expression of GRP78 were attenuated by K5 treatment and a specific inhibitor of ERK1/2. Neither JNK nor p38 was necessary for the K5-induced apoptosis (Supplementary Figure S4). Next we wanted to identify the receptor that mediated the downregulation of GRP78 by K5. As VDAC is a crucial regulator of apoptosis in many cancer types and mediates K5-induced apoptosis of endothelial cells, 17,18,35,36 VDAC was evaluated as the likely receptor of K5 in its regulation of GRP78. We used specific siRNAs to knockdown VDAC and found that the effects of K5 on GRP78 and p-ERK were abolished (Figures 6c and d), which highlighted the role of VDAC, as a critical receptor of K5, for both the antiangiogenesis and pro-apoptotic effects. GRP78 has also been identified as a receptor of K5, 19,20 and hence, whether GRP78 also acted as a receptor for K5 in our study needs to be investigated. Some studies also showed the co-localization of VDAC and GRP78 on the surface of 1-LN cells and endothelial cells. 30,37 It is possible that the downregulation of GRP78 by K5 results in the reduction of the cell surface GRP78, thus promoting the binding to VDAC, mediating the effects of antiangiogenesis, and pro-apoptosis.
In summary, the present findings provide in vivo and in vitro evidence for the antiangiogenic and the pro-apoptotic effects of K5 in gastric cancer (Figure 7). The mechanism underlying these effects, downregulation of GRP78 and HIF-1α by K5, indicates that K5 may be a promising agent to counter the drug resistance of antiangiogenic therapies that result in severe glucose and oxygen deprivation.   Subcutaneous implantation of SGC-7901 cells was performed by injecting 1 × 10 7 cells into the lower right flank of athymic nude mice. When tumour volume reached about 50 mm 3 , the mice were randomized into two groups and received intraperitoneal injection of K5 or PBS. Human K5 was obtained as described previously. 14,15 The K5 group received K5 injection at the dose of 2.5 mg/kg per mouse, every 3 days, until the total dose reached 12.5 mg/kg. Control group was treated with the same volume of PBS. The tumours were measured on alternate days. Tumour volumes were calculated by the equation: volume = (length × width 2 ) × 0.5. Tumours were collected 22 days after the injection of the cells and weighed.
Immune staining assay. Immunocytochemistry and immunohistochemical analysis were performed as described earlier. 38 Briefly, cells were seeded on gelatine-coated coverslips in six-well plates and treated with K5. The tissues were fixed in 4% paraformaldehyde solution and prepared as paraffin-embedded samples for immunohistochemical analysis. The samples were incubated with primary antibodies against CD34 from human epithelial cells (1 : 100, Dako, Copenhagen, Denmark), mouse CD34 (1 : 50, Abcam, Cambridge, UK), VEGF (1 : 100, Santa Cruz Biotechnology, Santa Cruz, CA, USA), HIF-1α (1 : 50, Santa Cruz), or GRP78 (1 : 100, Santa Cruz) at 4°C overnight. Avidin-biotin-peroxidase was used for detecting immunolabelled cells and the sections were counterstained with haematoxylin. Primary antibodies were omitted in the staining of negative controls. Tumour MVD was quantified by the Weidner's method. 39 The images were captured by an Olympus microscope (Tokyo, Japan) at the magnification of × 400. ImagePro Plus 6.0 (Bethesda, MD, USA) was used to analyse the results.
Quantitative RT-PCR. qRT-PCR was performed as described earlier. 40 Briefly, total RNA was isolated using Trizol (Invitrogen, CA, USA) and reverse-transcribed (Takara, Shiga, Japan). GRP78 cDNA was amplified using primers 5′-GGA GAAGACTTTGACCAGCG-3′ and 5′-CTTTGGAATTCGAGTCGAGC-3′. Amplification reactions were performed with an initial denaturation step at 94°C for 2 min, followed by 30 cycles consisting of a denaturation step at 94°C for 5 min, an annealing step at 55°C for 30 s, and an extension step at 72°C for 1 min. To estimate the efficiency of cDNA synthesis, β-actin was used as a control. The PCR products were resolved on 1% agarose gels.
Co-immunoprecipitation. For immunoprecipitations, 500 μg SGC-7901 cell lysates were incubated with 4 μg of anti-HIF-1α (Santa Cruz) over night at 4°C, followed by the addition of 20 μl Protein A or Sepharose G beads (Calbiochem) and incubated for 4 h at 4°C. Immunoprecipitates were washed four times with lysis buffer, eluted with sample buffer, and immunoblotted.
In situ detection of apoptosis by TUNEL staining. Frozen tumour sections (5 μm thick) were analysed by TUNEL staining, using an In Situ Cell Death Detection Kit (Fluorescein; Roche Diagnostics, Indianapolis, IN, USA). Apoptotic cells appear as green-stained cells. The total number of cells, indicated by DAPI staining, was calculated by counting the cells in five arbitrarily selected fields at × 400 magnification. The apoptotic index was calculated using the equation, apoptotic cells × 200/total number of cells.
Cell viability and apoptosis assays. SGC-7901 cells were seeded in 24well plates at a density of 2 × 10 4 cells per well and maintained in the growth medium until they reached 60% confluence. Cells were then treated with different concentrations (0, 80, 160, 320, 640, and 1280 nmol/l) of K5. Each treatment was set up in duplicate and incubated at 37°C for 72 h under normoxia (21% O 2 ) or hypoxia (1% O 2 ). Cell viability was measured by MTT (Sigma, St. Louis, MO, USA) assay, according to the manufacturer's protocol. Data represented absorbance and were expressed as the percentages of respective controls. For quantitative analysis of apoptosis, cells were treated as above for 48 h and then harvested for Annexin and propidium iodide staining, using the AnnexinV-FITC Apoptosis Detection Kit (Bender Med Systems, Vienna, Austria) according to the manufacturer's instructions. Cells treated with 25 μmol/l colchicine were used as positive control and those treated with PBS as negative control. The cells were subsequently counted using flow cytometry.


extfieldsofstudy: Medicine Biology

text:
Identification of miRNA, lncRNA and mRNA-associated ceRNA networks and potential biomarker for MELAS with mitochondrial DNA A3243G mutation

Researchers in the field of mitochondrial biology are increasingly unveiling of the complex mechanisms between mitochondrial dysfunction and noncoding RNAs (ncRNAs). However, roles of ncRNAs underlying mitochondrial myopathy remain unexplored. The aim of this study was to elucidate the regulating networks of dysregulated ncRNAs in Mitochondrial myopathy, Encephalopathy, Lactic Acidosis, and Stroke-like episodes (MELAS) with mitochondrial DNA (mtDNA) A3243G mutation, which might make contributions to the unveiling of the complex mechanisms underlying mitochondrial myopathy and, possibly, new tools applicable to clinical practice. Through high-throughput technology followed by quantitative real-time polymerase chain reaction (qRT-PCR) and bioinformatics analyses, for the first time, we found that the dysregulated muscle miRNAs and lncRNAs between 20 MELAS patients with mtDNA A3243G mutation and 20 controls formed complex regulation networks and participated in immune system, signal transduction, translation, muscle contraction and other pathways in discovery and training phase. Then, selected ncRNAs were validated in muscle and serum in independent validation cohorts by qRT-PCR. Finally, ROC curve analysis indicated reduced serum miR-27b-3p had the better diagnosis value than lactate and might serve as a novel, noninvasive biomarker for MELAS. Follow-up investigation is warranted to better understand roles of ncRNAs in mitochondrial myopathy pathogenesis.

understanding of mitochondrial epigenetics to mitochondrial dysfunction promises to deepen our understanding of mitochondrial disease pathogenesis, nominate new biomarkers, and motivate new therapeutic strategies.
The microRNAs and lncRNAs are non-coding RNAs that have been implicated as fine-tuning regulators controlling diverse biological processes at the level of post-transcriptional repression of their target genes in a sequence-specific manner 8 . At present, ncRNAs have become an integral component in understanding of the post-transcriptional mechanisms of disease. Nearly 2000 different miRNAs and more than 50,000 lncRNAs have been identified within the human genome with a constant increase in their number 9,10 . However, only a very limited number of nRNAs have biological function annotations. There have been accumulating evidences at remendous increase describing the involvement of ncRNAs regulating various biological and pathological processes in the past 10 years 11 . Researchers in the field of mitochondrial biology are increasingly incorporating these concepts into their studies and discovering a "gold mine" of regulatory interactions between mitochondrial function and ncRNAs biology 12 . An emerging body of research points to close interactions between miRNA and mitochondria 13 . Several studies have documented the localization of various miRNAs and their processing machinery in mitochondria, and that uncoupling of mitochondria decreased cellular miRNA functional efficiency 13,14 . It is noteworthy that miRNAs are the only class of relatively well-investigated ncRNAs in mitochondrial disorders until now, with previous studies profiling their dysregulation involved in metabolism and redox signaling, regulating mitochondrial gene expression, and association with occurrence and progression of varieties of mitochondrial disorders 14,15 . LncRNAs have also been implicated in mitochondrial function, such as previously mitochondrial-associated RMRP 16 , MDRL 17 , TGFB2-OT1 18 , CARL 19 , and UCA1 20 . But so far, the significance of ncRNAs in mitochondrial myopathy is sparsely characterized. Mitochondria-associated miRNAs, lncRNAs and their regulatory networks need to be specifically and thoroughly investigated in mitochondrial disease. In the present study we aimed to arouse the attention of ncRNAs regulation in studying mitochondrial myopathy, which might make contributions to the unveiling of the complex mechanisms underlying mitochondrial myopathy and, possibly, new tools applicable to clinical practice.
Through high-throughput technology, qRT-PCR and bioinformatics analyses, we performed, for the first time, one systematic and integrative analysis of the potential ncRNA-mRNA regulatory networks in muscle biopsies between MELAS patients with mtDNA A3243G mutation and controls. We found that the deregulated ncRNAs caused by mtDNA A3243G mutation formed complex regulation networks and participated in immune system, signal transduction, translation, muscle contraction and other pathways. Among these dysregulated ncRNAs, we found miR-27b-3p reduced in muscle and serum of MELAS patients with A3243G mutation or non-A3243G mutation, which was negatively associated with lactate, NMADS and the mutation load in muscle, had the better diagnosis value for MELAS than lactate. These results will enrich the genome-wide analysis of RNA molecules with potential cross-talk, provide the rationale for further confirmation studies in larger prospective cohorts and contribute to further systematic studies of mitochondrial disorders.

Results
Study population and design. The study population comprised a total of 74 MELAS patients with A3243G mutation, 22 MELAS patients with non-A3243G mutation and 69 controls without neuromuscular diseases. The demographic details are in Supplementary Table S1. This multistage, case-control study consisted of a preliminary step of high-throughput screening, secondary validation and quantification of candidate ncRNAs and mRNAs. In the primary screening, pooled muscle samples were from 20 MELAS patients with A3243G mutation and 20 normal controls for miRNA microarray analysis and lncRNA sequencing. The secondary screening of the most promising RNAs used the individual samples of the primary screening. The validation stage of selected ncRNAs was in a larger cohort of 54 MELAS patients with A3243G mutation and 49 age-and gender-matched normal controls. In addition, 34 MELAS with A3243G mutation from the validation cohort, 22 MELAS with non-A3243G mutation and 34 age-and sex-matched controls were employed to test the promising RNAs in serum. An overview of the study design is illustrated in Fig. 1.
Profiling the distinct miRNAs, lncRNAs and mRNAs caused by mtDNA A3243G mutation in muscle biopsies. Firstly, the genome wide miRNA, lncRNA and mRNA expression profiles in muscle biopsies were detected in discovery phase (including 20 MELAS patients with A3243G mutation and 20 matched controls). The criteria of corrected p-value < 0.05 and absolute fold change > 2 was used to identify the significantly dysregulated RNAs. The complete lists of dysregulated genes between two groups are shown in Supplementary Table S2. By microarry profiling, there were 31 and 41 miRNAs significantly up-and down-regulated in MELAS vs. controls, respectively ( Fig. 2A, Supplementary Table S2A). By sequencing profiling, 57 and 52 lncRNAs were identified significantly up-and down-regulated by A3243G mutation, respectively ( Fig. 2A, Supplementary Table S2B). Meanwhile, categorization of lncRNAs indicated that more antisense and intergenic than intronic and processed transcripts were dysregulated (Fig. 2B). Using the same sequencing profiling with lncRNAs, we found 167 and 162 mRNAs up-and down-regulated in muscle in MELAS patients vs. controls ( Fig. 2A, Supplementary Table S2C). In addition, mRNA profiles were compared to that of mitochondrial myopathy with TK2 mutation (GSE43698). When compared with normal subjects, differential gene expression profiling in muscle tissue was identified between TK2 and mtDNA A3243G mutation. Interestingly, 41 dysregulated mRNAs were common between two kinds of mitchondrial myopathy (Fig. 2C), while others seemed to be specific to the given phenotype. For example, compared patients with TK2 and mtDNA A3243G mutation to normal controls, 16 mRNAs, including EEF1A1, C3 and FMO2 were commonly upregulated; 25 mRNAs, including MYC, KLF9, PDE4B, ATP2A2, and HSPB6 were commonly down-regulated. This suggested that different kinds of mutation just caused partly common pathogenesis of mitochondrial myopathy. Technical validation of the dysregulated RNAs by qRT-PCR using a training set. The qRT-PCR was not only used to ensure the abnormal gene expression profiles detected by the initial high-throughout step not technical specific, but also to better define and select the meaningful ncRNAs for further analysis. MiRNAs/lncR-NAs that met the following criteria were selected for subsequent qRT-PCR validation. Firstly, meaningful ncRNAs should be dysregulated RNAs showing fold-change (log 2 ) > 1.5 or < − 1.5 (P < 0.01 and P ≠ 0) and having at least 30 copies/500 signal abundantly expressed in either group. Besides, the meaningful ncRNAs were selected from focus of the following network of miRNA-mRNA-lncRNA or lncRNA-miRNA-mRNA. Although seldomly selected, the mRNAs for subsequent validation should be potential targets of the selected miRNAs or lncRNAs. Corresponding to the above features, we screened out subsets of 6 miRNAs (including miR-6089, miR-27b-3p, miR-214-3p, miR-150-5p, let-7e-5p and miR-145-5p), 5 lncRNAs (including LINC01405, SNHG12, RP11-403P17.4, CTC-260E6.6 and RP11-357D18.1) and 6 mRNAs (including PDK4, CDKN1A, ATP2A2, SOD3, DDIT4 and EEF1A1) for qRT-PCR confirmation. Detailed characteristics of these RNAs are available in Supplementary Table S3. As shown in Fig. 2D, the findings were all consistent with results from our high-throughput data.
The miRNA, mRNA and lncRNA regulatory networks in MELAS with mtDNA A3243G mutation.
To explore biologic function of significantly dysregulated ncRNAs with fold change > 2 that were screened out from the discovery phase, firstly we obtained the potential target pairs between three kinds of RNAs by bioinformatics Figure 1. Flow-chat of the study design. First, differentially expressed RNAs in muscle were identified through high-throughout profiling. Then bioinformatics analysis was used to annotate the biological function and construct complex regulation networks between them. Furthermore, the selected meaningful RNAs were validated by qRT-PCR and then we searched their potential clinical implications in independent validation cohorts.
analyses. Because ncRNAs could induce repression or degradation of their targets, we just focused on target regulation with contrasting expression pattern. Based on that, three mtDNA A3243G-regulated interaction networks of miRNA-mRNA, lncRNA-mRNA and miRNA-lncRNA were constructed (Supplementary Figure S1 and Table S4). There were 61 miRNAs, 125 mRNAs and 259 potential pairs with contrasting expression within the miRNA-mRNA network. In regards to lncRNA-mRNA network, 973 potential lncRNA-mRNA pairs consisted of 34 lncRNAs and 90 mRNAs were constructed. In addition, the miRNA-lncRNA regulating network was formed by 77 miRNA-lncRNA pairs, which consisted of 7 up-regulated miRNAs with 6 down-regulated lncRNAs and 28 down-regulated miRNAs with 19 up-regulated lncRNAs. Finally, we combined all the three regulating networks together and constructed the whole target networks of miRNA-mRNA-lncRNA and lncRNA-miRNA-mRNA ( Fig. 3A   Gene ontology and KEGG pathway analyses. Functional annotation was performed by Gene Ontology (GO) and KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway analyses, to determine the biological significance of the whole dysregulated expressed mRNAs and two subsets of mRNAs which were potentially targeted by altered miRNAs/lncRNAs with negative expression between MELAS patients and controls. Classification and main enrichment items of GO and pathway analyses were indicated in Fig. 4 and Supplementary Table S5. In terms of the biological process, the common classification and significantly related categories of mRNAs, miRNAs, lncRNAs pointed to extracellular matrix organization, muscle contraction and protein binding (Fig. 4A, Supplementary Figure S2 and Table S5A). In addition, the significant GO items of the whole differential mRNAs still indicated the small molecule metabolic process, signal transduction and transcription (Fig. 4A,  Table S4. For the network, circular node represent miRNA, diamond node represent mRNA and triangle node represent lncRNA. Color variations of the nodes and lines represent the following: red node, up-regulated RNA; green node, down-regulated RNA; pink line, miRNA-mRNA network; blue line, miRNA-or mRNA-lncRNA network.
Scientific RepoRts | 7:41639 | DOI: 10.1038/srep41639 Supplementary Figure S2A and Table S5A). The significant GO classifications of miRNA-mRNA pairs were signal transduction, apoptotic process and development (Fig. 4A, Supplementary Figure S2B and Table S5A). The main associated GO classifications of lncRNA-mRNA pairs were small molecule metabolic process and muscle contraction (Fig. 4A, Supplementary Figure S2C and Table S5A). Correspondingly, all the three KEGG pathway analyses of mRNAs, miRNAs, lncRNAs indicated that the common significant items were immune system, infectious diseases, translation, signal transduction and cell communication. Specifically, KEGG pathway analyses of the whole differential mRNAs were associated with metabolic pathways, cytokine-cytokine receptor interaction, MAPK signaling pathway and so on (Fig. 4B, Supplementary Figure S2A and Table S5B). Subsequently, the pathway annotations of differential miRNAs were mainly involved in complement and coagulation cascades, RNA degradation and calcium signaling pathway (Fig. 4B, Supplementary Figure S2B and Table S5B). Pathway enrichment items of differential lncRNAs included metabolic pathways, muscle contraction, p53 signaling pathway and cytokine-cytokine receptor interaction (Fig. 4B, Supplementary Figure S2C and Table S5B).

Clinical validation of candidate miRNAs and lncRNAs using an independent MELAS patients and controls.
To test whether our findings in training set were replicable for the promising ncRNAs, an . Scatterplots of significant enrichment items of GO analyses (A) and pathway analyses (B) associated with abnormally expressed mRNAs, miRNA targets, lncRNA targets in muscle biopsies between MELAS patients and controls. The vertical axis is the GO or pathway term, and the horizontal axis is rich factor of GO or pathway enrichment. The P value denotes the significance of the GO or pathway item (the recommend cut-off of P value is 0.05).

Serum miR-27b-3p-based biomarker for MELAS diagnosis and association with clinical parameters.
To test whether reduced serum miR-27b-3p-based biomarker could distinguish MELAS patients from controls, ROC curve analyses were performed to compare the diagnostic value of serum miR-27b-3p with lactate. Serum lactate was up-regulated in MELAS patients vs. controls (Supplementary Figure S3). The AUC, adjusted p-value, and asymptotic 95% confidence interval of serum miR-27b-3p and lactate are shown in Supplementary Table S7. The AUC for lactate distinguishing the MELAS group from controls was 0.841 (95% CI: 0.762-0.919) (Fig. 6D). Serum miR-27b-3p showed higher diagnostic accuracy with AUC of 0.879 (95% CI: 0.806-0.952) than lactate (Fig. 6D). At the cut-off value of 0.93 for miR-27b-3p, the optimal sensitivity and specificity were 96.4% and 64.7%, respectively. Furthermore, correlation analyses were performed to examine the association of expression level of miR-27b-3p with clinical parameters, such as age, sex, BMI, NMDAS, disease duration, lactate and A3243G mutation load in muscle. Significant inverse correlation was found between miR-27b-3p and lactate (R = − 0.43, P < 0.001), NMDAS (R = − 0.644, P < 0.001) or A3243G mutation load in muscle (R = − 0.709, P = 0.002) (Fig. 6E,F and G). Serum miR-27b-3p appeared to decrease with the clinical severity of MELAS patients, and might serve as a useful tool to quantitatively monitor disease severity of MELAS patients.

Discussion
As the most commonly occurred pathogenic mutation locus of MELAS, mtDNA A3243G mutation in tRNALeu gene leads to pre-termination of transcription, compromising mitochondrial protein synthesis and ATP production, and activation of diverse retrograde signaling pathways 6 . Until recently, the clinical diagnosis and treament of patients with mitochondrial disorders remain a major challenge to clinicians due to the complex pathogenesis and the clinical variability 1 . In this study, we hypothesized that ncRNAs driven by mtDNA A3243G mutation could regulate the signaling pathways involved in mitochondrial myopathy and have firstly employed an integrated genomic approach to identify the dysregulated miRNAs and lncRNAs in muscle biopsies of a well characterized cohort of MELAS with mtDNA A3243G mutation (including 20 MELAS patients and 20 healthy controls in discovery and training phase, 54 MELAS patients and 49 controls in validation phase), which allowed us to identify ncRNAs of interest that may underlie MELAS pathogenesis for further studies.
The decisive factors involved in the secondary pathogenic cascades and compensatory adaptations in mitochondrial myopathy are still obscure. In current study, by coexpression profiling of mRNAs, miRNAs and lncRNAs, we identified the differential expressed RNAs in muscle biopsies of MELAS with mtDNA A3243G cytotoxicity. The patients were only restricted to MELAS with mtDNA A3243G mutation to reduce the potential for confounding effects due to sample heterogeneity, which would be expected with the inclusion of other ROC curve analyses of 6 miRNAs (A) and 5 lncRNAs (B) in muscle biopises for discriminating MELAS patients from controls. Serum and muscle miR-27b-3p were reduced in MELAS patients (C). It showed higher diagnostic accuracy than lactate for MELAS patients (D). Significant correlation was found between serum miR-27b-3p level and lactate (E), NMDAS (F) or A3243G mutation load in muscle (G). Expression level of serum miR-27b-3p was normalized to spiked-in cel-miR-39. Data are presented as the mean ± SD. *P < 0.05 versus control group.
Scientific RepoRts | 7:41639 | DOI: 10.1038/srep41639 kinds of mitochondrial mutation. Indeed the mitochondrial dysfunction associated with mtDNA A3243G cytotoxicity could activate a far greater retrograde signaling pathways than previously anticipated. We found that 329 mRNAs were significantly dysregulated in MELAS patients vs. controls. Interestingly, some of them were previously reported to be dysregulated in mitochondrial myopathy, including MYC, ATP2A2, ACTN3, AQP4, GLUL, TRIM63 and PFKFB3, which were partly consistent with transcriptome datasets from previous studies of mitochondrial genetic disorders [21][22][23][24] . This suggested that mitochondrial disease have disparate causes, and yet share common mechanistic themes 25 . In accordance with previous reports, further GO & Pathway analyses indicated these deregulated mRNAs were mainly involved in small molecule metabolic process, extracellular matrix organization, muscle contraction, protein binding, signal transduction, transcription and so on, indicating the importance of these signaling pathways in mitochondrial myopathy. As for the dysregulated ncRNAs, 72 miRNAs and 102 lncRNAs were found to be dysregulated in MELAS patients vs. controls. Interestingly, some miRNAs were also reported in gene expression profiling (GEP) datasets obtained from previously published oxidative stress cell models, such as miR-1-3p 26 , miR-15b-5p 27 , miR-29a-3p 28 , miR-128-3p 29 , miR-133a-3p 30 , miR-145-5p 31 , miR-150-5p 32 , miR-499a-5p 33 , which indicated that they might participate in the occurrence and development of MELAS through oxidative stress pathway. For example, some of these meaningful ncRNAs haven been investigated to be associated with mitochondrial function. Both miR-150-5p and let-7 family were thought to be mitochondria-enriched miRNAs 32 . The induction of miR-145-5p might promote cell apoptosis as its downstream target SOD2 was known to play the antioxidant defense 34 . However, miR-145-5p could suppresses ROS-induced Ca 2+ overload of cardiomyocytes by targeting CaMKIIδ 31 . In addition, PDE4B was found to be the potential target of miR-145-5p in our target analysis, and has been implicated in the suppression of apoptosis 35 . Over-expression of miR-145-5p significantly inhibited the H 2 O 2 -induced cellular apoptosis, ROS production, mitochondrial structure disruption as well as the activation of key signaling proteins in mitochondrial apoptotic pathway 36 . Previous studies have declared that miR-214-3p could regulate mitochondrial morphology and cell cycle by targeting Mitofusin2 37 , or Ca 2+ homeostasis by repression of Ncx1 38 . At present, more than 50,000 lncRNA genes have been identified. However, only a very limited number of lncRNAs have biological function annotations. Although some mitochondrial-related lncRNAs, including RMRP 16 , MDRL 17 , TGFB2-OT1 18 , CARL 19 and UCA1 20 , were already reported to regulate mitochondrial function, detailed participation mechanisms in mitochondrial myopathy of the dysregulated lncRNAs identified in this study have yet to be elucidated.
Integrated analyses of target prediction and functional enrichment based on coexpression networks may provide a clue for elucidating the precise underlying mechanisms of ncRNAs in pathogenesis caused by mtDNA A3243G mutation. We made interaction analyses and constructed three regulatory networks which consisted of miRNA-mRNA, lncRNA-mRNA or miRNA-lncRNA with contrast coexpresstion. In addition, the hypothesis of competing endogenous RNA (ceRNA) suggested large-scale regulatory networks between coding and noncoding RNAs across the transcriptome 39 . LncRNAs and mRNAs might interact through competition with the same miRNA, or miRNAs and lncRNAs could interact with each other, thus regulating function of the same mRNAs. So the mtDNA A3243G-regulated miRNA-mRNA-lncRNA and lncRNA-miRNA-mRNA networks in MELAS muscle were constructed. Because most of the ncRNAs in the coexpression network were not annotated yet, GO & Pathway analyses were performed to further annotate the biological function of these dysregulated ncRNAs. In accordance with analysis of the whole deregulated mRNAs, function analyses of potential targets of miRNA and lncRNA all showed that the mtDNA A3243G-regulated ncRNAs signaling pathways were enriched in immune response, signal transduction, apoptotic process and development, suggesting that they might contribute to MELAS pathogenesis through these signaling pathways. Finally, our study provided valuable clues for the identification of potential ceRNA pairs between miRNAs and mRNAs, lncRNAs and mRNAs or even miRNAs and lncRNAs caused by mtDNA A3243G mutation. It is very much worthy to perform further validation and functional examination in vivo and in vitro experiments to reveal the underlying mechanisms of these ncRNAs in MELAS.
Reduced serum miR-27b-3p showed potential as a biomarker for MELAS diagnosis. Although allowed for high-throughput analyses, the approach of the high-throughput assay with mixed samples in the first step of preliminary screening ignored the individual differences. So qRT-PCR validation was then followed to verify the meaningful ncRNAs selected from the regulatory network focus of MELAS. Consistent with results from the high-throughput data, we confirmed that muscle expressions of 4 ncRNAs (miR-6089, miR-27b-3p, miR-214-3p and LINC01405) were significantly down-regulated, and 7 ncRNAs (miR-150-5p, let-7e-5p, miR-145-5p, SNHG12, RP11-403P17.4, CTC-260E6.6 and RP11-357D18.1l) significantly up-regulated in MELAS muscle biopsies vs. controls, indicating that they might be acting as suppressors or accelerator in occurrence or progression of MELAS pathogenesis. Some of these meaningful ncRNAs haven been investigated to be associated with mitochondrial function. The oxidative stress-responsive miR-27b-3p has been shown to alleviate oxidative stress and inflammation 40 , regulate NO production 41 , improve mitochondrial function 42 and suppresses apoptosis 43,44 . On the other hand, miR-27b-3p might prevent the autophagic degradation of impaired mitochondria 45 . The dual effect of miR-27b-3p is possibly tissue-or cellular context-dependent, which is consistent with other miRNAs 46 . Apart from many targets including PINK1, FOXO1, PPARγ and Apaf-1 investigated before, MKNK2, GREM1 and EEF1A1 that might be potential targets of miR-27b-3p were revealed in the regulatory network of MELAS pathogenesis. Anyhow the role of miR-27b-3p in MELAS pathogenesis is still needed to elucidate. In addition, the following ROC curve analyses showed the highest sensitivity and specificity of muscle miR-27b-3p in differentiating MELAS patients from controls among the selected ncRNAs. So we evaluated its serum expression and confirmed its reduction not only in muscle and serum of MELAS with A3243G but also non-A3243G mutation. ROC curve analyses further indicated that serum miR-27b-3p had the better diagnostic value for MELAS that lactate. Correlation analyses also showed significant inverse correlation between miR-27b-3p and lactate or NMDAS. Muscle tissue might be an important source of serum miR-27b-3p because of its relevance with A3243G mutation load in muscle. Although elevated lactate level in serum and CSF may be helpful in diagnosis of mitochondrial diseases, it still lacks consistency as some patients with confirmed diagnosis don't show increased lactate level. The results in this study indicated that reduced serum miR-27b-3p could be a new potential biomarker applicable for MELAS diagnosis. However, it is still needed to expand the sample size in future studies.
Taken together, to our knowledge, this is the first study to construct coexpression networks among a series of miRNAs, mRNAs, and lncRNAs downstream of mtDNA A3243G mutation in mitochondrial myopathy. And we identified reduced serum miR-27b-3p for MELAS diagnosis with high sensitivity and specificity. However, sufficient sensitivity and specificity of the miR-27b-3p should be determined in the further well designed clinical studies. These findings may have substantial clinical significance or implications for illustrating the pathogenesis of mitochondrial myopathy, and follow-up investigation is warranted to better understand the detailed participation mechanisms of miRNA and lncRNA in mitochondrial myopathy pathogenesis.

Materials and Methods
Characteristics of individuals and muscle biopsy. The diagnosis of MELAS was based on the criteria of Newcastle Pediatric Mitochondrial Disease scale NPMDS and Sanger sequencing 47 . All muscle biopsy specimens were retrospectively collected from patients with mitochondrial myopathy and from subjects who were ultimately deemed to be free of neuromuscular diseases matched for age and gender selected as non-diseased controls. All procedures were approved by the Ethics committee of Qilu Hospital affiliated to Shandong University (Jinan, China), and carried out in accordance with the principles of Helsinki Declaration. All collecting samples of muscle biopsies were informed consent of patients at the Neuromuscular Center at Qilu Hospital.
MicroRNA microarray Assay. Total RNA (including miRNAs) from each frozen muscle was extracted using Animal Tissue Kit(LCS, TRK-1002) according to the manufacturer's instructions. The total RNA quantity and purity was analysed with Bioanalyzer 2100 and RNA 6000 Nano Lab Chip Kit (Agilent, CA, USA) with RIN number > 8.0. Qualified total RNA was further purified by RNeasy mini kit (QIAGEN, GmbH, Germany) and RNase-Free DNase Set (QIAGEN, GmbH, Germany). Total RNA from the same stage were mixed in equal amounts into two pooled samples: MELAS and Control. Microarray assay (miRHumanViruses-21b array) was performed by the service provider (LC Sciences). The assay started from 4 to 8 μ g total RNA sample which were 3′ -extended with a poly(A) tail using poly(A) polymerase. An oligonucleotide tag was then ligated to the poly(A) tail for later fluorescent dye staining. Hybridization was performed overnight on a μ Paraflo microfluidic chip using a micro-circulation pump (Atactic Technologies) 48 . On the microfluidic chip, each detection probe consisted of a chemically modified nucleotide coding segment complementary to target microRNA (from miRBase, http://www.mirbase.org/) or other RNA (control or customer defined sequences) and a spacer segment of polyethylene glycol to extend the coding segment away from the substrate. The detection probes were made by in situ synthesis using PGR (photogenerated reagent) chemistry. The hybridization melting temperatures were balanced by chemical modifications of the detection probes. Hybridization used 100 L 6xSSPE buffer (0.90 M NaCl, 60 mM Na2HPO4, 6 mM EDTA, pH 6.8) containing 25% formamide at 34 °C After RNA hybridization, tag-conjugating Cy3 dye were circulated through the microfluidic chip for dye staining. Fluorescence images were collected using a laser scanner (GenePix 4000B, Molecular Device) and digitized using Array-Pro image analysis software (Media Cybernetics). Data was analyzed by first subtracting the background and then normalizing the signals using a LOWESS filter (Locally-weighted Regression) 49 . To identify statistically significance for differentially expressed miRNAs, a criterion of signal > 500, |fold change| ≥ 1.0 and FDR < 0.01 was used. Microarray data was deposited in GEO database (GSE89059).
LncRNA library construction and Illumina sequencing. Total RNA was collected as described above. First, approximately 5 ug of total RNA was used to deplete ribosomal RNA according to the manufacturer's instructions of the Ribo-Zero ™ Magnetic Kit (Epicentre). Following purification, the RNA is fragmented into small pieces using divalent cations under elevated temperature. Then the cleaved RNA fragments were reverse-transcribed to create the final cDNA library in accordance with the protocol for the RNA-Seq sample preparation kit (Illumina, San Diego, USA). And then we performed the single-end sequencing (50 nt) on an Illumina Hiseq2500 sequencer at the LC Biotech (Hangzhou, China) following the vendor's recommended protocol. The sequence results were obtained as the FPKM (fragment per kilobase of exons per million reads) for each transcript and were deposited in GEO database (GSE89065).

Confirmation and quantification of miRNAs, lncRNAs and mRNAs by qRT-PCR.
To validate the microarray and sequencing datasets, qRT-PCR were performed for most promising differential miRNAs, lncR-NAs and mRNAs. For muscle miRNA and lncRNA/mRNA confirmation, U6 and GAPDH were used as the endogenous reference gene, respectively. Total muscle RNA for qRT-PCR was collected as described above. For serum miR-27b-3p, synthetic cel-miR-39 was added to each denatured sample (25 fmol in a 5 μ l total volume) as stable control for normalization in serum. The serum miRNA was isolated from 400 μ L serum sample for each individual using the Ambion mirVana ™ PARIS ™ Kit (Ambion). For muscle lncRNA/mRNA, cDNA was synthesized from total RNA using the M-MLV RTase cDNA Synthesis kit (Takara, Dalian, China) according to the manufacturer's instructions. For muscle miRNA, cDNA was reverse transcribed from total RNA in a final volume of 20 μ l using One Step PrimeScript miRNA cDNA Synthesis Kit (Takara, Dalian, China). Quantitative real-time PCR was carried out on an ABI Prism 7900HT (Applied Biosystems) with SYBR Premix Ex Taq TM II (Takara, Japan) according to the manufacturer's protocol. The primer sequences are listed in the Supplementary Table S7. All reactions were tested in triplicate, and differences of threshold cycles between target genes and house-keeping genes (U6 in muscle miRNA, GAPDH in muscle lncRNA/mRNA, and cel-miR-39 in serum miR-27b-3p) were calculated using the 2 −ΔΔCT method and analyzed by paired t tests to evaluate statistical significance relative to Scientific RepoRts | 7:41639 | DOI: 10.1038/srep41639 controls. Correct qPCR product size was verified by agarose gel electrophoresis, and melting curve measurements were made after every experiment.
The construction of miRNA, lncRNA and mRNA regulatory networks. Integration of ncRNA and mRNA expression profiling may allow the identification of functional links between dysregulated ncRNAs and their targets. MiRNA target genes were predicted by TargetScan, miRecords, TarBase, and miRecords software. The putative miRNA binding sites of the lncRNAs sequences were performed using the RegRNA software (http://regrna2.mbc.nctu.edu.tw/). LncRNAs were predicted via two independent algorithms: cis-or trans-regulatory effects 50 . We selected the dysregulated target RNA with an inverse correlation of expression with the respective differential ncRNA. The miRNA-mRNA-lncRNA and lncRNA-miRNA-mRNA interaction networks were built with the Cytoscape Software (version 3.1.0).
Functional annotation of differentially expressed miRNAs and lncRNAs. Ingenuity pathways analysis (IPA) software was used to identify signaling pathways regulated by mtDNA A3243G mutation mining from our miRNAs/lncRNAs datasets. Functional annotation was performed by Gene Ontology (GO) and Pathway analysis to determine the biological significance of the differentially expressed mRNAs targeted by miRNA/lncRNA. Dysregulated mRNAs were input into the Database for Annotation, Visualization and Integrated Discovery (DAVID; http://david.abcc.ncifcrf.gov/) v6.7, which utilized GO to identify the molecular function represented in the gene profile and the Kyoto Encyclopedia of Genes and Genomes (KEGG; http://www. genome.jp/kegg/) to analyze the potential functions of these genes in pathways. The GO analysis was divided into molecular function, biological process and cellular component. Finally, results of the enrichment analysis were ranked by the significantly enriched GO term and KEGG pathway according to hypergeometric distribution. A corrected p-value < 0.05 was the threshold for statistically significant correlation. As the enrichment increases, the corresponding function is more specific, which helps us to identify GOs and pathways with more concrete functional description in MELAS.

Statistical analysis.
Clinical characteristics were compared using χ 2 test of independence for qualitative variables, t-test of quantitative variables with normal distribution, the non-parametric Kruskall-Wallis test or the Mann-Whitney U test of quantitative variables with skewed distribution. The threshold value we used to screen differentially expressed lncRNAs/miRNAs/mRNAs was a fold change > 1.5(P < 0.05). Expression levels of lncRNAs/miRNAs were compared using the Kruskall-Wallis test or the Mann-Whitney U test. ROC curves and AUC were established to evaluate the diagnostic value of lncRNAs and miRNAs for differentiating melas patients with healthy controls. In ROC analysis, the normalized expression level of RNAs ( 2 −ΔΔCT ) was selected as the test variable for the selected lncRNAs and miRNAs. All data were expressed as mean ± standard deviation (mean ± SD). A p value of less than 0.05 was considered statistically significant. All analyses were performed by SPSS 17.0 software (SPSS, Chicago, IL, USA), R (http://www.r-project.org/, version 2.15.3) and Bioconductor 30, or Graphpad Prism (version 5.0; Graphpad software).


text:
Molecular Docking of Phytoligands to the viral protein receptor P. monodon Rab7

The development of shrimp aquaculture has been severely affected by viral diseases resulting in a huge economic burden to the industry. White spot disease (WSD) has caused severe mortality in farmed shrimp in many countries. Globally aquaculture industries face huge economic losses due to rapid spread of White Spot Syndrome Virus (WSSV) disease that can cause 100% mortality in a short period of 3-10 days. In the present study in order to prevent the spread of WSSV disease in shrimps, the receptor, PmRab7 has been chosen as the drug target. Due to the absence of a precise 3D structure of the target, homology-modeling approach was employed to obtain the structure that was validated later. This structure was then used as a template to screen selective phytomolecules as potential antiviral agents and their docking results with the target are analyzed based on their energy scores. Identification of the drug-like molecule obtained from the docking analysis would be used to optimize to a candidate drug. This is expected to play a role of the inhibitor that blocks the binding of the viral protein to the receptor, duly preventing the WSSV disease.

Background:
Viral diseases strongly affect the shrimp aquaculture industry resulting in a heavy loss to the economy, of all the known viruses in the field of aquaculture, such as, Yellow head virus (YHV), Infectious myonecrosis virus (IMNV), Taura syndrome virus (TSV), White spot syndrome virus (WSSV), and Infectious hypodermal and hematopoietic necrosis virus (IHHNV). WSSV exists to be the most devastating shrimp pathogen in cultured shrimp, globally [1-3]. The White spot disease (WSD), resulting from the infection of WSSV has caused severe mortality in farmed shrimps in many countries [4][5][6]. There is no treatment currently available to control the spread of the disease due to lack of effective therapeutics.
The WSSV genome has been completely sequenced for four isolates, and more than 50 virus-encoded proteins have been identified as structural proteins [7]. Among these structural proteins, the envelope proteins are extremely important because they are believed to be the first molecules to interact with the host and, consequently, play critical roles in cell targeting as well as in triggering host defences. The precise mechanism of entry of WSSV into P. monodon is yet unknown and it is speculated that some receptor proteins viz. viral attachment protein,PmCPB, βintegrin and PmRab7 are involved in directing WSSV into the host. The viral envelope and its structural proteins form the first and most important component of virus to directly come into contact with the shrimp [8,9].
The virus entry into the host cells occurs through the interactions of viral envelope proteins with the host cellular receptors. These interactions are so specific that inhibitors that could prevent these interactions show promising anti-viral activity. In case of WSSV, VP28 is one of the major envelope proteins that are identified as the interacting partner with the host surface cellular protein P. monodon Rab7 [8,10]. We have been working on an aim at studying the structural details on the residue specific interactions between the WSSV VP28 -Shrimp PmRab7 proteins by employing a wide range of biophysical, especially NMR, and biochemical methods. The outcome of this study identifies the binding interface and specific residues involved in the viral-host protein-protein interactions [11]. This will help us to design new drug like molecules mimicking the envelope protein binding interface or the receptor protein binding groove, such that the molecules that we design from the knowledge of the atomic structure of VP28 / PmRab7 prevent and block the viral entry onto the shrimp host, thus controlling the WSSV infection. The scheme representing the design of inhibitors for prevention of WSSV entry into the shrimps is shown in Figure 1. Further, the high-resolution structural characterization of the proteinprotein interfaces deliver us the essential knowledge to design inhibitors with high accuracy.
In line with the above hypothesis, here in the present study, in order to prevent the spread of WSSV disease in shrimps, the receptor, PmRab7 has been chosen as the drug target. The viral protein VP28 could also be considered as another target, as the reverse approach through which the viral infection can be combated too. The use of the designed / synthesized molecule in the experiments to be done in ponds would be difficult and is also not an eco friendly approach. Hence, use of recombinant pmRab7 is considered to be the ideal choice. But the three dimensional atomic structure of the target protein PmRab7 is unavailable in Protein Data Bank. Hence, its 3D model has been constructed using homology-modeling approach.
Plant-derived natural products play a significant role by being a lead molecule in the development of drug candidates. Herbal extracts represent the primary form of health care for a major proportion of the world population and are the important sources of single-molecule drug leads [12]. A prominent example is the anti-malarial activity of Artimisiaannua discovered by Professor Tu, of China Academy of Chinese Traditional Medicine, recipient of the 2015 Nobel Prize for Physiology and Medicine [13,14]. Phytomolecules like Baicalein, Luteolin, Quercetin and Kaempferol are potential antiviral agents against a wide range of important viruses including Dengue, HIV, H5N1 influenza A viruses, Coxsackie virus, CHIKV and Japanese encephalitis virus [15]. The quest for identification of the leads for various drug targets led us to investigate the phytomolecules as antiWSSV leads and the above flavonoids were considered for our present study.
The contribution by homology modeling and computational biology has accelerated the pace of drug discovery. It is used in the biopharmaceutical industry to discover and develop new lead compounds. By this route, one can visualize the possibilities of binding of potential small molecules as ligands / inhibitors. The analyses of the different docked conformations are carried out using the scores / energies based on their binding affinities, as parameters to evaluate the ideal ligand. Therefore, we have designed the present study using computational approaches to discover the potential of these flavonoids targeting PmRab7 and preventing the complex formation as depicted in Figure 1. Though Figure 1 is displayed with two ligands, each one encompassing the interface region of VP28 and Rab7 respectively, the idea is to deal with them independently. In the present study, we concentrate on the surface of Rab7 interacting with VP28. Essentially, we have to design a small molecule that is mimicking this interface so that it can act as a good inhibitor.

Methodology: Protein Preparation and Homology Modeling:
In the isolation and characterization of white spot syndrome virus (WSSV)-binding proteins from shrimp, out of the three membrane-associated molecules identified, a 25-kDa protein was identified that had a primary structure with high homology to the small GTP-binding protein Rab7, named it Penaeus monodon Rab7 (PmRab7). This is the first Rab homologue from crustaceans to be identified and characterized. It is shown that the GTP-binding protein Rab7 may be a receptor for VP28 envelope protein of WSSV in shrimp. Rab proteins are known to be regulators of vesicle budding and fusion events and represent a family of over 30 proteins that are localized on the surfaces of distinct membrane-enclosed compartments of exocytic and endocytic pathways. They are found in all eukaryotes, including yeasts, plants, insects, and mammals.
PmRab7 has the four conserved GTP-binding regions of the small G protein superfamily (G1 and G3 to G5), as well as an effector site (G2). These five regions are characteristic of Rab proteins. The sequence analysis suggests that PmRab7 may be an active GTPase that could cycle between the GDP-and GTP-bound states.
The 3D structure of the PmRab7 protein is unavailable in protein data bank. The protein sequence of Penaeusmonodon Rab7 (Accession No.ABB70064.1) was taken from National Centre for Biotechnology Information (NCBI). The protein sequence was used to develop homology-modeled structure. The threedimensional coordinates of PmRab7 were constructed using Schrödinger [16], a comparative protein-modeling program. It computes a model based on the alignment of the sequence to be modeled with known related 3D structures. The homology modelled 3D structure of PmRab7 in its complex form with the ligand GDP has been obtained from Schrodinger and the model has structural similarity with the rat Rab7. Water molecules, ligands and other heteroatoms were removed from the model. Hydrogen atoms were added to the model using CHARMm force field. Energy minimization was performed by using conjugate gradient method with an RMS gradient of 0.01kcal/Å mol on Accelyrs-Discovery studio client [17].

Model Validation and Structural Characterization:
Structure verification programs such as PROCHECK, Ramachandran Plot and VERIFY3D were used from the Structure Analysis and Verification Server (http://nihserver.mbi.ucla.edu/SAVES/) to evaluate the validity of the homology modeled PmRab7. These programs validate the homology-modeled structure by checking the quality of the modeled structure, its resolution and refinement. PROCHECK, the structure verification program checks the stereo-chemical quality of the protein structure and determines the quality of the predicted structure by assessing parameters such as bond lengths and angles, planarity of the peptide bonds, geometry of the hydrogen bonds, and side chain conformations of protein structures as a function of atomic resolution. Ramachandran Plot visualizes energetically allowed regions for backbone dihedral angles ψ against φ of amino acid residues in protein structure. The Verify3D assess protein structures using three-dimensional profiles. It analyzes the compatibility of an atomic model (3D) with its own amino acid sequence (1D) by assigning a structural class based on its location and environment (alpha, beta, loop, polar, apolar etc). Then a score for each of the 20 amino acids in this structural class is obtained from a database generated from good structures. The scores range from -1 (bad score) to +1 (good score).

Ligand Preparation:
Among the plant-derived flavonoids, the potential phytochemical ligands namely Kaempferol ,Luteolin, Baicalein, and Quercetin were chosen based on the literature survey as potent compounds possessing antiviral activity. The molecules were docked against PmRab7. Their basic structure is a skeleton of diphenylpropane, namely, two benzene rings (two extreme rings) linked by a threecarbon chain that forms a closed pyran ring (heterocyclic ring containing oxygen, the middle ring) with benzenic first ring. These phytochemicals were retrieved from the pubchem database and the chemical structures were generated using SMILES (Simplified Molecular Input Line Entry Specification) by using the Discovery Studio 2.4 version. Also considered in this analysis is GDP (guanosinediphosphate), the ligand that is being complexed in the crystal structure of rat Rab7. GDP consists of three moieties: The guanine ring, which is a nitrogen-containing heterocycle; the ribose sugar ring; and the diphosphate group. The structural details namely, the chemical name, PubChem ID and Molecular formula of the selected phytochemicals along with GDP are given in Table 1.

Drug likeliness prediction:
Drug like properties of the ligands were predicted by using the Discovery Studio 2.4 version. The Lipinski's rule helps in distinguishing drug-like and nondrug-like properties and predicts high probability of success or failure due to drug likeliness of the molecules. The Lipsinki's filter helps in early preclinical assessment, thereby avoiding costly late stage preclinical and clinical failures. Lipinski's rules state that ideal drug molecules possess 5 hydrogen bond donors, not more than 10 hydrogen bond acceptors, Molecular weight not more than 500 and LogP not more than 5 [18].

Molecular Docking:
The grid-based molecular docking method is used here using the program Cdocker (Accelrys) that employs CHARMmforce field. The targetis held rigid while the ligands are allowed to be flexible during the refinement. Since the ligand, GDP was already complexed in the crystal structure; the same binding site info is considered. Hence, it is possible, however, to specify the ligand placement in the active site using a binding site sphere. For this purpose, the GDP ligand present in the active site of the target protein was used to generate the sphere around the active site. Then the prepared ligands are docked to the active site using default parameters. The results of the docking enabled the ranking of the docked conformation of the ligands according to their Cdocker energy values.

Analyses of the ligand binding sites:
The docking poses were ranked according to their docking energies. The scoring function in Cdockerwas used to predict the binding affinity of one ligand to the target molecule. In addition to the structural information, each record includes the Cdocker energy reported as negative value, where the higher value indicates a more favorable binding. This enables the energy to be used like a score. This score includes internal ligand strain energy and receptor-ligand interaction energy, and is used to sort the poses of each input ligand. The molecular visualizations of the docked complexes were analyzedusing the Discovery Studio 2.4 version [17].

Results:
In this study, theshrimp receptor protein, Pm Rab7, was considered as the target protein towards the lead identification. There are no 3D atomic structures available for the target. Hence, the aim of this study consisted of two parts. First is to build the homology model of the receptor PmRab7 and then to validate and characterize the constructed model. Secondly, use the constructed model as the target to dock with the selected Phytoligands. The homology model was constructed and the 3D structure of shrimp receptor protein PmRab7 was designed by using molecular modeling software Schrödinger. The coordinates of the resultant model were used as the template for the target. The ribbon diagram of the homology modeled PmRab7 is as shown in Figure 2.
As the model has been derived from the structural homology principles, based mainly on Rat Rab7of RattusNorvegicus, the resultant modeled protein has been found to be complexed with the ligand GDP and has been considered for the docking study. Rab7 is a small GTPase and act like a switch which is turned on and off by GTP and GDP molecules. The active site that has been occupied by the ligand GDP was used for the docking of the phytochemical ligands. Four phytochemical ligands belonging to flavonoids class of compounds were considered as the ligands for this computational study. The drug likeliness of the selected ligands along with that of GDP has been presented in Table 2.
The table provides the details of the molecular weight, number of donor and acceptor atoms and AlogP values. The data would reveal the likelihood of the ligands that could form a potential lead for further analysis.  Table 3 lists the results of the docking analysis by providing the Cdocker energy for all the selected ligands with PmRab7 as the target. Higher the negative energy corresponds to the stable binding with the target. The Cdocker energy for GDP was also considered and it is proved to be a better complex with PmRab7  Table 4. Further, the intermolecular interactions between these ligands with PmRab7 have been shown as a 2D representation in Figure 3 (a) for GDP, (b) for Quercetin.

Discussion:
With the non-availability of the atomic structure for the target Pm Rab 7, a homology model was constructed. This model resembles the rat Rab7 with the exclusion of about 40 residues at the Cterminus. The model consists of five helices, a double, and a triple stranded beta sheets (Figure2).
It is globular in conformation with the turns protruding outside the globularity. The structure Analysis and Verification Server have validated the success of the constructed model, which consist of the individual programs such as PROCHECK, Ramachandran Plot and VERIFY3D. According to PROCHECK, that is based on the stereochemical quality, the structure is refined to a 1.5 A resolution. The Ramachandran plot and the plot statistics confirm the good quality of the structure. Over all, the model has good quality with more than 96.8 % of residues under the most favoured and additional allowed regions in the Ramachandran Map. Only two residues are in the disallowed region. Hence, the model is considered as the best with the dihedral angles of the maximum number of residues in the core region and a minimum most in the disallowed region of the Ramachandran Map. From Verify3D, it is observed that nearly 92% of the residues from the homology model of PmRab7 are compatible with its own amino acid sequence, indicating the high quality of the model.
The results of this docking study with the good quality homology modeled Pm Rab7 were analyzed through four different routes, namely, structural visualization and comparison with GDP, the drug like property analysis, analysis of Cdocker energy and finally, the hydrogen bonding and other ligand interactions with PmRab7. Simple visualization of the ligand structures given in Table 1 indicates clearly that all the ligands belong to a single cluster and have structural similarity among themselves havingbenzo-γpyrone derivatives with phenolic and pyrane rings. The variation in the molecular weights for the ligands is by a mere oxygen and one could expect a consistency of the derivatives that could bind with PmRab7. .
The analyses of the Drug Likeliness of the selected phytochemicals indicate that all the ligands satisfy completely the Lipinski's rules (Table 2) and found as potential candidates to be the leads [16]. Surprisingly, the ligand GDP doesn't satisfy the Lipinski's rules of acceptor and donor atoms. Though GDP fits nicely into the catalytic site of PmRab7, it could not be a proper candidate towards the making of a drug. Thus, all the ligands chosen for this study are suitable to be the leads that have the potential to form the drug. Though all four phytoligands were similar in structure visualization, drug likeliness and high negative value of Cdocker energy, we have considered Quercetin only for analysing the intermolecular interactions between Quercetin and PmRab7, due to the reason that Quercetin is found to be the best among all in all respects. In support of the previous analyses and results (Figure 3 and Table 4), Quercetin interacted with PmRab7 by making four H-bonds among them. These H-bond interactions with Gly18, Gly20, Gln67 & Lys126 were also present between GDP and PmRab7. We were interested to find out the common motifs / residues in PmRab7 that interact with all the ligands. GDP is found to interact with 9 residues (Gly18, Val19, Gly20, Lys 21, Thr 22, Ser 23, Gln 36, Gln67& Lys126) of PmRab7.
Quercetin has all four residues in common in GDP. The superiority of Quercetin is observed from the hydrogen bonds formed with the residues of PmRab7 to its core ring structures, whereas in GDP, the hydrogen bonds are only to the extended diphosphate region. More hydrogen bonding cause the lowest Cdocker energy in GDP, but with no interaction with the core ring group, it may not be the proper lead molecule. Thus, this analysis strongly supports the earlier observation that Quercetin is best among the four and might act as the potential inhibitor against WSSV. This would be verified and validated by conducting "wet lab" experiments and can be taken up further for confirming as a potential drug candidate.   In order to validate our claim, we are working on the following, as part of our wet lab studies based on our new proposal on the protein-protein interactions between VP28 and PmRab7. We are independently working on the expression and purification of these proteins by recombinant methods. These proteins would be labelled with N15 and C13 to have samples for Multidimensional NMR studies.
We also plan to demonstrate the protein crystallization of Rab7 for the X-ray analysis and then use this model for searching potential inhibitors. Protein complex crystallization experiments would also be attempted by our structural biology collaborators to determine the VP28 and PmRab7 complex structure. Since the model showed that Quercetin was very likely to interact with PmRab7, we plan to validate the interaction easily using the following techniques, isothermal titration calorimetry and surface plasmon resonance.
We also plan to test if the interaction between VP28 and PmRab7 is disrupted by adding of Quercetin in vitro, following the GST pull-down assay conditions provided by Sritunyalucksana et al.2006 [9]. We could also conduct some WSSV challenge tests with potential inhibitors, to test whether Quercetin can reduce the mortality caused by WSSV.

Conclusion:
The White Spot Disease occurring due to WSSV is a major threat to shrimp farming as it can cause complete mortality of the infected shrimps within few days. There is an urgent need to find effective inexpensive anti-WSSV drugs from natural sources to prevent the infection.
The target for this structural and computational host pathogen interaction studies is considered as the shrimp receptor Pm Rab7. But, the major hurdle in having an ideally usable target for anti WSSV drug is the 3D structure of PmRab7.Hence, homology modeling was used to build a 3D structure of PmRab7 and it was successfully validated by the standard procedures.
Plants have derived effective antiviral drugs for some of the viruses of both plant and animal origin. Working along this direction, the present research work is based on the screening of the Phytochemicals that could act as the leads towards the development of the drugs. Hence docking studies were performed that depicted more affinity of binding of Quercetin with PmRab7. Out of the Phytochemicals that were screened, Quercetin is predicted to be the suitable lead to proceed further. Hence use of this compound during maintenance of shrimp population in ponds, tanks or natural reservoirs might decrease the interactions between PmRab7 and the viral protein.

Reduction in the interactions
between the shrimp receptor protein PmRab7 and Viral envelop protein VP28 can further block one of the major route of entry of the virus into the shrimps. Further wet lab experimental validation as mentioned above can pave way for the generation of potential drug that can inhibit the entry of WSSV into the shrimp by this route, which in turn, can reduce the WSSV occurrence in the shrimp Penaeus monodon.


text:
ATP serves an anti-inflammatory role by enhancing β-defensin-2 response in acute pneumonia of rat

The aim of the current study was to evaluate the effect of ATP on the expression of rat β-defensin-2 (rBD-2) in a time-dependent manner, as well as its therapeutic value in an acute pneumonia rat model. A total of 30 rats as a treatment group and 30 as a control group were treated with the same dose of ATP and normal saline, respectively, lung tissues were isolated from rat and expression of rBD-2 mRNA was assessed with reverse transcription-quantitative polymerase chain reaction (RT-qPCR) at 12, 24 and 36 h following treatment. Rats were divided into five groups: The control group treated with normal saline, the Pseudomonas aeruginosa (PA) infected group, group treated with ATP, group treated with cephalosporins, and the group treated with both ATP and cephalosporins. At 24 h following treatment, rat serum and lung tissues were collected for assessment of histological changes, and alterations to expression of the rBD-2 protein by immunohistochemistry, expression of tumor necrosis factor (TNF)-α and interleukin (IL)-6 proteins by ELISA. RT-qPCR results indicated that the expression of rBD-2 mRNA was upregulated in response to ATP stimulation in lung tissues of rat, reaching its highest peak at 24 h. Immunohistochemistry demonstrated that ATP treatment enhanced the expression of rBD-2 protein in rat lungs. Ceftazidime and ATP protected lungs from infection of PA and reduced the pathological damage of the lung. Overexpression of rBD-2 by ATP led to decreased protein expression of TNF-α and IL-6 in lung tissues and serum. ATP upregulates the expression of rBD-2 and serves an anti-inflammatory role in the acute pneumonia of a rat model.

ATP is a well-known energy source inside cells, it also functions as an essential messenger between cells. Activation of P2 receptors, P2X and P2Y, by ATP may initiate the release of sodium and calcium ions, producing different cellular effects (18). ATP has a wide range of functions in the human body, including the development of pain, the heart and blood vessel disorders, secretion of digestive enzymes, as well as tumor-killing effects (18). In addition, ATP serves an important role in inflammation and infection. ATP activation of the P2X7 receptor results in the processing and secretion of mature interleukin (IL)-1β, which is a potent cytokine capable of inducing inflammation (19,20). Moreover, active ATP release and adenosine receptor signaling are required for activation of the NLR pyrin domain containing 3 (NLRP3) inflammasome that is a major component of the innate immune system (21). Nonetheless, little is known regarding the relationship between BD2 and ATP, as well as their roles in the acute pneumonia of rats.
In the present study, the authors evaluated the impact of ATP on the expression of rat BD-2 (rBD-2) in a time-dependent manner, as well as their protective effect in the acute pneumonia of rat caused by Pseudomonas aeruginosa (PA) (Xinjiang Uygur Autonomous Region Food and Drug Administration, Urumqi, China). The completion of the study provides some novel findings on rBD-2 in acute pneumonia which will be useful to guide further investigation and targeted therapy for this disease.

Materials and methods
Animal housing and treatment. Given the homology between rBD-2 and hBD-2, such as 6 cysteine and other lysine and arginine residues (22), the rat was chosen as a model to study the function of hBD-2 in vivo. BALB/c rats (4-6 weeks old, male/female ratio, 1:1; 20-30 g), provided by Medical School of Tongji University (Shanghai, China), were housed under controlled conditions. Procedures involving animals and their care were approved in accordance with the guidelines of the Experimental Animal Ethics Committee of Medical School of Tongji University (Shanghai, China; no. 2010-1-1). A total of 60 rats were randomly divided into treatment and control groups with 30 rats in each. A control group was administrated with normal saline at 250 µl/100 g body weight, treatment group was injected with ATP at 0.12 mg/10 g body weight through the tail vein. Rats were sacrificed at 12, 24 and 36 h, respectively after treatment, and the lung tissues were isolated from each rat and stored in liquid nitrogen.
Establishment of acute pneumonia model of rat, animal grouping and treatment. Rats were randomly divided into four treatment groups and one control group, with eight rats each. Rats were anesthetized with ether, 1 cm skin incisions were made along the middle of the mouse's neck, and muscle tissues were bluntly dissected exposing the trachea and thyroid. A total of 0.2 ml PA (6x10 8 CFU/ml, Xinjiang Uygur Autonomous Region Food and Drug Administration) was injected into each rat trachea in treatment groups, while, the control group was injected with the same dose of sterile phosphate buffered saline, then incisions were sutured. After three days, treatment groups were divided into four different groups: PA infected group, group treated with ATP at 16 mg/kg or cephalosporins at 0.4 g/kg through the tail vein, group treated with both ATP at 16 mg/kg and cephalosporins at 0.4 g/kg through tail vein. Rats were sacrificed 24 h following treatment, lung tissues and serum were collected and stored in liquid nitrogen for subsequent experiments.
Pathological examination. The grouping and treatment were the same as described above, lung tissues were isolated from rats 24 h following ATP and cephalosporin intervention. Fresh lung specimens were fixed in 4% paraformaldehyde, embedded in paraffin, cut and heated, and subjected to hematoxylin and eosin (H&E) staining. The morphological changes of H&E-stained tissues were analyzed by Leica digital microscopy.
Immunohistochemistry. Formalin-fixed, paraffin-embedded tissue sections of lung tissues were obtained from treatment and control groups as described above in animal grouping and treatment. The expression of rBD-2 was detected following the immunohistochemistry protocol of Sigma-Aldrich; Merck KGaA (Darmstadt, Germany). The sections were incubated with primary rabbit polyclonal antibody rBD-2 (YY3068R; Santa Cruz Biotechnology, Inc., Dallas, TX, USA) at 1:200 dilution overnight at 4˚C and secondary antibody [SP9000D; DAB Detection kit (Streptavidin-biotin); ZSGB-BIO, Beijing, China] for another 2 h at room temperature. Finally, visualization of the sections was conducted and images were taken with a Leica confocal laser scanning microscope (Leica Microsystems, Wetzlar, Germany).

Determination of IL-6 and tumor necrosis factor (TNF)-α expression in lung tissues and peripheral blood.
At 24 h following ATP and cephalosporin intervention, rats from control and treatment groups were sacrificed and lung tissues were isolated, centrifuged at 2,045 x g for 20 min at room temperature followed by collection of supernatants.
The peripheral blood samples were also collected. IL-6 and TNF-α expression was evaluated using ELISA assay kits in accordance with the instructions (Invitrogen; Thermo Fisher Scientific, Inc.).
Statistical analyses. Data were presented as mean ± standard deviation. Student's t-test or SNK-q test was used to determine the statistical significance between each two groups, differences among multiple groups were analyzed using an analysis of variance (ANOVA). The SNK-q test was used for multiple comparisons. P<0.05 was considered to indicate a statistically significant difference.

ATP enhanced expression of rBD-2 mRNA in rat lungs.
Expression of rBD-2 in lung tissues of rats in response to ATP has not yet been investigated. Therefore, the authors determined induction of rBD-2 expression by ATP in the lung tissues of rat. RT-qPCR was conducted to analyze the expression pattern of rBD-2 following ATP stimulation.  Fig. 1). It increased greatly at 12 h, reached highest peak at 24 h and gradually declined by ~20% at 36 h following ATP treatment, however, significant differences of the expression of rBD-2 were observed between each control and treatment group at 12, 24 and 36 h. Therefore, the authors confirmed that ATP could upregulate rBD-2 mRNA expression in rat lungs.
ATP increases expression of rBD-2 protein in rat lungs. As expression of rBD-2 mRNA was upregulated by ATP in rat lungs, the next question concerned whether ATP promotes the production of rBD-2 protein in vivo. In the present study, the authors reported that the expression of the rBD-2 protein was slightly increased in the PA infected group, an enhanced expression of the rBD-2 protein was observed in the ATP or cephalosporin treatment group and the cephalosporin and ATP treatment group, in comparison with that in PA infected group and control group, suggesting ATP could induce rBD-2 expression at both mRNA and protein levels in rat lungs (Fig. 2).
ATP alleviates histological changes in rat lungs. Microcopy indicated that alveolar and interstitial lung structures were normal, alveolar inflatable was good, no exudation of red blood cells was observed, part of bronchial and alveolar walls were infiltrated with a small amount of inflammatory cells in the lung tissues of the control group. The PA infected group presented alveolar wall edema, infiltration of a large number of inflammatory cells (neutrophils and lymphocytes) on bronchial and alveolar walls, red blood cells and cellulose exudation in the alveolar space. ATP or ceftazidime intervention groups presented an improved lung edema, less inflammatory cell infiltration, thicker partial alveolar wall and a small amount of red blood cells as compared with the PA infected group. The results demonstrated that ceftazidime, ATP and ceftazidime + ATP treatment groups improved the histological destruction of lung tissues caused by PA (Fig. 3).
ATP decreases expression of TNF-α and IL-6 levels in rat lungs. Following this, the authors next investigated mechanisms by which enhancement of rBD-2 by ATP presented a beneficial effect on the PA infection in rat lungs. In order to clarify whether rBD-2 regulated the secretion of proinflammatory cytokines, IL-6 and TNF-α expression was analyzed by ELISA. Significant differences of IL-6 and TNF-α expression were identified in both lung and serum among five groups . Moreover, IL-6 and TNF-α expression in the rat lung and serum was reduced in ceftazidime, ATP and ceftazidime + ATP treatment groups relative to the PA-infected group (P<0.05 for all cases, SNK-q test). No significant difference was observed for TNF-α expression in serum and IL-6 in lungs among ceftazidime, ATP and ceftazidime + ATP treatment groups (P>0.05 for all cases, SNK-q test). However, the ceftazidime + ATP treatment group presented a decreased expression of TNF-α in lungs and IL-6 in serum relative to ceftazidime or ATP treatment groups (P<0.05 for all cases, SNK-q test; Fig. 4).

Discussion
There has been an increasing body of evidence that ATP serves a protective role in the immune system. ATP promotes a wide range of epithelial responses through the activation of P2X and P2Y receptors in response to infectious stimuli and injuries. For instance, ATP mediates the activation of the NADPH oxidase DUOX1 through the P2Y receptor, phospholipase C and intracellular calcium signaling in the airway epithelium, leading to the activation of ERK1/2 and NF-κB pathways, IL-8 release and the production of H 2 O 2 in response to bacterial stimulations or tissue injury (24,25). In addition, ATP modulates the activities of the NLRP3 inflammasome during bacterial infection by activation of the purinergic receptor P2X7, resulting in IL-1β maturation and secretion of chemokines (26,27). Here, believed to be for the first time, the authors reported that ATP upregulated the expression of rBD-2 mRNA and enhanced the production of rBD-2 protein in response to the bacterial stimulation of PA in rat lungs. A possible mechanism underlying this process is that ATP activates the NF-κB pathway, which is critical to the production of rBD-2 (11,28) through either the P2X or P2Y receptor, however, further studies are still required in the future. In addition, enhanced rBD-2 expression by ATP reduced inflammatory damage was caused by PA in lung tissues, indicating that rBD-2 ameliorates the infection of PA. These findings are consistent with studies conducted by Wu et al (12) and Hu et al (29). The primary reasons may lie in the antibacterial function of increased rBD-2 and its regulation of inflammatory factors during bacterial infection and inflammation. Therefore, the authors analyzed the secretion of cytokines, TNF-α and IL-6 following stimulation of ATP and ceftazidime. It was identified that PA infection significantly increased pro-inflammatory factors, such as TNF-α and IL-6 in the rats, and overexpression of rBD-2 by ATP inhibited the increase of these inflammatory factors. These results are in line with a previous report by Wu et al (12), indicating that rBD-2 accelerated the weakening process of inflammation.
In summary, these results demonstrated that ATP induces rBD-2 mRNA and rBD-2 protein expression in rat lungs and    . Changes of TNF-α and IL-6 protein production measured by ELISA in rat lung and serum. Data are presented as mean ± standard deviation. TNF-α, tumor necrosis factor-α; IL, interleukin. decreases the expression of proinflammatory cytokines TNF-α and IL-6, suggesting ATP has protective effects against infection of PA. These results may provide a better understanding of the therapeutic values of ATP and hBD-2 in infectious diseases.