The juvenility-associated long noncoding RNA Gm14230 maintains cellular juvenescence

Young animals possess unique traits, such as a higher growth rate, quicker wound healing and faster learning, and also show a higher capability for tissue plasticity and regeneration. These juvenile traits are helpful for the emergence of the mature individual and could provide a therapeutic strategy for intractable diseases. However, the biological and molecular bases behind the regulation of the juvenile traits have been poorly investigated. We previously identified juvenility-associated genes (JAGs) that are highly expressed selectively in juvenile cardiomyocytes and hepatocytes in mice (Jam et al., 2018). Here, we analyzed the noncoding elements of the transcriptome to obtain a comprehensive landscape of the transcriptome regulating cellular juvenescence.

lncRNAs are transcripts of more than 200 nucleotides that do not encode for proteins (Kopp and Mendell, 2018; Mercer et al., 2009; Ponting et al., 2009; Rinn and Chang, 2012; Wilusz et al., 2009). Preceding studies have reported that the genome is pervasively transcribed and contains numerous transcriptional units with unknown function (Bertone et al., 2004; Hangauer et al., 2013; Kapranov et al., 2007a,b; Xu et al., 2009). The current understanding of the expression patterns and functions of the lncRNA molecular group is incomplete, although recent studies have identified lncRNAs in the human transcriptome with accurate 5′ ends (Hon et al., 2017).

During cellular senescence, the lncRNA MIR31HG de-represses the transcription of p16^INK4A (also known as CDKN2A), thus mediating B-RAF-induced cellular senescence (Vanderschuren et al., 1997). The PANDAR lncRNA interacts with scaffold-attachment-factor A (SAFA, also known as HNRNPU) and polycomb repressive complexes (PRCs), thus suppressing senescence in proliferating cells (Puvvula et al., 2014).

Numerous lncRNAs function as ‘scaffolds’ for chromatin-modifying complexes (Khalil et al., 2009). The Xist lncRNA functions by recruiting polycomb repressive complex 2 (PRC2) to one of the X chromosomes in cells from females, repressing its transcription and dose-compensating X chromosome gene expression (Brockdorff et al., 1991; Zhao et al., 2008). Hotair, which is encoded by the HoxC locus, represses HoxD genes by recruiting PRC2 to the HoxD locus (Rinn et al., 2007). Braveheart (Bvht) interacts with Suz12, a PRC2 component, during cardiomyocyte differentiation, thereby regulating gene expression (Klattenhoff et al., 2013). PANDAR recruits PRC1 and PRC2 to specific gene loci, including CDKN1A (Puvvula et al., 2014). Fendrr binds the PRC2 and TrxG/MLL complexes, thereby regulating gene function during heart development (Grote et al., 2013). Airn recruits the histone methyltransferase G9a (also known as Ehmt2) to a target gene locus for gene silencing (Braidotti et al., 2004; Nagano et al., 2008).

Although the molecular machinery by which lncRNAs perform their functions have been investigated as described above, knowledge regarding the age-dependent functions of lncRNAs, particularly their functions during the juvenile phase, is limited. Here, we established a lncRNA analysis bioinformatics pipeline and identified juvenile-selective lncRNAs. These lncRNAs, which we call juvenility-associated lncRNAs (JAlncs), constitute a transcriptional network responsible for the juvenile properties of a cell. Furthermore, Gm14230 (ENSMUSG00000087179), which is an evolutionarily conserved JAlnc, maintains the juvenile properties of cells, and the loss of Gm14230 promotes cellular senescence. Thus, our lncRNA analysis identified the noncoding elements in the transcriptome network controlling the juvenile properties of cells.

To identify lncRNAs that are highly expressed selectively in juvenile cells, we established a pipeline to analyze the nucleotide sequence of the noncoding transcripts extracted from mouse hepatocytes and cardiomyocytes at different times during postnatal development (Fig. 1A). The original dataset was generated in a previous work of ours (Jam et al., 2018) and deposited in the DNA Data Bank of Japan (DDBJ) under the accession number DRA007101 (https://ddbj.nig.ac.jp/DRASearch/submission?acc=DRA007101). Thus, we analyzed the lncRNAs in the transcriptome of hepatocytes and cardiomyocytes isolated from juvenile [postnatal day (P)1 and P7] and adult (P56) mice (n=3, Fig. 1A). P1 and P7 correspond to an early and late neonatal phase in mice, respectively. Two points (P1 and P7) rather than one were included to capture consistent alterations in comparison to the adult phase (P56). The Ensembl database (version 91.38) (Zerbino et al., 2018) covers 5034 lncRNAs in the mouse transcriptome. We defined the juvenility-associated lncRNAs (JAlncs) as lncRNAs that were more highly expressed in juvenile (both at P1 and P7) cells than adult (P56) cells by more than 2-fold. Consequently, we identified 327 JAlncs in the hepatocytes (hepato-JAlncs) and 262 JAlncs in the cardiomyocytes (cardio-JAlncs) (Fig. 1A).

Then, we analyzed the components of the hepato- and cardio-JAlncs. The hepato-JAlncs included functionally well annotated lncRNAs, such as H19, Lockd and Snhg5 (Fig. 1B). The adult (P56)-expressed lncRNAs included Lncbate1 (Fig. 1B). The cardio-JAlncs also included functionally annotated lncRNAs, such as H19, Xist, Jpx and Crnde (Fig. 1C). The results of the previous RNA-seq analysis (Jam et al., 2018) were validated by quantitative PCR (qPCR) (Fig. 1D). Thus, we identified JAlncs in juvenile hepatocytes and cardiomyocytes in mice.

To identify lncRNAs that play an essential role in the juvenile physiology, we focused on the JAlncs common (shared) between the hepato- and the cardio-JAlncs. While numerous lncRNAs are expressed in a tissue-specific manner (Cabili et al., 2015; Washietl et al., 2014), others including MALAT1 are broadly expressed (Ji et al., 2003; Nakagawa et al., 2012), indicating their important function for all types of cells. A total of 62 lncRNAs were identified as the common JAlncs by retrieving the overlapping lncRNAs between the hepato- and cardio-JAlncs (Fig. 2A,B). The expression data of the 62 common JAlncs are shown in Table S1. The common JAlncs included H19, Lockd, Snhg5, Crnde, Meg3 and Rian (Fig. 2C).

The common JAlncs also included numerous lncRNAs with unknown functions. With the aim to discover novel lncRNAs and their functions to broaden our understanding of the molecular machinery of the juvenile transcriptome, we focused on the JAlncs with an unknown function. Based on the hypothesis that functionally important lncRNAs are evolutionarily better conserved, we investigated interspecies conservation of the common JAlncs. A BLAST search identified a homologous lncRNA in human for 17 common JAlncs (Table S2), of which nine were novel. We compared a secondary structure of the conserved JAlncs, because a lncRNA sequence is not generally well conserved. The predicted secondary structures of lncRNAs were compared using RNAforester (Höchsmann et al., 2003, 2004). We identified Gm14230 (Fig. 2D) as having the highest similarity to a human homolog in the secondary structure among the conserved JAlncs (Table S2). An open reading frame analysis with ORFfinder (https://www.ncbi.nlm.nih.gov/orffinder/) showed Gm14230 had a potential to code a peptide consisting of up to 70 amino acids. The gene encoding Gm14230 is adjacent to the Tgif2 gene in the mouse genome (Fig. 2E). In the human genome, the DLGAP4-AS1 lncRNA gene is adjacent to TGIF2 (Fig. 2E). The sequence of Gm14230 is similar to that of the DLGAP4-AS1 variant DLGAP4-AS1-205 (62.4% match, Fig. S1A). DLGAP4-AS1-205 is not considered an antisense transcript of DLGAP4 because they do not overlap, and DLGAP4-AS1-205 has exon–intron patterns that are distinct from those of other DLGAP4-AS1 variants (201, 202, 203 and 204) (Fig. 2E). Similar to Gm14230, DLGAP4-AS1-205 has not been functionally annotated.

The proximity of Gm14230 to Tgif2 suggests that these two genes are coordinately regulated. The coordinated regulation of a lncRNA and a protein-coding gene is a frequently observed phenomenon (Hon et al., 2017; Sigova et al., 2013). An analysis using the ChIP-Atlas (http://chip-atlas.org), which is a comprehensive database of previous ChIP-seq findings, revealed dense histone modifications associated with the active promoter (H3K4me3 and H3K9ac) upstream of Gm14230 in the first intron of Tgif2 in the mouse liver, heart and embryonic fibroblasts (Fig. S1B). Therefore, Gm14230 is likely regulated by a bidirectional promoter present in the first intron of Tgif2 that regulates both Gm14230 and Tgif2. Thus, we identified the common JAlncs and selected Gm14230 for further functional analysis.

Based on the hypothesis that a variant of Gm14230 might exist, we cloned the Gm14230 cDNA from NIH3T3 cells. Sequence analysis in Gm14230 revealed an adenosine (A) to guanosine (G) mutation in half (5 of 10) of the clones (Fig. S2A). This nucleotide alteration was also observed in the mouse neonatal hepatocytes (Fig. S2B). The recurrent substitution implied that Gm14230 is a substrate for adenosine to inosine (A-to-I) RNA editing (Gott and Emeson, 2000; Kawahara et al., 2007; Nishikura, 2010; Sommer et al., 1991; Zheng et al., 2016). Because inosine (I) base-pairs with cytosine (C), I is recognized as G by sequencing after cDNA synthesis. To investigate whether the nucleotide change is due to A-to-I editing, we performed cyanoethylation on the extracted RNA (Suzuki et al., 2015). Acrylonitrile treatment induces the cyanoethylation of inosine bases (ce¹I), and ce¹I does not base-pair with C and prevents cDNA synthesis. The acrylonitrile treatment of the RNA extracted from the NIH3T3 cells resulted in a decrease of the G peak (Fig. S2C), corroborating the idea that Gm14230 undergoes A-to-I editing. We further hypothesized the RNA editing of Gm14230 might have a functional relevance. We investigated an age-dependency of the RNA editing of Gm14230 and found the A-to-I editing was more frequent in cardiomyocytes derived from older mice (Fig. S2D). These results suggest that the function of Gm14230 is modified through RNA editing in the age-dependent manner.

The subcellular localization of lncRNAs has been systematically investigated in human cells (Mas-Ponte et al., 2017) but not in mouse cells. To investigate the subcellular localization of the Gm14230 lncRNA, we performed RNA-FISH and found Gm14230 localized both in the cytoplasm and nucleus (Fig. S3A). For a more quantitative analysis, we performed subcellular fractionation followed by qPCR analysis. This analysis showed that approximately two thirds of Gm14230 localized in the cytoplasm, and one third in the nucleus (Fig. S3B). Serving as control, primary (pri)-microRNAs, which are precursors of microRNAs serving as substrates for the so-called microprocessor complex in the nucleus (Gregory et al., 2004), were detected almost exclusively in the nuclear fraction, and Actb mRNA was detected mainly in the cytoplasmic fraction (Fig. S3B). We also tested the hypothesis that different cell densities could affect localization, but the Gm14230 subcellular localization was not significantly altered by the cell densities (Fig. S3B). Then, we examined whether cytoplasmic Gm14230 returns to the nucleus after cytoplasmic distribution following transcription. RNAi operating in the cytoplasm (Zeng and Cullen, 2002), if Gm14230 cannot translocate from the cytoplasm back to the nucleus, therefore RNAi cannot decrease the Gm14230 levels in the nucleus. The siRNA-mediated knockdown of Gm14230 leads to a decrease in Gm14230 in the nucleus, implying that the Gm14230 lncRNA shuttles between the cytoplasm and nucleus (Fig. S3C).

We further evaluated whether the fibroblast aging affected the subcellular localization of Gm14230. We induced cellular senescence by using the genotoxic agent zeocin (Fig. S3D; Mao et al., 2016; Robles and Adami, 1998). Gm14230 expression was repressed in the senescent cells, and the decrease was more significant in the cytoplasm than the nucleus, suggesting a cytoplasmic machinery that downregulates Gm14230 (Fig. S3E,F). These findings suggest that the amount of nuclear Gm14230 is regulated under the control of equilibrium that is affected by cellular aging.

To determine the biological relevance of Gm14230, we repressed Gm14230 by using siRNA in mouse NIH3T3 cells (Fig. 3A). In this experiment, we utilized three different siRNAs to distinguish any possible off-target effects. The Gm14230 knockdown using each of the three siRNAs significantly suppressed cell growth (Fig. 3B,C). The knockdown of a different JAlnc, AC153370.2, did not affect cell growth suggesting that the effect on Gm14230 is specific (Fig. S4A–C). Gm14230 knockdown decreased the number of Ki67-positive cells, suggesting that it induced cell cycle exit (Fig. S4D). In addition, Gm14230 knockdown did not increase the number of cells positive for active caspase 3 (Fig. S4E) or cleaved PARP (Fig. S4F), indicating that the growth suppression was not accompanied by apoptosis. We also noticed that the Gm14230 knockdown induced changes in the cell shape, such as a flattened and extended cytoplasm, reminiscent of senescent cells (Serrano et al., 1997). β-galactosidase staining, which is associated with senescence (Serrano et al., 1997), showed significantly increased positivity in cells transfected with siRNA against Gm14230 (Fig. 3D,E). Similar findings were observed with DLGAP4-AS1-205 in human ONS-76 cells (Fig. S4G–I). Furthermore, Gm14230 knockdown induced growth impairment (Fig. S4J–L) and cellular senescence (Fig. S4M,N) in Hepa1-6 cells, suggesting that Gm14230 has a functional relevance of in the hepatocytes in which Gm14230 was identified. Moreover, Gm14230 knockdown led to the induction of the senescence-associated proteins P16 (Alcorta et al., 1996), P27 (also known as CDKN1B) (Collado et al., 2000) and the marker γH2Ax (Fagagna et al., 2003) in NIH3T3 cells (Fig. 3F). Gm14230 knockdown also suppressed the protein expression of Srsf7, which is a juvenility-associated gene (JAG) that is highly expressed in juvenile cells (Jam et al., 2018) (Fig. 3F). Taken together, these observations indicate that the loss of Gm14230 leads to cellular growth impairment and cell senescence.

Subsequently, we sought to determine the mechanism by which Gm14230 regulates cell growth. Gm14230-adjacent Tgif2 belongs to a family of genes encoding TALE-homeodomain proteins, including Tgif1, Tgif2 and Tgif2lx1. Tgif2 functions as a transcriptional repressor and suppresses transcription initiated by TGF-β (Melhuish et al., 2001). Gm14230 knockdown enhanced Tgif2 mRNA expression in the NIH3T3 cells, suggesting that Gm14230 functions as a negative regulator of Tgif2 expression (Fig. 4A). To assess the regulation of Tgif2 by Gm14230 further, we examined an expression of Evi5l, a transcriptional target of Tgif2 (Anderson et al., 2017). Gm14230 knockdown significantly downregulated Evi5l, while the effect was abolished by simultaneous knockdown of Gm14230 and Tgif2 (Fig. 4A). This observation indicates that Gm14230 regulates Evi5l through Tgif2.

Next, we aimed to identify the partner required for Gm14230 to repress Tgif2. lncRNAs, such as Xist (Brockdorff et al., 1991; Zhao et al., 2008), Hotair (Rinn et al., 2007), Bvht (Klattenhoff et al., 2013), Pandar (Puvvula et al., 2014), Fendrr (Grote et al., 2013) and Airn (Braidotti et al., 2004; Nagano et al., 2008), interact with the PRC complex and induce gene silencing. Thus, the interactions between Gm14230 and Ezh2, which is a component of the PRC2 complex, were assessed. RNA immunoprecipitation and subsequent qPCR (RIP-qPCR) analyses revealed that Gm14230 bound the Ezh2 protein in NIH3T3 cells (Fig. 4B). Furthermore, the re-analysis of the RNA immunoprecipitation-sequencing (RIP-seq) that comprehensively identified the Ezh2-bound transcriptome (Zhao et al., 2010), revealed interactions between Gm14230 and Ezh2 in mouse embryonic stem cells (Fig. S5A).

To directly assess the role of Gm14230 in Ezh2 recruitment to Tgif2 locus, we performed chromatin immunoprecipitation (ChIP) analysis for Ezh2 protein in normal and Gm14230-knockdown NIH3T3 cells. Among the evolutionarily conserved sequence near the Tgif2 transcription start site (Fig. 4C), two loci in the Tgif2 first intron showed significant enrichment when immunoprecipitated by Ezh2 (Fig. 4D). siRNA-mediated depletion of Gm14230 abrogated the occupation of Tgif2 loci by Ezh2, indicating the pivotal role of Gm14230 for the recruitment of Ezh2 at these loci (Fig. 4D). H1foo, which is expressed in an oocyte-specific manner (Teranishi et al., 2004), showed Ezh2 occupation regardless of Gm14230 expression (Fig. 4D). Thus, Gm14230 recruits Ezh2 to the Tgif2 regulatory elements, thereby repressing Tgif2 transcription (Fig. 4E).

We then investigated the biological relevance of Tgif2 expression on cell growth. Compared to the EGFP control, overexpression of Tgif2 significantly suppressed the growth of NIH3T3 cells (Fig. 5A,B). The growth suppression was accompanied by decreased levels of Ki67 (Fig. S5B) but no active caspase 3 (Fig. S5C) or cleaved PARP (Fig. S5D) was detected, excluding the involvement of apoptosis. Moreover, the Tgif2-transfected cells exhibited an extended cytoplasm and positive staining for senescence-associated β-galactosidase (SA-β-gal) activity, indicating that enhanced Tgif2 expression leads to cell senescence (Fig. 5C,D). The Tgif2-transfected cells expressed the senescence-associated proteins P16, P27 and γH2Ax (Fig. 5E).

To investigate the role of Gm14230–Tgif2 axis in senescence induction, we performed simultaneous knockdown of Gm14230 and Tgif2. A significant rescue of the cell growth impairment that was seen with Gm14230 knockdown alone was observed (Fig. 5F,G), thus confirming the function of Tgif2 in combination with Gm14230. In addition, the simultaneous knockdown decreased the positive staining for SA-β-gal activity induced by Gm14230 knockdown (Fig. 5H,I), supporting the notion that Gm14230 regulates the cellular juvenescence by repressing Tgif2. Knockdown of Evi5l, a downstream target of Gm14230-Tgif2 axis, neither impaired cell growth (Fig. S5E–G) nor induced cell senescence (Fig. S5H), suggesting that Evi5l is not used by Gm14230–Tgif2 axis to regulate the cell senescence. Collectively, Gm14230 regulates cell growth and titrates Tgif2 expression to avoid undergoing cell senescence (Fig. 6). Taken together, the analysis of the observations from this study has revealed the essential roles of JAlncs, which are lncRNAs that are highly expressed selectively in juvenile cells. Moreover, Gm14230 is an indispensable factor that regulates cell growth and cellular juvenescence.

In this study, we performed a comprehensive lncRNA analysis in differentially aged cardiomyocytes and hepatocytes to identify lncRNAs that are highly expressed selectively in juvenile cells. Among the juvenility-associated lncRNAs (JAlncs), we identified Gm14230 as an essential factor for the maintenance of cellular juvenescence. The cellular juvenile properties, including the higher proliferation, differentiation and maturation abilities, have a tendency to be lost during cell senescence. Molecules involved in cellular juvenescence might increase the ability of cells to proliferate and regenerate, and, hence, counteract cellular senescence. As demonstrated by the depletion of the JAlnc Gm14230 causing premature senescence in the various cell types, a role of the JAlncs is to maintain the juvenile state to counteract senescence.

Since the discovery that the genome is pervasively transcribed, numerous lncRNAs have been functionally investigated. Among JAlncs, Lockd is an enhancer RNA transcribed from an enhancer of the Cdkn1A gene (Paralkar et al., 2016), Snhg5 promotes cancer cell proliferation (Damas et al., 2016) and Crnde is associated with cancer progression (Ellis et al., 2012). Adult-expressed lncRNAs include Lncbate1, which is associated with mitochondrial biogenesis (Alvarez-Dominguez et al., 2015).

lncRNAs have the potential to serve as therapeutic targets; specifically targeting a lncRNA via an oligonucleotide, such as siRNA or antisense oligonucleotide (ASO), is possible and has been clinically applied for the treatment of spinal muscular atrophy (Lim and Hertel, 2001) and Duchenne muscular dystrophy (van Deutekom, 2001). Whenever JAlncs might be associated with intractable human diseases, they might as well become a therapeutic target via an oligonucleotide.

Gm14230 expression is age dependent and is virtually non-existent in adult cells in mouse hepatocytes and cardiomyocytes. As shown by the ChIP-atlas analysis of histone modifications, Gm14230 transcription can be regulated by the promoter present in the first intron of Tgif2. This regulatory element is also likely to regulate the transcription of Tgif2 because a region with a similar histone modification does not exist around the Tgif2 transcription start site (TSS). Divergent or bidirectional promoters are frequently observed to regulate lncRNA and protein-coding gene pairs (Hon et al., 2017; Sigova et al., 2013). Elucidating the mechanism by which the divergent promoter of Gm14230 and Tgif2 is silenced in an age-dependent manner may reveal the mechanistic machinery involved in the loss of juvenility and initiation of aging.

We found that the subcellular localization of Gm14230 is maintained at an equilibrium between the nucleus and the cytoplasm (Fig. S3). The perturbation of the equilibrium occurred in the premature senescent cells, resulting in a lower Gm14230 expression in the cytoplasm, implying that a cytoplasmic machinery downregulates Gm14230. One aspect of the cytoplasmic machinery might involve microRNAs. lncRNAs are targets for miRNA-mediated degradation, as in the case of Hottip lncRNA, which is targeted by miR-192 and 204 (Ge et al., 2015). Gm14230 contains multiple potential target sites for miRNAs (Fig. S6A), and the expression of Gm14230 is enhanced in Dicer-deficient mouse embryonic stem cells (Zheng et al., 2014), Dicer-deficient spermatocytes (Modzelewski et al., 2015) and Ago2-deficient embryonic stem cells (Ngondo et al., 2018; Fig. S6B). These findings implied the cytoplasmic expression of Gm14230 might be regulated by microRNAs.

Our analysis also discovered Gm14230 is a substrate for the A-to-I RNA editing. A biological role for the RNA editing of lncRNAs is not fully understood. An alteration of the sequence is considered to change an affinity for binding with other molecular species such as DNA (Batista and Chang, 2013), RNA (Jalali et al., 2013), protein (Ferrè et al., 2016) and others, thereby influencing the downstream molecular function of a lncRNA. Further elucidation of the effects of RNA editing might reveal an as yet unknown role for the RNA editing of lncRNAs.

The Gm14230 knockdown led to enhanced Tgif2, indicating that the Gm14230 lncRNA represses Tgif2 mRNA expression. The repression of a gene by a lncRNA is also observed with other lncRNAs, such as Xist (Brockdorff et al., 1991; Zhao et al., 2008), Hotair (Rinn et al., 2007), Bvht (Klattenhoff et al., 2013), Pandar (Puvvula et al., 2014), Fendrr (Grote et al., 2013) and Airn (Braidotti et al., 2004; Nagano et al., 2008). We showed that Gm14230 interacted with Ezh2, which is a subunit of the PRC2 chromatin remodeling complex. Thus, Gm14230 likely recruits PRC2 to repress Tgif2 expression. Tgif2 is a homeodomain-containing transcriptional repressor that suppresses Evi5l transcription. As shown by the knockdown of Gm14230 and Tgif2, Gm14230 upregulates Evi5l through Tgif2. Gm14230, Tgif2 and Evi5l are expressed more highly in juvenile than in adult and Gm14230 titrates the expression of Tgif2 and Evi5l. In the absence of Gm14230, Tgif2 is transcribed in excess, leading to cell cycle arrest and cellular senescence.

The following question remains: how does Gm14230 maintain cellular juvenescence? Gm14230 and Tgif2 are co-regulated by an element in the first intron of Tgif2. Tgif2 suppresses cellular growth and eventually induces cell senescence. Tgif2 overexpression is observed in ovarian cancer cell lines (Imoto et al., 2000), suggesting that Tgif2 may function as an oncogene. Similar oncogene-induced senescence has been observed in genes, such as Ras (Serrano et al., 1997), BRAF (Wajapeyee et al., 2008), Akt (Chen et al., 2005) and E2F (Denchi et al., 2005). Gm14230 is believed to titrate the magnitude of Tgif2 expression. Gm14230 knockdown also promotes excess Tgif2 expression, leading to cell senescence. Thus, Gm14230 represses Tgif2 expression to maintain cellular juvenescence and prevent premature senescence.

In this study, we identified JAlncs, which are lncRNAs that are selectively highly expressed in juvenile hepatocytes and cardiomyocytes in mice. Among the 62 common JAlncs, Gm14230 is evolutionarily conserved and is essential for the maintenance of cellular juvenescence. The loss of Gm14230 leads to cell growth impairment and cell senescence. Therefore, our findings reveal the previously unrecognized role of lncRNAs in the transcriptional network that governs the physiology of juvenile cells.

In the previous study (Jam et al., 2018), we performed RNA-seq analysis with the hepatocytes and the cardiomyocytes isolated from the C57BL/6N male mice aged at postnatal day 1 (P1), 7 (P7) and 56 (P56) in triplicate. All animal experiments were approved by the institutional animal care and use committee at Tokyo medical and dental university. All experiments were performed in accordance with the relevant guidelines and regulations. We analyzed the RNA-seq data as below.

The lncRNA analysis was performed using our previously reported dataset, which has been deposited to DNA Data Bank of Japan (DDBJ) under the accession number DRA007101 (https://ddbj.nig.ac.jp/DRASearch/submission?acc=DRA007101). The analyses were executed using BioLinux8 (Field et al., 2006). The Fastq files were processed with Fastx trimmer (FASTX-Toolkit, http://hannonlab.cshl.edu/fastx_toolkit/) to remove the first base from all the reads because the first base that was nearest to the annealing site of the sequence primer showed a frequent sequencing error. From the 3′-terminus, the characters with a quality score lower than 20 were trimmed using the Fastq quality trimmer. The reads were removed if their length was shorter than 30 characters after trimming. Reads with more than 80% of characters harboring a quality score higher than 20 were extracted for further analysis. Additionally, reads were removed using the Fastq quality filter if 10% or more of the characters had a quality score less than 30. The Fastq files were processed as described above and mapped to GRCm38 using Tophat2 (Kim et al., 2013). The GTF file of the lncRNA analysis was generated from the Ensembl GTF files by extracting the genes with a gene or transcript biotype of ‘lincRNA’, ‘bidirectional promoter lncRNA’ or ‘macro lncRNA’. The differential gene expression was analyzed using Cuffdiff (Trapnell et al., 2010) and the lncRNA GTF file.

To identify the juvenility-associated lncRNAs (JAlncs), we processed the fragments per kilobase of transcript per million fragments sequenced (FPKM) data as follows. First, lncRNAs that showed 0 FPKM for all samples were removed. Next, those lncRNAs that were detected in one sample or less among the P1 triplicates were removed. The lncRNAs with fold changes of greater than 2 at both P1 and P7 compared to those at P56 were included. Some lncRNAs had specific expression in juvenile phases and expression in P56 was not detected at all (FPKM=0). In those cases, a fold change is impossible to calculate. To include the lncRNAs that were expressed in juvenile phases but not at all in the adult phase, those lncRNAs were included if two or three of the three samples for both P1 and P7 showed more than 0 FPKM.

Heatmaps were generated using MeV (http://mev.tm4.org). The expression was visualized with the Integrative Genomics Viewer (http://software.broadinstitute.org/software/igv/). Open reading frame analysis was performed with ORFfinder (https://www.ncbi.nlm.nih.gov/orffinder/). The evolutionary conservation between Gm14230 and DLGAP4-AS1-205 was assessed using ClustalW (Thompson et al., 1994). To investigate the regulatory elements associated with gene expression, the ChIP-atlas database was utilized (http://chip-atlas.org). microRNA target sites in Gm14230 were searched using LncBase Predicted v.2 (http://carolina.imis.athena-innovation.gr/diana_tools/web/index.php) and RegRNA2.0 (http://regrna2.mbc.nctu.edu.tw/). For the quantification of Gm14230 in microRNA-deficient conditions, we reanalyzed the deposited data of RNA-seq performed in Dicer-knockout (KO) mouse embryonic stem (ES) cells (GSE55338), Dicer-KO spermatocytes (GSE63166) and Ago2-KO ES cells (GSE80454). We downloaded the fastq files from the Sequence Read Archive (SRA) and calculated the expression levels as transcripts per million (TPM) using Salmon (Patro et al., 2017).

The interspecies conservation of lncRNA between mouse and human were analyzed as below. Sequences of the mouse common JAlncs were obtained from the Ensembl database. Sequences of the human lncRNAs were obtained from Ensembl BioMart (https://www.ensembl.org/biomart). A database of human lncRNAs were generated using the makeblastdb command of BLAST (version 2.7.1). The sequence similarity between mouse and human lncRNAs was assessed by running blastn, where a mouse JAlnc sequence was queried in the generated human lncRNA database. To quantitatively compare the extent of conservation, E-values and bit scores were used. The interspecies similarity in the secondary structure of lncRNAs were assessed using RNAforester, a tool that compares RNA secondary structures via forest alignment (Höchsmann et al., 2003, 2004, https://bibiserv.cebitec.uni-bielefeld.de/rnaforester). RNAforester calculates similarity between a pair of RNA structures and calculates the relative scores that are upper-bounded by 1, which is the score for equal structures.

The Hepa1-6 cell line was purchased from The European Collection of Authenticated Cell Cultures (ECACC). The NIH3T3 and ONS-76 cell lines were purchased from Japanese Collection of Research Bioresources (JCRB) Cell Bank. NIH3T3 and ONS-76 cells were used in this study because these cell lines exhibited hallmarks of the cellular senescence prominently among tested cell lines. Hepa1-6 and NIH3T3 cell lines were cultured in Dulbecco's modified Eagle's medium (DMEM) containing 10% fetal bovine serum (FBS) and penicillin/streptomycin. ONS-76 cells were cultured in RPMI containing 10% fetal bovine serum (FBS) and penicillin/streptomycin. In the knockdown experiments, Lipofectamine RNAiMAX (Invitrogen) was used to transfect the siRNA duplexes at 10 nM: Gm14230-siRNA1 (5′-CGAAGACCACUUGGGAUAAdTdT-3′), Gm14230-siRNA2 (5′-GGUGCUAAAGGACCAGUUGdTdT-3′), Gm14230-siRNA3 (5′-GCUCAGUCGGCUCAGAGUAdTdT-3′), AC153370.2 (5′-CUGCUGACGUCUGCCGGCAdTdT-3′), DLGAP4-AS1-205 (5′-CGGUGCGAUUGGUCCCGUUdTdT-3′), Evi5l (Sigma-Aldrich, Mission siRNA SASI_Mm01_00020067) and siRNA negative control (Applied Biosystems, AM4611). The siRNA negative control was designed to have no significant sequence similarity with mouse or human transcript sequences, and used for the purpose of establishing a baseline and eliminating the systemic effects produced by transfection.

pCAGIP-EGFP has been previously described (Mori et al., 2009). pcDNA4-TO-HA-Tgif2 was generated by PCR cloning human Tgif2. For the plasmid transfection, Lipofectamine 2000 (Invitrogen) was used according to the manufacturer's instructions. The cell images were obtained using EVOS FL (Thermo Fisher Scientific). The number of cells per field was counted manually.

For the staining of SA-β-Gal activity, 7 days after the transfection, the cells were fixed with 3% (v/v) formaldehyde in PBS, washed with PBS twice, and incubated with the staining solution [5 mM K₃Fe(CN)₆, 2 mM MgCl₂, 150 mM NaCl, 30 mM citric acid/phosphate buffer, 5 mM K₄Fe(CN)₆, and 1 mg/ml X-Gal] for 16 h at 37°C.

For the zeocin-induced cellular senescence model, cells were treated with zeocin (100 μg/ml) for 48 h and were continued to be cultured in fresh medium (Robles and Adami, 1998; Mao et al., 2016).

The total RNA was extracted using TRIzol reagent (Thermo Fisher Scientific). The extracted RNA was quantified using a NanoDrop Lite Spectrophotometer (Thermo Fisher Scientific) and reverse transcribed using a high-capacity RNA-to-cDNA kit (Applied Biosystems) and a PCR thermal cycler (Dice; Takara) according to the manufacturer's instructions. qPCR was performed using a LightCycler 480 SYBR Green I Master Kit on a LightCycler 480 instrument (Roche) using the reverse transcribed cDNA as a template. The specificity and quality of the qPCR amplification was assessed by performing a melting curve analysis. The data were normalized to mouse Polr2a or Tubb5 as indicated in the figure legends. The primers used for the qPCR are listed in Table S3.

The sequence of Gm14230 in the NIH3T3 cells and mouse neonatal hepatocytes (at P1) was analyzed by cloning the Gm14230 cDNA. The PCR products were cloned into a T-vector using a Mighty TA-cloning Kit (Takara). The sequences were analyzed for each of the colonies. The RNA editing analysis was performed using cyanoethylation as previously described (Suzuki et al., 2015). Briefly, 4 µl of the RNA sample adjusted to 0.125 µg/µl was mixed with 30 µl of the CE solution (50% ethanol and 1.1 M triethylammonium acetate, pH 8.6) and 4 µl of 15.2 M acrylonitrile and incubated at 70°C for 60 min. The reaction was quenched by the addition of 162 µl of RNase-free water on ice. The RNA was extracted by ethanol precipitation and reverse transcribed as described in the ‘Gene expression analysis’ section. The sequences of the synthesized cDNAs were analyzed by direct sequencing.

For the subcellular fractionation of the cultured cells, the cells were trypsinized, pelleted by centrifugation (300 g for 2 min), resuspended in lysis buffer (50 mM Tris-HCl, pH 8.0, 140 mM NaCl, 1.5 mM MgCl₂, and 0.5% NP-40) and incubated at 4°C for 5 min. Then, the lysates were centrifuged at 300 g for 2 min. The supernatants were collected as the cytoplasmic fractions. To wash the nuclear fraction, chilled PBS was added to the pellet, followed by centrifugation at 300 g for 2 min. The supernatants were discarded, and the pellet remained as the nuclear fraction. For the RNA extraction, 1 ml of TRIzol reagent (Thermo Fisher Scientific) was added to both fractions. The RNA was extracted as described in the ‘Gene expression analysis’ section.

Fluorescent in situ hybridization for Gm14230 in NIH3T3 cells was performed as described previously (Sone et al., 2007). More specifically, Gm14230 cDNA was cloned and PCR-amplified using the forward primer: 5′-GAGTGTTTCCGGAGCGCCTAC-3′ and the T7 site-adding reverse primer: 5′-TTAATACGACTCACTATAGG TCACCCAGCATTCAGTCCGC-3′. A repetitive sequence was searched by RepeatMasker (www.repeatmasker.org/). The PCR product was used as a template for RNA synthesis. Digoxigenin-labeled RNA probes were prepared by using RNA-labeling mixture (Roche) and T7 RNA polymerase (Roche). NIH3T3 cells were fixed with 4% paraformaldehyde for 10 min at room temperature (RT) followed by a PBS wash. Cells were then permeabilized with PBS containing 0.5% Tween 20 for 5 min. After incubation in prehybridization buffer (50% formamide, 1× Denhardt's solution, 2× SSC, 10 mM EDTA, 100 μg/ml yeast tRNA, 0.01% Tween 20) for 2 h at 55°C, cells were incubated with synthesized RNA probes diluted in hybridization buffer (50% formamide, 1× Denhardt's solution, 2× SSC, 10 mM EDTA, 100 μg/ml yeast tRNA, 0.01% Tween 20, 5% dextran sulfate) overnight at 55°C. The hybridized probes were detected with standard immunohistochemical procedures using a mouse monoclonal anti-DIG primary antibody (1:200, 21H8, #ab420, Abcam) and Alexa Fluor 488-conjugated anti-mouse secondary antibody (ThermoFisher Scientific). Images were acquired using Nikon C1si.

The cells were harvested in RIPA lysis buffer containing 25 mM Tris-HCl (pH 7.6), 150 mM NaCl, 1% NP-40, 1% sodium deoxycholate and 0.1% SDS. The protein samples (10 µg per lane) were separated by SDS-PAGE and transferred to a Hybond-P PVDF membrane (GE Healthcare). Western blotting was performed using the antibodies against the following proteins: (1) P16 ARC (1:1000, Abcam, ab51243), (2) P27 Kip1 (1:1000, Cell Signaling, #3698), (3) γH2Ax (1:1000, Cell Signaling, #9718), (4) Srsf7 (1:1000, Thermo Fisher Scientific, PA5-39482) and (5) β-actin (1:2000, Cell Signaling, #5125). Secondary antibodies conjugated with horseradish peroxidase (Thermo Fisher Scientific, anti-mouse-IgG, 32430 and anti-rabbit-IgG 32460) were used at a 1:1000 dilution. The immunoreactive bands were detected using Chemi-Lumi One L or Chemi-Lumi One Ultra (Nacalai Tesque). Uncropped blots were provided in Fig. S6C.

NIH3T3 cells transfected with siRNA or a plasmid vector were fixed with 4% paraformaldehyde for 10 min at RT, permeabilized with 0.1% Triton X-100 for 2 min at RT, blocked with 2% FBS, and incubated with antibody against Ki67 (1:400, Cell Signaling, #12202), active caspase 3 (1:200, Abcam, ab2302), cleaved PARP (1:200, Cell Signaling, #9544) and vimentin (1:80, Santa Cruz Biotechnology, sc-6260) at 4°C overnight. After washing PBS, cells were incubated with anti-mouse-IgG conjugated to Alexa Fluor 488 and anti-rabbit-IgG conjugated to Alexa Fluor 546 (1:1000, Invitrogen). Nuclei were stained with DAPI (Invitrogen).

Crosslinking was performed by adding 3.3% formaldehyde to the medium, followed by incubation at RT for 10 min. The reaction was stopped by the addition of 6.5 μM glycine and further incubated for 5 min at RT. After a few washes for several minutes each with cold PBS, the cells were lysed with NETN buffer [20 mM Tris-Cl, pH 8.0, 1 mM EDTA, 100 mM NaCl, and 0.5% (v/v) Nonidet P-40] and collected by scraping. The cell lysates were treated with 10 units of Turbo DNase at 37°C for 10 min. After centrifugation at 20,000 g for 5 min at 4°C, a part of the supernatant was collected as the input sample. The remaining supernatant was incubated with the anti-Ezh2 antibody (Cell Signaling, 5246S) or control rabbit IgG (Santa Cruz Biotechnology) conjugated to Dynabeads Protein A (10002D, Life Technologies) overnight at 4°C. The Dynabeads complexes were collected using a magnet stand (Dynal, Thermo Fisher Scientific) and washed with 170 µl CLIP wash buffer (1× PBS without Mg²⁺ or Ca²⁺, 0.1% SDS, and 0.5% NP-40) twice for 10 min and then with 170 µl high-salt wash buffer (5× PBS without Mg²⁺ or Ca²⁺, 0.1% SDS, and 0.5% NP-40) twice for 10 min. Then, the beads were resuspended with RIPA buffer (25 mM Tris-HCl, pH 7.6, 150 mM NaCl, 1% NP-40, 1% sodium deoxycholate and 0.1% SDS). The crosslinks were reversed by incubating the beads with proteinase K at 70°C for 3 h. The supernatants were collected as the RIP samples, and the RNA was extracted using TRIzol as described in the ‘Gene expression analysis’ section. The qPCR quantification was performed based on standard curves. The relative enrichment by the Ezh2 antibody was calculated by normalizing to the amounts of RNA pulled down by the control IgG.

Chromatin immunoprecipitation was performed as previously described (Nimura et al., 2009) with NIH3T3 cells transfected with Gm14230 siRNA_2 or negative control siRNA at 20 nM. Primer sequences are listed in Table S4. The data and the image for the interspecies conservation at the mouse Tgif2 locus was obtained from the UCSC genome browser (https://genome.ucsc.edu/).

The mean±s.d. is presented for all quantified data. The statistical significance between two experimental groups is indicated by an asterisk, and the comparisons were performed using Student's t-tests. P-values less than 0.05 were considered significant.

We would like to thank all the laboratory members from the Molecular Neuroscience Research Center (MNRC) for their helpful discussions and sincere cooperation. This study was also supported by the Central Research Laboratory at Shiga University of Medical Science (SUMS).