Long non-coding RNA GAS5 controls human embryonic stem cell self-renewal by maintaining NODAL signalling

Long non-coding RNAs (lncRNAs) are known players in the regulatory circuitry of the self-renewal in human embryonic stem cells (hESCs). However, most hESC-specific lncRNAs remain uncharacterized. Here we demonstrate that growth-arrest-specific transcript 5 (GAS5), a known tumour suppressor and growth arrest-related lncRNA, is highly expressed and directly regulated by pluripotency factors OCT4 and SOX2 in hESCs. Phenotypic analysis shows that GAS5 knockdown significantly impairs hESC self-renewal, but its overexpression significantly promotes hESC self-renewal. Using RNA sequencing and functional analysis, we demonstrate that GAS5 maintains NODAL signalling by protecting NODAL expression from miRNA-mediated degradation. Therefore, we propose that the above pluripotency factors, GAS5 and NODAL form a feed-forward signalling loop that maintains hESC self-renewal. As this regulatory function of GAS5 is stem cell specific, our findings also indicate that the functions of lncRNAs may vary in different cell types due to competing endogenous mechanisms.

mbryonic stem cells (ESCs), which are derived from the inner cell mass, are pluripotent cells that possess unlimited proliferation potential and the ability to differentiate 1 . Pluripotency is tightly controlled by core transcription factors, signalling pathways and other regulators. Among several early attempts to reveal the signalling pathways that control ESC pluripotency, members of the transforming growth factor-b (TGFb) superfamily were found to be crucial for the maintenance of the undifferentiated state 2,3 . Two signalling branches are involved in the self-renewal process in ESCs, the Nodal/Activin branch and the bone morphogenetic protein (BMP) branch. Downstream of these signalling pathways, NODAL/ACTIVIN signalling activates intracellular Smad2/3, whereas BMP signalling primarily activates Smad1/ 5/8 (ref. 4). However, compared with the BMP signalling, the regulation and function of NODAL/ACTIVIN signalling in hESCs is less been elucidated.
Non-coding RNAs (ncRNAs) were recently found to be important players in cell development, metabolism, differentiation and homoeostasis 5 . Of these, microRNAs (miRNAs) and long ncRNAs (4200 nucleotides, lncRNAs) are believed to play major regulatory roles in all multicellular organisms 6,7 . The roles of lncRNAs in hESCs are largely unclear; however, emerging evidence indicates that lncRNAs also play an essential role in regulating hESC-specific processes [8][9][10] . Many studies have shown that lncRNAs activate transcription, facilitate epigenetic modification and participate in posttranscriptional regulation in hESCs 11 . A recent work performed by our group showed that lncRNA-ROR functions as a sponge to protect the core transcription factors from miRNA binding 9 . However, the specific roles of lncRNAs in pluripotency regulation are still largely unknown.
In this study, taking advantage of high-throughput RNA sequencing technology, we investigate a set of highly expressed lncRNAs and identifiy that growth-arrest-specific transcript 5 (GAS5) correlates with hESC self-renewal. We show that GAS5 can increases OCT4, NANOG and SOX2 expression, and promotes the self-renewal of hESCs. We also show that GAS5 expression is directly controlled by the pluripotency factors OCT4 and SOX2, thus forming a circuit that promotes pluripotency. Through mechanism studies, we found that GAS5 attenuates miRNAs targeting the pluripotency-related TGFb receptor family ligand NODAL, thus maintaining NODAL expression and promoting hESC self-renewal and pluripotency. Taken together, these findings demonstrate a new pluripotency regulatory circuit that functions via a miRNA competitive mechanism mediated by the lncRNA GAS5.

Results
GAS5 is highly expressed and controlled by OCT4/SOX2. To identify lncRNAs that affect the pluripotency of hESCs, we first searched for highly expressed candidate lncRNAs via highthroughput RNA sequencing in two hESC lines (H1 and X-01). Of the annotated lncRNAs that we identified, the highly expressed GAS5 drew our interest (Fig. 1a). This formerly identified tumour suppressor and growth arrest lncRNA is extremely enriched in hESCs. We overexpressed some of the highly expressed lncRNAs and found only GAS5 remarkably and dose-dependently increased pluripotency-related OCT4, NANOG and SOX2 expression (t-test, Po0.01 Fig. 1b left panel; Supplementary Fig. 1a), while GAS5 silencing significant reduced the expression of these pluripotency genes in hESCs (t-test, Po0.01 Fig. 1b right panel; Supplementary Fig. 1b).
Of the lncRNAs reported in hESCs [12][13][14][15] , GAS5 is the most highly expressed, with an expression level similar to that of OCT4 (Fig. 1c). We performed northern blotting and quantitative PCR (qPCR) analyses to determine GAS5's expression in other cell types; however, results showed the expression in hESCs remains the highest (Fig. 1d,e). To demonstrate the dynamic changes in GAS5 expression during both embryonic development and hESC differentiation, we used fluorescence in situ hybridization (FISH). We found that GAS5 levels increased along with human embryo cleavage and that GAS5 was abundantly expressed in the cytoplasm of hESCs (Fig. 1f,g). However, its expression decreased rapidly under various differentiation conditions ( Fig. 1g and Supplementary Fig. 1d), which is similar to the previous reports demonstrating that GAS5 expression changes during mouse embryonic development 16,17 and localization of GAS5 in HeLa cells 18 . We thus hypothesized that pluripotency factors are related to the GAS5 expression. To test this hypothesis, we overexpressed OCT4, NANOG and SOX2 individually, and found that GAS5 was significantly upregulated only when OCT4 and SOX2 were overexpressed ( Supplementary Fig. 1e). We performed a chromatin immunoprecipitation (ChIP) assay along with an electrophoretic mobility shift assay (EMSA) assay to confirm the binding of OCT4 and SOX2 to the GAS5 promoter at sites specified by the JASPAR database and corresponding histone modifications as specified by ENCODE ( Fig. 1h; Supplementary  Fig. 1f; Supplementary Data 1). The results confirmed the binding of OCT4 and SOX2 to the region extending from À 400 bp to the TSS (Transcription start site), where the histone markers of transcriptional activation were most highly expressed ( Fig. 1I; Supplementary Fig. 1g-i). These findings indicate that GAS5 transcription is directly controlled by OCT4 and SOX2, and may explain why GAS5 is more highly expressed in hESCs than in other cell types. GAS5 is essential for hESC self-renewal. To determine the role of GAS5 in hESC self-renewal, we constructed lentiviral vectors expressing either GAS5 or GAS5 short hairpin RNA ( Fig. 2a; Supplementary Fig. 2a). Knockdown of GAS5 (Lenti-shGAS5) in hESCs resulted in more spontaneous differentiation (Fig. 2a, lower panels), while overexpression of GAS5 resulted in the formation of denser and larger colonies (Lenti-GAS5, Supplementary Fig. 2b). This observation was confirmed via colony formation assay (Fig. 2b), where more and larger colonies formed in Lenti-GAS5 hESCs than in those that received an empty vector (negative control, NC) group, and fewer colonies formed in Lenti-shGAS5 hESCs. By examining pluripotency gene expression, we found that Lenti-shGAS5 cells exhibited decreased pluripotency gene expression, while Lenti-GAS5 cells exhibited increased pluripotency gene expression (Fig. 2c). This effect was partially reversed by knockdown of GAS5 in Lenti-GAS5 hESCs or overexpression of GAS5 in Lenti-shGAS5 cells ( Supplementary Fig. 2c).
We next examined the impact of GAS5 on the hESC cell cycle. The proportion of G2/M phase cells was decreased in Lenti-shGAS5 cells, and the proportion of G1/G0 phase cells was significantly increased (t-test, Po0.01). Contrasting findings were noted in Lenti-GAS5 cells, indicating that GAS5 expression influences the hESC cell cycle (Fig. 2d,e; Supplementary Fig. 2d). We also employed an EdU (5-ethynyl-2 0 -deoxyuridine) incorporation assay to quantitatively assess cell proliferation. As expected, Lenti-GAS5 cells showed greater EdU incorporation than NC cells, whereas Lenti-shGAS5 cells showed less incorporation than NC cells (Fig. 2f). In addition, we noted that Cyclin and CDK gene expression was upregulated and that negative cell cycle regulator expression was downregulated in Lenti-GAS5 cells. Contrasting findings were observed in Lenti-shGAS5 cells (Fig. 2g,h). Moreover, we performed a cell cycle synchronization assay to exclude the influence of hESC heterogeneity ( Supplementary Fig. 2d). The results indicated that Lenti-shGAS5 cells recovered more slowly, while Lenti-GAS5 hESCs recovered more quickly. Together, these results confirmed the essential role of GAS5 in sustaining hESC self-renewal.
To determine the function of GAS5 in differentiation, we employed an fetal bovine serum (FBS)-induced differentiation model. We found that Lenti-GAS5 hESCs expressed fewer differentiation markers than control cells and Lenti-shGAS5 cells ( Supplementary Fig. 2f). Alkaline phosphatase (AP) staining  revealed better preservation of pluripotency in Lenti-GAS5 cells than in other cells ( Supplementary Fig. 2g). Embryo body formation assays also showed that Lenti-GAS5 cells formed denser embryoid body (EB) spheres, while Lenti-shGAS5 cells formed smaller and fewer EB spheres ( Supplementary Fig. 2b, lower panels). These results are indicative of the essential role of GAS5 in maintaining hESC pluripotency and self-renewal. However, GAS5 exerted growth-arrest effects in HEK-293T cells (Fig. 2d-f), which have been observed in other cell lines [19][20][21][22] , indicating that the abovementioned effects of GAS5 may be hESC specific.
NODAL-SMAD signalling is activated by GAS5 in hESCs. To determine the mechanism underlying the function of GAS5 in hESC self-renewal, we analysed global gene expression via RNA sequencing. We identified 1698 upregulated genes (fold change Z2) and 2,963 downregulated genes (fold change Z2) in the GAS5 knockdown group on GAS5 overexpression ( Fig. 3a; Supplementary Fig. 3a). Here we first analysed the GAS5 co-expressed genes (overlapped region in Fig.3a) using gene ontology analysis, we found that many of the differentially expressed genes are related to SMAD signalling ( Fig. 3b; Supplementary Data 2) and that many of the genes downregulated in Lenti-shGAS5 hESCs are related primarily to the cell cycle ( Supplementary Fig. 3b). Gene set enrichment analysis (GSEA) and pathway analyses also demonstrated enrichment of TGFb signalling in Lenti-GAS5 hESCs ( Fig. 3c; Supplementary  Fig. 3c). In contrast, the changes in global gene expression in Lenti-shGAS5 hESCs were negatively enriched with respect to TGFb signalling and SMAD-binding sites (Fig. 3c). However, the analysis results of GAS5 inversely correlated genes (downregulated in GAS5 overexpression and upregulated in GAS5 knockdown, Supplementary Fig. 3d,e) showed none pluripotency-related terms are enriched. Thus indicates GAS5 may played a role in TGFb-SMAD signalling in hESCs.
TGFb signalling involves various ligands and receptors 23,24 . To identify specific genes affected by GAS5, we screened various TGFb signalling-related genes and found a significant correlation between NODAL and GAS5 expression ( Fig. 3d; Supplementary  Fig. 3f). We noted elevated NODAL protein expression and increased SMAD2/3 phosphorylation on GAS5 overexpression, as well as contrasting results on GAS5 knockdown, and no significant changes in the expression of other key signalling genes ( Fig. 3e; Supplementary Fig. 3g). Because NODAL can be secreted, we measured net NODAL production in serum free, TGFb/NODAL/Activin-free medium E7. As expected, Lenti-GAS5 cells produced more NODAL in the supernatant than control cells, whereas GAS5 knockdown cells produced less than half of the normal amount (Fig. 3f).
To find additional evidence indicating that the function of GAS5 is related to NODAL signalling, we chose three hESC-related or embryonic carcinoma-related microarray data sets from the GEO database [25][26][27] that were treated with rhNODAL, a NODAL inhibitor or Activin (a cytokine that shares receptors with NODAL) and analysed gene enrichment in the TGFb and SMAD-binding subsets via GSEA. The rhNODALand Activin-treated samples exhibited enrichment in the TGFb and SMAD-binding-related subsets similar to that exhibited in Lenti-GAS5 cells, while NODAL inhibition (SB431542, an inhibitor of NODAL receptors) produced results more closely resembling those produced in Lenti-shGAS5 hESCs ( Fig. 3g; Supplementary Fig. 3h). Together, these data provide evidence that GAS5 activates NODAL signalling by promoting NODAL secretion and that activation of NODAL signalling results in phenotypes similar to that of GAS5.
NODAL is required for GAS5 to promote hESC self-renewal. Given the contrasting roles of NODAL in embryonic development, we first tested the contributions of NODAL to hESC self-renewal. We examined the cell cycle and EdU incorporation of hESCs treated with rhNODAL or SB431542, a NODAL receptor inhibitor. Similar to GAS5 overexpression, rhNODAL treatment increased hESC self-renewal, while SB431542 exerted contrasting effects (Fig. 4a,b). The effects of SB431542 were partially rescued by GAS5 transfection in both EdU and colony formation assays (Fig. 4b, lower right panel; Supplementary Fig. 4a). qPCR analysis showed that NODAL increases expression of the pluripotency factors OCT4 and NANOG in a dose-dependent manner ( Supplementary Fig. 4b); however, NODAL had little effect on SOX2 expression, which was consistent with of previous findings 28,29 .
We next tested these effects in an FBS-induced hESC differentiation model. By adding different doses of SB431542 to GAS5-overexpressing, FBS-induced hESCs, we found that GAS5 overexpression improved AP staining (Fig. 4c) and pluripotency gene expression during differentiation ( Supplementary Fig. 4c), which was reversed by NODAL inhibition in a dose-dependent manner (Fig. 4c).
Using a luciferase reporter containing multiple SMAD responsive elements (pGMSMAD-Lu), we confirmed its downstream activation in both GAS5-overexpressing and rhNO-DAL-treated cells. In contrast, GAS5 knockdown and ALK4/7 inhibition reduced luciferase activity (Fig. 4d). We next tested whether NODAL is essential for GAS5 function. We first tested the effects of SB431542 treatment in GAS5-overexpressing hESCs and found that the messenger RNA (mRNA) and protein levels of the pluripotency factors were decreased, except for those of NODAL (Fig. 4e). In GAS5 knockdown cells, we found that after transfected with different GAS5 short interfering RNAs (siRNAs). NC represents a scramble siRNA as negative control, GAS5-A/B/C/D represents different GAS5 knockdown siRNAs. **Po0.01, t-test, n ¼ 3. (c) RPKM value of several known lncRNAs and genes from the RNA sequencing were shown. (d) Northern blot analysis (n ¼ 2) using digoxingenin-labelled GAS5 anti-sense probe, the lower band indicates GAS5 transcript. Each lane is loaded with 2 mg total RNAs and a blot labelled with b-actin anti-sense probe served as loading control. (e) Relative RNA level of GAS5 in different cell types compared with HEK-293T. **Po0.01, t-test, n ¼ 3. (f) FISH analysis (n ¼ 2) of GAS5 in different human embryo stages using the digoxingenin-labelled GAS5 antisense probe. A sense strand probe was used as negative control probe. Scale bar, 50 mm. (g) FISH analysis of GAS5 in hESCs of different culture conditions (n ¼ 3). MTesR1 medium represents undifferentiated state, while À bFGF and MEF medium were used to induce hESCs differentiation. Scale bars, 100 mm. (h) A scheme of primers and probes used to detect the binding of pluripotency genes to the predicted promoter of GAS5. Lower part shows the coordinate active histone modification state from ENCODE ChIP-sequencing data. (i) ChIP analysis of the binding of OCT4 and SOX2 to GAS5 promoter using the indicated antibody with the primers in h. The analysis of NANOG promoter was shown as positive control. Data were first normalized to input then compared with IgG groups. **Po0.01, t-test, n ¼ 3. Error bars represent s.d. of the indicated experiment replicates. RNA level of b-actin served as internal reference for qPCR. See also Supplementary Fig. 1  R b 1 P 1 9 P 2 1 P 5 7 P 2 7 P 1 8 Fold change in expression     rhNODAL incorporation rescued the expression of OCT4 and NANOG (Fig. 4f). Taken together, these results indicate that GAS5 function is largely dependent on NODAL signalling. NODAL is an autocrine cytokine. Thus, we tested whether GAS5 affects NODAL secretion in hESCs via co-culture assay. First, we seeded 'donor' hESCs stably expressing Lenti-GAS5, Lenti-shGAS5 or control virus into the lower chamber and allowed the cells to reach at least 60% confluence to ensure sufficient NODAL secretion. Then, we seeded 'receiver' hESCs into the upper chamber and co-cultured the cells for an additional 3 days. The E7 medium contained no extra NODAL or Activin ( Relative SMAD-reporter activities co-culture with GAS5-overexpressing cells increased NODAL and pluripotency factor mRNA levels compared with the cells infected with the control vector. ALK receptor inhibition by SB431542 attenuated these effects (Fig. 4g). NODAL protein level elevation in the supernatant was also confirmed by enzyme-linked immunosorbent assay (ELISA; Fig. 4h). Collectively, these data demonstrate that NODAL secretion is required for GAS5 to enhance hESC self-renewal.
NODAL is post-transcriptionally regulated by GAS5. Several reports have demonstrated the involvement of GAS5 in the mammalian target of rapamycin (mTOR) and glucocorticoid receptor (GR) pathways 18,30 . Therefore, we aimed to determine whether the function of GAS5 in hESC self-renewal is related to these pathways. Because the expression of GAS5 can be elevated by mTOR inhibition 30 , we first tested whether mTOR inhibition by the mTOR inhibitor rapamycin would increase the expression of pluripotency genes. The results showed that pluripotency gene expression and NODAL levels did not significantly decrease until day 3 of treatment (t-test Supplementary Fig. 5a), which is inconsistent with the effects of GAS5. GR function can be inhibited by GAS5, so we performed GR knockdown and activation and used dexamethasone, a GR agonist, to analyse the effects. The results showed that neither GR knockdown nor stimulation significantly affected pluripotency gene expression or NODAL expression (t-test. Supplementary Fig. 5b,c). In addition, previously identified GR-binding element mutations 18 in the GAS5 transcript did not prohibit NODAL and NANOG upregulation in GAS5-overexpressing hESCs (Fig. 5a). On the basis of these results, we hypothesized that other GAS5-related mechanisms modulate NODAL expression.
We searched for the functional segment of GAS5 by constructing several vectors encoding different truncated fragments of GAS5 (Fig. 5a, left panel). qPCR results showed that the GAS5 transcript increases NODAL and Nanog expression in a length-dependent manner (Fig. 5a, right panel) and that the segment from 1-251 bp has no function. Next, we investigated the possibility that GAS5 directly regulates NODAL transcription. We tested the luciferase activity of the NODAL promoter in both hESCs and HEK-293T cells and noted no significant changes in luciferase activity, when GAS5 was overexpressed or silenced (t-test Fig. 5b; Supplementary  Fig. 5d). We also analysed the asymmetric enhancer region of Nodal to determine whether GAS5 promotes NODAL expression by transcriptional regulation (Supplementary Fig. 5e). The results showed that GAS5 modulation does not alter asymmetric enhancer reporter luciferase activity in hESCs. We also assessed GAS5 transcript's localization in either GAS-overexpressing or GAS5 knockdown hESCs and found that the levels of cytoplasmic GAS5 transcripts were significantly affected (t-test, Po0.01 Fig. 5c). Taken together, we hypothesized that indicating that cytoplasmic GAS5 may play an important role in modulating NODAL expression.
The regulatory mechanisms of lncRNA in the cell cytoplasm involve RNA-RNA or RNA-protein interactions. To examine direct RNA-RNA interactions, we investigated the possibility of sequence homology between the GAS5 and NODAL transcripts and found none ( Supplementary Fig. 5f). To examine RNAprotein interactions, we performed an RNA immunoprecipitation (RIP) assay using an MS2-binding protein (MS2bp) system, in which tagged MS2bp specifically binds RNA containing MS2-binding sequences (MS2bs; Fig. 5d). We constructed an MS2bs-GAS5 vector and used MS2bp to pulldown the transcript. ARTICLE SDS-polyacrylamide gel electrophoresis (SDS-PAGE) and silver staining were performed to identify a differential protein band (Fig. 5e, indicated by a solid arrow), and subsequent liquid chromatography-mass spectrometry/mass spectrometry proteomic analyses were performed to analyse the content of this band (Fig. 5e). By analysing these data, we found that AGO2, a key component of the miRNA-mediated silencing complex, is more enriched in MS2-GAS5 pulldown products than other control products (Fig. 5f)  differentially stained band, we performed label-free highthroughput proteomics analysis using MS2-GAS5 pulldown products. Gene ontology analysis of all predicted pulldown proteins showed enrichment of translation regulation-related and miRNA-mediated silencing-related proteins (Supplementary Fig. 5g; Supplementary Data 3), supporting the idea that an AGO2-and miRNA-based mechanism may underlie the function of GAS5. RNA FISH analysis combined with immunofluorescence using a specific antibody also showed that GAS5 and AGO2 co-localized in the cytoplasm of hESCs (Fig. 5g) and that GAS5 does not directly bind to its downstream effectors NODAL and SMAD2/3 (Fig. 5h). These data suggest that GAS5 binds to AGO2 and may be associated with miRNA-mediated silencing.
Previous reports have demonstrated that lncRNA and other transcripts serve as 'sponges' that regulate target transcripts by competing with post-transcriptional elements. We therefore employed AGO2-RIP to test whether differential GAS5 expression affects the binding of AGO2 to NODAL transcripts. The results showed that on GAS5 overexpression, AGO2-NODAL binding was significantly decreased, while on GAS5 silencing, this binding was significantly increased (t-test, Po0.01 Fig. 5i). These findings suggest that GAS5 regulates NODAL by competing with AGO2 post transcriptionally.
GAS5-interacting miRNAs regulate NODAL expression. The binding of AGO2 to a specific transcript requires specific miRNA. Therefore, we first used lentivirus-mediated knockdown of Dicer expression and found global miRNA expression reductions in hESCs compared with control cells (Supplementary Fig. 6a). Simultaneously, the GAS5 knockdown induced decreases in NODAL expression were attenuated in Dicer knockdown hESCs (Fig. 6c), which indicated a key role of miRNA-mediated regulation in this process. To identify GAS5-and NODALinteracting miRNAs, we analysed the miRNome via highthroughput sequencing and constructed a network, in which the abundances of both target genes and miRNAs were considered. To scale down the network, we included only the top 30 miRNAs that were predicted to target GAS5 (Fig. 6a). We found that seven miRNAs in this network also targeted NODAL (Fig. 6b).
To identify the functions of these candidate miRNAs, we transfected their mimics or inhibitors directly and analysed the changes in GAS5 and NODAL expression in hESCs. The results showed that miR-2467-5p, -3200-3p and -let-7a/e-5p overexpression significantly reduced both NODAL and GAS5 RNA levels (t-test, Po0.01 Fig. 6d; Supplementary Fig. 6b), whereas their inhibition caused RNA upregulation ( Supplementary  Fig. 6c). Validation of target sites were done using luciferase reporter system, including NODAL wild-type and site-mutated constructs ( Fig. 6e; Supplementary Data 4). We found that these miRNAs significantly suppressed the luciferase activities of the wild-type reporters (t-test, Po0.01), but not those of the site-mutated reporters, except for a reporter featuring mutations at site 2 of miR-2467 in NODAL (Fig. 6f). To validate the direct binding of the predicted miRNAs to GAS5 transcripts, we performed an MS2-mediated RIP analysis using constructs containing either wild-type GAS5 or mutated GAS5 (Fig. 6e) fused with MS2bs. qPCR analysis showed that miR-2467-5p, miR-3200-3p and let-7a/e-5p were significantly enriched in the wild-type GAS5 and NODAL transcripts, and that no significant enrichment occurred in the relative site-mutated GAS5 transcripts or in a GAS5 transcript, in which all miRNA-binding sites were mutated (Fig. 6g). The binding of miR-145-5p to OCT4 3 0 -untranslated region (UTR) served as a positive control. We again performed an Ago2-mediated pulldown (Fig. 6h) and an MS2-mediated pulldown (Fig. 6i) to determine whether all mutated GAS5 sites still affected NODAL expression. As expected, the mutated GAS5 transcript did not compete with Ago2 to bind the NODAL transcript. These data provide evidence of the direct binding of the candidate miRNAs to GAS5 and NODAL mRNA.
How GAS5 affects these miRNAs remains unknown. Here we first performed northern blotting and confirmed their existence ( Supplementary Fig. 6e). Using qPCR analysis, we noted downregulation of the validated miRNAs in conjunction with wild-type GAS5 overexpression, but not GAS5mut overexpression (t-test, Po0.01 Supplementary Fig. 6f), indicating that miRNA degradation mechanism may be associated with this interaction. To clarify this point, we explored the turnover of these RNAs using actinomycin D, an inhibitor of de novo RNA transcription. The turnover of validated miRNAs was accelerated on GAS5 expression ( Supplementary Fig. 6g,h). At the same time, overexpression of miR-2467-5p, miR-3200-3p or let-7a/e-5p also promoted GAS5 turnover (Supplementary Fig. 6i). These data suggest that the GAS5 transcript post-transcriptionally regulates the levels of miR-2467-5p, miR-3200-3p and let-7a/e-5p, and promotes their degradation.
MiR-2467/3200/Let7e are essential for GAS5 function. To investigate the functions of these miRNAs in hESC self-renewal, we employed colony formation assays and cell cycle analysis. The results showed that overexpression of these miRNAs significantly disrupted hESC self-renewal (Fig. 7a,b). Increased expression of the pluripotency factors OCT4 and NANOG, but not SOX2, was facilitated via inhibition of these miRNAs ( Supplementary  Fig. 7a). Unlike the wild-type transcript, the GAS5-mut transcript exerted no effects (Fig. 7a,b; Supplementary Fig. 7b). We found that miR-2467-5p, miR-3200-3p and let-7a/e-5p did not directly target these pluripotency factors, using dual luciferase assays ( Supplementary Fig. 7c). Because miR-2467-5p, miR-3200-3p and let-7a/e-5p facilitated expression of similar phenotypes resembling NODAL inhibition, we investigated whether NODAL rescues the effects of these miRNAs. We found that rhNODAL rescues pluripotency gene downregulation caused by individual miRNA overexpression (Fig. 7c), suggesting that NODAL is the critical target of these miRNAs in hESCs. To determine whether these miRNAs are essential for GAS5-mediated NODAL expression, we transfected these miRNAs mimics into GAS5-overexpressing hESCs, and found that both pluripotency gene expression and NODAL expression were attenuated compared with control cells (Fig. 7d). Western blot analysis confirmed this finding (Fig. 7e). On the basis of these data, we concluded that miR-2467-5p, miR-3200-3p and let-7a/e-5p play important roles in GAS5-mediated NODAL regulation.
The biological significance of miR-2467-5p, miR-3200-3p, let-7a/e-5p and their regulatory transcripts is unknown. Therefore, we detected the changes in their endogenous expression during hESC differentiation. Real-time analysis of both undifferentiated (H1 and H9) and differentiated depletion of basic fibroblast growth factor ( À bFGF) and mouse embryonic fibroblasts (MEF) induction) cells showed that miR-2467-5p, miR-3200-3p and let-7a/e-5p expression was significantly upregulated during differentiation (t-test, Po0.01 Fig. 7f), accompanied by gradual reductions in GAS5, NODAL signalling-related gene and pluripotency gene expression (Fig. 7g), indicating the existence of an inverse relationship between the expression of these miRNAs and GAS5-mediated NODAL signalling.
We validated the function and mechanism of GAS5 in hESCs and demonstrated the existence of a regulatory loop, in which pluripotency-driven GAS5 promoted NODAL signalling by competing with NODAL-targeting miRNAs. However, to gain additional in vivo evidence, we performed teratoma formation analysis. By injecting stably transfected hESCs into the backs of nude mice, we found that GAS5-overexpressing hESCs formed larger teratomas than control cells, while mutant GAS5-overexpressng cells did not form significantly larger tumours. GAS5 knockdown resulted in reduced tumour size (Fig. 7h). We also found that GAS5 increased NODAL expression, while reducing the expression of its targeting miRNAs in these tumours (Fig. 7i).
Specific interactomes are essential for GAS5 function. To explore the function of GAS5 in other cells, we first evaluated the expression levels of the above miRNAs and NODAL in other cell lines, and found that although these interacting miRNAs are differentially expressed in different cell lines (Fig. 8a), the expression of NODAL was consistently lower in other cell lines (Fig. 8b). Furthermore, we found that GAS5 did not promote pluripotency and NODAL expression in cells other than hESCs and umbilical cord mesenchymal stem cells (uMSCs) (Fig. 8d).
As previous studies have shown that NODAL is not expressed in normal tissues 28 and that specific ratios of competing endogenous RNA (ceRNA)/target RNA are critical for the formation of functional competing endogenous pairs 31,32 , we hypothesized that lack of NODAL expression may be the cause of the loss of GAS5 function in other cells. To verify this hypothesis, we evaluated NODAL secretion levels in these cell types and found that only hESCs and uMSCs secrete detectable levels of NODAL (Fig. 8d), and that only in uMSCs does GAS5 induce NODAL section (Fig. 8e). We measured the copy number changes of these miRNAs and NODAL after transfecting different cells with different concentrations of in vitro-synthesized full-length GAS5 transcripts. The results showed that NODAL expression is initially upregulated after transfection and increases in the presence of higher GAS5 concentrations, whereas the expression of their interacting miRNAs decreases continuously (Fig. 8f). Similar expression patterns were observed in uMSCs (Fig. 8g). Although miRNA expression was gradually downregulated, NODAL expression did not change due to its low or non-existent basal expression level in HEK-293T cells (Fig. 8h), which supported our hypothesis. Here we demonstrated that a cell-specific GAS5 regulation circuit fine tunes NODAL expression via interacting miRNAs, thus promoting hESC selfrenewal and pluripotency (Fig. 8i).

Discussion
LncRNAs are newly identified players in pluripotency and self-renewal regulation in hESCs. Many studies have examined the functions of long ncRNAs in both ESC lineage commitment and the self-renewal 12,33 . GAS5 was identified decades ago 34 ; however, its function and working mechanism were not identified until recently in tumour suppression 21,[35][36][37][38] , growth arrest 18,22,38,39 and affected by rapamycin 40,41 . Our finding that GAS5 is highly expressed in hESCs appears to be arbitrary given its significant effect on tumour suppression and proliferation inhibition. However, its role in development was indicated decades ago 16,17,34,39 ; among the GAS gene family members, GAS5 42 and GAS6 17 were found to be upregulated in ESCs, and only GAS5 42 was highly expressed in hESCs. Although previous studies have suggested a relationship between GAS5 and hESC self-renewal maintenance, no experiments regarding the function of GAS5 in hESCs have been conducted. To address this idea, we investigated GAS5 expression in hESCs and found that GAS5 maintains hESC pluripotency and self-renewal by promoting NODAL signalling. Because previous reports showed that GAS5 promoted apoptosis, and growth inhibition of various tumour cells 18,38,43,44 and normal cells 40,45 , it is interesting that GAS5 also promotes hESC proliferation, which may due to its modulation to NODAL signalling.
In identifying the signalling pathways that control ESC pluripotency, the members of the TGFb superfamily were crucial for maintenance of the undifferentiated state 2,46 . Unlike mESCs, NODAL is important for hESC self-renewal 47 . NODAL is highly expressed in hESCs and is rapidly downregulated during differentiation 48 . Its overexpression inhibits hESC differentiation into neuro-ectoderm and maintains pluripotency markers 49 . Previous reports 29,50,51 have also shown that the NODAL signalling pathway directly controls NANOG expression in hESCs and mouse epiblast stem cells (mEpiSCs), thereby affecting hESC pluripotency and self-renewal. Intriguingly, an important feature of the NODAL pathway is its complete lack of activity in normal tissues 52 , making its role in vertebrate developmental biology more mysterious. In our study, we showed that GAS5 regulates NODAL expression through a competing endogenous network controlled by miRNAs. However, we also tested this interaction in other cell types and did not observe elevated NODAL expression on GAS5 overexpression, suggesting that GAS5 functions differently in other cell types than in hESCs. This phenomenon may be partly explained by the complete lack of NODAL expression in differentiated cells; accordingly, we did not detect NODAL secretion in the supernatants of these cells.
Another possible reason for the unique functions of GAS5 may also be tied to its involvement in competing endogenous networks. Highly expressed site-containing RNAs, whether naturally occurring or delivered as research reagents, can act as 'sponges' to titrate miRNAs away from other normal targets [53][54][55][56][57] . Although many studies have identified competing endogenous networks in various cell types, controversy remains regarding the biological significance of these networks in primary cells 32 .  A major concern is the quantity of these 'sponges' in competing endogenous pairs. However, global RNA sequencing showed that GAS5 is the most abundant lncRNA and it is among the top 50 expressed RNAs in hESCs, with an RPKM (reads per kilobase of transcript per million reads mapped) similar to that of OCT4, which is the most highly expressed pluripotency gene in hESCs.
In contrast, NODAL expression is 20 times lower than that of GAS5, and all known interacting miRNAs have copy numbers that are higher than NODAL, but lower than GAS5. Our findings indicate that GAS5-mediated promotion of hESC self-renewal via NODAL expression elevation may depend on certain miRNomes or target ceRNA expression to form this regulatory loop. However, confirming this hypothesis may require additional studies comparing miRNome and transcriptome changes in various cell types.
Here we found that GAS5 not only sponges, but also downregulates and increases the degradation of its target miRNAs. This could be an intriguing issue concerning ceRNA or may be all miRNA-related mechanism. Traditionally, binding of miRNAs to target sites causes translational repression or RNA degradation. Perfect matches appear more likely to cause RNA degradation, and seed region matches appear more likely to cause translational repression. However, such criteria are also changing quickly that can be inferred from previous reports 58, 59 , for miRNA causes deadenylation and decapping of the target RNA even in an imperfect binding motif. So it is better to understand that miRNA's binding to target RNA could cause not only the degradation of target, but also the instability of miRNA itself [60][61][62] , which indicates a mutual degradation mechanism exists. Nevertheless, such mechanism is still controversial and needs further investigation to fully reveal this phenomenon.
Since miRNAs are critical in forming ceRNA network, it is hard to fully classify and standardize ceRNA due to the diverse functions of miRNAs. However, whether miRNA and its target degrade or not is certainly not a required factor in defining ceRNA. The fundamental factor in forming a ceRNA pair is that both RNAs share mutual miRNA and one can affect the expression of other target. Since the quantity strongly affects the possible function of ceRNA, it is accepted that ceRNA must have greater expression than the miRNAs and the other functional targets. Still, to make certain such regulation is caused by ceRNA mechanism, one must exclude other possibilities. Our proposed GAS5-miRNA-NODAL model fits these criteria for both GAS5 and NODAL is targeted by these miRNAs, we found and GAS5 is much higher in expression than the others. However, to fully reveal the whole ceRNA network in ESCs and to further define the underling mechanism in such mutual degradation warrant further study.

Methods
Informed consent. This study was approved by the Ethics Committee of Second Military Medical University. All embryos were obtained with written informed consent signed by the donors voluntarily for research on human early embryonic development mechanisms with no financial payment.
Cell culture and embryo collection. Human ESC lines H1, H9 (purchased from WiCell Research Institute, ID: WAe001-A for H1 and WAe009-A for H9) and X-01 (ref. 63); provided by Professor Lei Xiao from Shanghai Institutes for Biological Sciences, ID: SDSHES-01) were routinely maintained on feeder-free cultures with commercial available maintaining medium as instructed. E8, E7 or mTeSR1 (Stemcell Technologies) medium were used for feeder-free cultures as indicated respectively with daily change. Cells were passaged using Gentle Cell Dissociation Reagent (Stemcell Technologies) when 70% clonal confluence was reached. Differentiation by forming EB suspension was carried out in a low attachment six-well plate (Corning) with hESCs medium. Undirected differentiation of hESCs was performed in a feeder-free condition using DMEM supplied with 10% fFBS, 1 mM glutamine, 1% nonessential amino acids (all form Invitrogen) and 0.1 mM b-mercaptoethanol (Sigma-Aldrich). Another alternative differentiation method of feeder-free hESCs involved the use of conditioned hESC medium deprived of bFGF. Cells were routinely tested for mycoplasma contamination, and no contamination was found during the experiments.
Embryos were donated for this study with the informed consent of healthy couples, who already had a healthy baby with no genetic disorders from the same in vitro fertilization cycle and wished to donate the remaining cryopreserved embryos. The cleavage stage embryos were collected using fine needle aspiration after fertilization and were stored in liquid nitrogen storage tank. On usage, the embryos were thawed rapidly by taking straws from the liquid nitrogen storage tank, and embryos of different stages were randomly collected for further analysis by the embryologists at the clinical in vitro fertilization lab of the Center for Reproductive Medicine, Renai Hospital, Shanghai Jiao Tong University. All human-related procedures were carried out in accordance with guidelines approved by the Ethics Committee for Research on Human Subjects of Second Military Medical University. Samples were used and reported anonymously with the experimental findings.
Immunofluorescence and FISH. Samples were fixed with 4% paraformaldehyde. For the detection of proteins, samples were blocked with 5% bovine serum albumin, than incubated with rabbit polyclonal anti-Ago2 (1:100, ab32381, Abcam). For the detection of GAS5 and DANCR, samples were further hybridized with commercial available SuperHyb Solution (Tiandz Inc. No.130906-10) containing respective RNA probes. RNA probes were used and labelled with digoxigenin (DIG)-UTP (Roche) using the mMESSAGE T7 Ultra In Vitro Transcription Kit (Ambion) in accordance with the manufacturer's directions.
Lentiviral transduction and selection in hESCs. Lentivirus production and tittering were done by GenePharma Corp, all lentiviruses contain a puromycin selection tag. For lentivirus transduction in hESCs, the lentiviruses were added the second day of hESC passage at multiplicity of infection 10, and continued to be added in the medium with daily change for another 2 days. To generate stable colonies, puromycin was added on the fifth day of passage at the concentration of 100 ng ml À 1 for selection, and hESCs colonies that survived were passaged and analysed by either quantitative PCR with reverse transcription (qRT-PCR) or fluorescence microscopy.
Vector construction and small RNA synthesis. The complementary DNA (cDNA) encoding lncRNA-GAS5, LOC100506647, LINC00458 or lncRNA-DANCR was PCR-amplified by the Pfu Ultra II Fusion HS DNA Polymerase (Stratagene) and subcloned into the XhoI and EcoRI sites of pcDNA3.1 vector (Invitrogen). The pcDNA3.1-GAS5 with point mutations in the DBD mimic sequence, miR-2467, miR-3200 and let-7a/e response elements was synthesized using a QuikChange Site-Directed Mutagenesis kit (Stratagene). Vector pcDNA 3.0-MS2 (12X) was double digested with NotI and XhoI, and the GAS5 wild-type and mutated fragment was synthesized and subcloned into this construct. The3-0 UTR of NODAL were amplified using PCR and subcloned into the pMir-Report vector (Promega) for luciferase reporter assay using the one-step-directed cloning kit (Novoprotein, Shanghai, China). The 3 0 -UTR of NODAL mRNA containing the miRNA-binding sites were mutated and subcloned into pMir-Report.
SiRNAs, miRNA mimics, miRNA inhibitors and other oligoes used in vector construction were all designed and synthesized by GenePharma Corp. RNAiMAX (Invitrogen) was used for small RNA transfection. SMAD luciferase reporter (pGMSMAD-Lu) was purchased from Genomeditech (USA) where multiple SMAD responsive element was inserted before the luciferase gene.
Luciferase reporter transfection and dual luciferase assay. Luciferase reporter transfection and dual luciferase assay was performed as following. In brief, HEK-293T cells were plated in 96-well plates and transfected using Lipofectamine 2000 (Invitrogen) with 50 ng pMir-Report vector (carrying firefly luciferase, Promega) inserted with indicated target sequence, empty vector was served as a control, and a PRL-TK vector (carrying Renilla luciferase, Promega) was transfected simultaneously to serve as internal control. Twenty four hours after transfection, cells were lysed and subjected to luciferase assay according to the manufacturer's protocol. For luciferase reporter assay in hESCs, cells were seeded in 48 well, while other procedures were accordingly. The details of construction of reporter vectors and transfections are shown in the Supplementary Experimental Procedures.
ChIP and RIP assay. ChIP assays were performed in accordance with the manufacturer's instructions of the EZ-Magna ChIP A/G Kit (Millipore, Billerica, MA, USA). The MS2bp/MS2bs-based RIP assay was performed using the EZ-Magna RIP Kit (Millipore) in accordance with the manufacturer's instructions. MS2bp/MS2bs transfected hESCs were collected and lysed using the kit and an anti-GFP antibody (5 mg per assay, Abcam, # ab290) were ligated to the magnetic beads for RNA pulldown.
FACS and cell cycle analysis. FACS and flow cytometry analysis were performed according to the standard protocol. For cell cycle analysis, each 1 ml of cell suspension (1-5 Â 10 5 cells) was incubated with the (propidium iodide (PI), Cell Signalling Technology). For EdU analysis, 1-5 Â 10 5 cells from each sample were processed with EdU (EdU Detection Kit, Ribobio), in accordance with the manufacturer's instructions.
Cell cycle synchronization assay. The cells were treated with NOCODAZOL (Beyotime, Shanghai, China) at concentration of 100 ng ml À 1 for 12 h when 80% confluent reached. Two hours before experiment, the cells are washed in phosphate-buffered saline for two times and changed for new medium to release the NOCODAZOL. Cells were collected at different time points and fixed according to standard protocol 64 . Cell cycle were analysed using BD FACs-Calibur (BD Biosciences).
Total RNA isolation and RT-qPCR. Total RNAs from cultured cell lines were extracted with Trizol (Invitrogen). For RT-PCR, after treatment with DNase I (Ambion, DNA-free kit), the complementary DNA was transcribed with Revertra Ace (TOYOBO). For qPCR, the relative expression of different sets of genes was quantified to GAPDH mRNA.
Northern blotting. Northern blots were performed with DIG-labelled anti-sense probes (DIG Northern Starter Kit, Roche). In brief, total RNA isolated from various cells using Trizol (Invitrogen) was separated by 12% denaturing RNA PAGE and transferred to a positively charged nylon membrane (GE Healthcare) by wet electro-transfer blotting. Membranes were crosslinked by ultraviolet for 3 min, prehybridized for 1 h at 55°C and hybridized overnight at 60°C with GAS5 probes complementary to GAS5. After hybridization, membranes were washed twice 10 min with 2 Â saline-sodium citrate buffer (SSC) and once 10 min with 0.2 Â SSC. AP conjugated anti-DIG antibody is used for secondary labelling. CDP-STAR (Roche) was used for exposure using an Illuminator Chemiluminescent Detection System (Stratagene). For miRNA detection, total miRNAs were extracted using miRcute miRNA extraction kit (Tiangen), in vitro-transcribed anti-sense miRNA probes were generated following instruction, the hybridization temperature was adjusted to 55°C overnight. Uncropped images were shown in Supplementary Fig. 9.
NODAL ELISA. NODAL levels in the normal culture medium collected for 24 h from different treatment groups were detected with the Human NODAL ELISA Kit (Cusabio, Wuhan, China) according to the manufacturer's instructions. The ELISA kit is tested using the recombinant full-length human NODAL (Abnova, CA, USA) to assess the reliability of the kit.
Teratoma formation assay. All procedures relating to animal subjects were performed under ARRIVE guideline and institutionally approved protocol deemed in accordance with the guidelines of the Institute of Laboratory Animal Resources, the Second Military Medicine University. Twelve 6-week-old male NOD/SCID mice were obtained from Laboratory Animal Resources. All animals were housed in a specific pathogen-free environment with 12-h photoperiods and ad libitum access to standard chow and water. In vivo pluripotency of the treated hESCs were tested as previously described 46 . In brief, H1 hESCs were collected with collagenase IV and were implanted underneath the back skin of the mice. The mice were randomly divided into each group without blinding to receive implantation. Teratoma growth was determined by weekly observation and palpation. Mice were sacrificed killed 6 weeks after implantation. Teratomas were washed weighed and photographed, than part of the tumour were minced for RNA extraction while other part were preserved for future use. The size of the tumours generated in this experiment were within the limit of those allowed in the ethical guidelines of our institution. A standard ruler was used to define the size of the teratomas (0.5 mm per minor mark).
High-throughput RNA-sequencing and computational analysis. For RNA sequencing, we aligned 100-bp paired-end sequencing reads to the human genome (hg19) using Tophat2/Bowtie2 allowing for five mismatches. We identified read-pair mappings to gene structures derived from RefSeq using the summarize Overlaps function with mode Intersect Strict (Genomic Ranges, Bioconductor). Using these raw counts, we identified genes expressed across the sample groups. We removed genes with a count o10 across all samples prior to statistical analysis. The differential analysis was carried out using edgeR, applying TMM (trimmed Mean of M-values) library normalization and a 0.001 false discovery rate to select expressed transcripts. In addition, we calculated gene-level RPKM values for the same gene set using Cufflinks.
The likelihood of the binding of a mature miRNA to a certain gene's mRNA was evaluated using the miRanda miRNA Target Detection Software (August 2010 Release, http://www.microrna.org/microrna/getDownloads.do).
We used GSEA v2.0 to perform GSEA on various gene signatures. Gene sets were either obtained from the MSigDB database v4.0 or from published gene signatures. Statistical significance was assessed by comparing the enrichment score to enrichment results generated from 1,000 random permutations of the gene set to obtain P values (nominal P value).
For miRNA: mRNA interacting network construction, we first ranked the top upregulated mRNAs and top downregulated miRNAs after GAS5 overexpression, and excluded low expression mRNAs (o25 normalized counts) and miRNAs (o25 normalized counts). MiRanda algorithm was used to predict miRNAbinding sites and only miRNAs that binds to GAS5 transcript were selected to form the network.
Isolation of cytoplasm and nuclear RNA. Cytoplasmic and nuclear RNA were isolated and purified using the PARIS Kit (Ambion, CA, USA) according to the manufacturer's instructions. Extracted RNAs were put to reverse transcription immediately for qPCR analysis.
Statistical analysis. Each experiment was repeated two or three times or more as mentioned in each figure legend. Data are presented as mean±s.d., unless stated otherwise. Student's t-test (unpaired, two tailed) was used to compare two groups for independent samples. Po0.05 was taken to indicate statistical significance. No statistical method was used to predetermine sample size for these experiments. All observed data were included without exclusion criteria. Data collection was conducted randomly. Uncropped scans of the most important blots are shown in Supplementary Fig. 8b in the Supplementary Information Section, all oligos used in the article are listed in the Supplementary Data 5.
Data availability. Data for mRNA expression of either GAS5 overexpressed, knockdown or NC hESCs, and miRNA expression of GAS5 overexpressed or NC hESCs have been deposited in the Gene Expression Omnibus database under accession code GSE66993. The authors declare that all data supporting the findings of this study are available within the article and its Supplementary Information Files or from the corresponding author on reasonable request.