Senescence induction universally activates transposable element expression

ABSTRACT

Senescent cells constitutively secrete inflammatory cytokines, known as the senescence-associated secretory phenotype (SASP). Previous work has implicated SASP in immune-mediated clearance of senescent cells; however, its regulation remains unknown. Our recent transcriptome profiling study has shown that human senescent human stem and progenitors (s-HSPCs) robustly express genomic transposable elements (TEs) and pathways of inflammation. Furthermore, hypomethylating agents have been previously shown to induce expression of TEs and activate the dsRNA recognition pathway and downstream interferon-stimulated genes, leading to immune mediated cell death. Therefore, to examine whether activation of TEs occurred universally, independent of their modality of senescence induction, we performed transcriptomic analysis in artificially-induced senescent cell-lines and observed a robust activation of TEs. Hence we propose that the expression of TEs might play a role in immune mediated clearance of senescent cells.

KEYWORDS:

• Senescence
• transposable elements
• ERV
• dsRNA
• inflammation
• ISG
• SASP

Introduction

A complex cellular process contributes to the aging phenotype, characterized by progressive loss of organ function [1–3]. Cellular senescence, defined as a state of permanent cell cycle arrest, plays an important role in aging [4–6]. Stressors such as telomere shortening, mitochondrial dysfunction, oxidative stress, and DNA damage induce senescence [7]. Senescent cells are known to be actively cleared by the immune system [8–11]. Furthermore, they constitutively secrete cytokines to generate a SASP (senescence associated secretory phenotype) wherein, components of SASP in turn, enable the clearance of damaged senescent cells (“senescence surveillance”) [12], indicative of an important link between immunity, senescence clearance and longevity. In addition, SASP is known to exhibit several pleiotropic effects in tumorigenesis [13] – dependent on a cell growth arrest or pro-oncogenic milieu- based secretome [14–16]. In this regard, although the role of SASP and its implication in senescence surveillance are well studied, the regulation of SASP – specifically its upstream activators, are yet to be identified.

Approximately 50% of the genome comprises of TEs [17]. In steady state, DNA methylation, chromatin modification, RNA editing and other mechanisms repress their activation [18]. Recent works demonstrating derepression of TEs with hypomethylating agents (HMA), resulted in activation of the ds-RNA recognition pathway, leading to tumor clearance by immune cells (“viral mimicry”)[19,20]. Furthermore, we have shown that human s-HSPCs in vivo show activated expression of TEs along with immune pathways [21]. Thus, it is plausible that the expression of TEs in senescent cells foretell the immune system for their orderly clearance in order to maintain homeostasis.

Recent studies investigating the functional consequences of senescent cell clearance, demonstrated: rejuvenation of aging hematopoietic and muscle stem cells [22], improvement in renal function and frailty [23] and an increase in life span [24–26]. Similar growing interests in the development of therapeutics targeting senescent cells (senolytics) to improve organ function, underscores the importance to further our understanding of the molecular determinants/underpinnings/basis of senescence and its homoeostasis.

We therefore hypothesized that senescence induction led to activation of TEs and immunity related genes. Furthermore, we queried whether various modalities of senescence induction activated similar TE profiles. For this we characterized the transcriptome of cultured cells [27] induced to senescence states by higher cell-passages and chemotherapy agent, oxidative stress, HMA.

Materials and methods

The RNA-seq experimental data was obtained from Purcell et. al [27]., and the published dataset were obtained from GEO (Accession No.: GSE60340). RNA-seq and gene expression array data for oncogene induced senescence (OIS) was obtained from Pawlikowski et. al [28]. Transcript abundance was quantified by pseudoalignment using Kallisto [29]. Sequenced reads were mapped to a reference transcriptome prepared using Arkas [30], comprised of ENSEMBL human (GRCh38-build 81) coding and non-coding transcripts, repetitive elements (RepBase issue 21-03). Arkas managed the preparation of indices, quantification of transcript abundance, and annotation of Kallisto results by leveraging S4-summarized experiment structure.

TE transcript expression for each experimental group using upper-quartile normalization was used for principal component analysis (PCA) [31]. The relative log-expression of all TE transcript expression relative to each corresponding experimental group was depicted in a bar-plot. The original manuscript [27] labeled the lowest low passage samples as ‘LLP’, and in this manuscript please note that the ‘LLP’ experimental group was labeled as ‘Low Passage’.

TE expression counts less than 2 across samples were filtered from differential expression analysis. Differential expression used upper-quartile normalized counts comparing expression between treatment and control determined using ‘edgeR’[32]. The linear model contrasted Adriamycin, 5-Azacitidine, Oxidative (H202 treated), and natural, replicative, senescence against lowest low passage (original nomenclature ‘LLP’). The statistically significant differences in normalized log2 (CPM) (likelihood ratio test, adjusted p-value <0.05, FDR adjustment method: Benjamini-Hochberg procedure) were reported as differentially expressed. The statistically significantly dysregulated TE were subset, and their log₂ transcript-per-million (TPM) values were depicted in each figure. The differential expressed (FDR < 0.05, adjustment method: Benjamini-Hochberg) TE transcript biotypes were tabulated, and the linear model coefficients as log-fold change were used to determine the estimated dysregulation. The expected TE expression was estimated by compounding the linear model coefficients with the corresponding significant differential expression of TE (DE TEs) normalized observed expression values.

The negative binomial GLM (‘edgeR’) was constructed using a design matrix, which identified all the experimental groups. For each experimental group, the fitted negative binomial GLM was contrasted and tested using the likelihood ratio test (LRT), and differentially expressed TEs were identified using FDR filtering alpha threshold of 0.05. After each LRT identified significantly DE TEs per experimental group, a set-analysis was performed. The Venn diagram identified 33 common significantly DE TEs activated by the 4 different induced senescent experiments.

The viral recognition pathway gene set analysis used canonical interferon stimulated genes (ISGs) involved with viral pathogenesis [33]. The enrichment of the ISGs were calculated using QuSAGE [34] and a FDR cutoff of 0.05 using Benjamini-Hochberg adjusted p-values. The QuSAGE estimated 95% confidence intervals were used to identify genes that had significant pathway enrichment activity with a p-value alpha threshold of 0.05.

Results

Induction of senescence activates TE expression

Senescence can be induced by various stimuli including genotoxic agents and oxidative, replicative, and oncogenic stresses [5]. We have previously reported that primary senescent hematopoietic progenitors exhibited elevated expression of TEs [21]. In order to examine if the activation of TE expression was a universal consequence of senescence induction, we analyzed transcriptome data generated from fibroblasts with Li-Fraumeni syndrome (LFS fibroblasts; heterozygous for p53), in which senescence was induced through replicative stress, 5-azacytidine, oxidative stress, or Adriamycin [27]. This data-set was previously published, analyzing only the coding regions, which revealed activation of various inflammatory genes, including those involved in SASP.

Prior studies have shown that HMA induces TE expression in cell lines [19,20,35]. We analyzed differential expression of TEs in LFS fibroblasts with and without HMA treatment and observed significant upregulation of TEs in HMA treated cells (Figure 1, Supplement data 1); thereby confirming findings from previous studies implicating the role of HMA in activating TE expression, as well as support the possibility of an interplay between senescence and TE- this could prove beneficial in cancer therapeutics.

Figure 1. Differential expression of transposable elements following 5-Azacitidine treatment. (a) Differential expression of TE following senescence induced by 5-Azacitidine. The x-axis is the experimental groups, and the y-axis are the significantly dysregulated TE. (b) Differentially expressed TE biotype dysregulation following 5-Azacitidine treatment. The significantly DE TEs biotypes were identified and the modeled expression estimates are summarized for each TE transcript biotype. (c) 5-Azacitidine expected TE expression group summary. The two groups expected expression was determined by compounding the dysregulation coefficients of significantly DE TEs with the observed normalized corresponding DETE expression. The y-axis denotes the expected expression predicted by the linear model constructed in A.

We next analyzed the transcriptome of LFS fibroblast cells treated with other inducers of senescence and consistently observed a similar pattern of significant upregulation of TEs, in all models of senescence, compared to actively dividing cells (Figure 2(a) and 2(b), Supplement Figure 1 and Supplement data 1).

Figure 2. Transposable element overall expression of senescence induction groups. (a) Principal component analysis (PCA) of TE expression following treatment of 5-Azacitidine, Replicative senescence, oxidative (H202) treatment, and Adriamycin. (b) Relative log-expression of TE expression following treatments that induce senescence.

Replicative senescence significantly dysregulated 79 TEs, of which 72 (~91.1%) were up-regulated (Figure 3(a), Supplement Data 2) and 58 (~80.6%) had a minimum log fold change of at least 1. The predominant TE biotypes significantly differentially expresses were ERV1 (30 dysregulated, 27 up-regulated, 3 down-regulated) and L1 with 18 up-regulated (Figure 3(b)) (Supplement Data 2).

Figure 3. Differential expression of transposable elements following replicative senescence. (a) Differential expression of TE following senescence induced by replicative cell passages. Labeling is similar to Figure 2(a). (b) Differentially expressed TE biotype dysregulation following replicative senescence. Labeling is similar to Figure 2B. (c) Replicative senescent expected TE expression group summary. Labeling is similar to Figure 2(c).

5-Azacitidine treated LFS had the most differentially expressed TEs (DE TEs), with 169 significantly differentially expressed repeats (q-value <0.05, adjusted method: BH), of which 163 (~96.5%) transcripts were upregulated (Figure 2(a)) (Supplemental Data 1). The predominant dysregulated TE biotypes were ERV1 (85 TEs dysregulated, 82 up-regulated, 3 down-regulated), and L1 with 27 TEs significantly up-regulated (Figure 2(b)).

Adriamycin induced senescence and activated 114 TE transcripts, 107 (93.9%) were up-regulated (Figure 4(a), Supplement Data 3). The 114 DE TEs consisted of predominantly ERV1 (49 dysregulated, 45 up-regulated, 4 down-regulated) and L1 (27 dysregulated, 26 up-regulated, and 1 down-regulated) (Figure 4(b)).

Figure 4. Differential expression of transposable elements following Adriamycin treatment. (a) Differential expression of TE following senescence induced by Adriamycin. Labeling is similar to Figure 2(a). (b) Differentially expressed TE biotype dysregulation following Adriamycin chemotherapy induced senescence. Labeling is similar to Figure 2(b). (c) Chemotherapy induced senescent expected TE expression group summary. Labeling is similar to Figure 2(c).

Oxidative stress induced senescence activated 121 TEs, 114 (~94.2%) (Figure 5(a), Supplement Data 4). The predominantly dysregulated TE biotypes included ERV1 (50 dysregulated, 47 up-regulated, 3 down-regulated) and L1 (33 up-regulated) (Figure 5(b)).

Figure 5. Differential expression of transposable elements following senescence induction by oxidative stress (H202 treatment). (a) Differential expression of TE following senescence induced by oxidative stress. Labeling is similar to Figure 2(a). (b) Differentially expressed TE biotype dysregulation following oxidative stress induced senescence. Labeling is similar to Figure 2(b). (c) Stress induced senescent expected TE expression group summary. Labeling is similar to Figure 2(c).

In addition to ERV1 and L1, which were the predominant up-regulated biotypes in all models of senescence induction, ERVK and ERVL were also up-regulated; although to a lesser degree. We also identified 33 differentially expressed TEs that were commonly dysregulated in all 4 experimental models (Figure 6(a)) and they consistently comprised of ERV1 and L1 (Figure 6(b)). Collectively, these data suggest that, this model system of senescence induction leads to the upregulation of TEs.

Figure 6. Senescence induced differentially expressed TE transcript per each experimental modality. (a) Set analysis of the DETE identified from the 4 modeled groups. The overlaps in the Venn Diagram indicate the number of DETE significantly dysregulated in the corresponding experimental group. (b) The differential expression of the 33 commonly TE dysregulated across all senescent modalities. The x-axis depicts the 4 experimental groups. The y-axis depicts the 33 DETE shared across groups.

Senescence induction co activates genes within the viral recognition pathway

Previous studies demonstrating derepression of TEs in the presence of HMAs, have shown induction of dsRNA with downstream activation of a suicidal interferon response [19,20,35]. Although we observed a similar TE dysregulation profile within each model of senescence induction, downstream activation of interferon stimulated genes (ISGs) involved in regulation of viral pathogenesis, demonstrated variability (Figure 7, Supplemental Data 5).

Figure 7. Viral recognition gene set analysis in senescence model. (a) Viral recognition gene pathway enrichment following 5-Azacitidine treatment. The x-axis depicts, in order, genes involved in viral recognition such as sensors, adaptors, signal transducers, transcription factors, interferon genes, and downstream interferon stimulated genes. For each gene in this pathway, there exists a 95% confidence interval also depicted. The global pathway activation is depicted as a dashed blue line, and the corresponding confidence interval is a gray solid bar. The y-axis depicts the log-fold change enrichment score. (b) Viral recognition gene pathway following replicative senescence. The x-axis and y-axis are identical to (a). (c) Viral recognition gene pathway following adriamycin chemotherapy induced senescence. The x-axis and y-axis are identical to (a). Viral recognition gene pathway following oxidative stress. The x-axis and y-axis are identical to (a).

Upstream of ISG activation, the general innate immune response to viral infection includes replicating viral genome pattern recognition sensors (RIG-1/DDX58, MDA5, IFIH1, TLR3, IFI16), adaptor protein initiators of cascade signaling (MAVS, MYD88), and signal transducers involved with the activation of transcription factors (TANK, TBK1). The activated transcription factors such as IRF3, IRF7 lead to the activation of IFNB1, where IFNB1 drives the production of Type I/III IFNs and subsequent downstream ISGs [33]. The pattern recognition pathway was up-regulated significantly in Azacitidine (log-fold change 2.2159, adjusted p-value: 2.29e-04), replicative senescence (log-fold change: 1.6080, adjusted p-value: 0.0211), and Adriamycin chemotherapy treatment (log-fold change: 1.4116, adjusted p-value: 0.0038). However oxidative stress did not show any significantly dysregulated the pattern recognition response despite the significant up-regulation of TE (log-fold change: 1.8198, adjusted p-value: 0.5623) (Figure 7, Supplemental Data 5).

Similarly, pattern recognition sensors (RIG-1/DDX58, MDA5/IFIH1) and IFNB1 were significantly up-regulated in Azacitidine, replicative and Adriamycin mediated senescence. Consistent with a prior report [27], we too observed significant up-regulation of IRF7 with Azacitidine and replicative senescence. Next, further downstream, interferon stimulated genes such as ISG15, ISG20, IFITM1/2, RSAD2, BST2, OAS1, and IFIT3 were all significantly enriched in Azacitidine, replicative senescence and Adriamycin (p-value <0.05) (Figure 7). Thus, different modalities of senescence induction depict similarities in dysregulated TE profiles and downstream activation of the innate immune pathway in all modalities of senescence except for senescence induced by oxidative stress.

Since oncogenes are known to induce senescence in cancer, we used RNA sequencing data from BRAF600VE oncogene induced senescence (OIS) melanocyte cells from Pawlikowski et. al [28]. and observed a similar upregulation of TE along with genes in the viral recognition pathway (Supplemental Figure 2(a) and 2(b)). Gene expression array showed a statistically significant upregulation of the viral recognition pathway upon OIS in melanocytes compared with uninfected melanocytes (Supplemental Figure 2(c)).

Discussion

Prior studies have shown increased expression of repetitive elements and inflammatory response genes with aging[21,36–38]. This hitherto, is the first comprehensive analysis looking at the landscape of TE expression across various senescence models, demonstrating significant overlap in their TE expression profiles.

Elimination of senescent cells have previously shown increase in longevity[24–26]- partly mediated by immune cells [8–11]. In this regard, evidence suggests that TEs mediated activation of innate immune response facilitate clearance of cancer cells [19,20]. More directly, we recently showed that leukemic stem cells have repressed expression of TEs along with immune pathways [39], indicative of a survival mechanism to evade an innate immune response. Identification and targeting of similar oncogenic-addictions could prove essential therapeutic strategies, as recently shown [40].

Cytosolic DNA sensor has been shown to play a major in immune mediated clearance of senescent cells [41] and is being extensively investigated for cancer therapy. In contrast, TE-activated dsRNA recognition pathway, is known to mediate tumor clearance [14,17,19–21], although its role in senescence clearance, remains unexplored. This highlights the importance of intracellular DNA and RNA sensing pathways in benefitting organismal health by promoting clearance of cancerous and senescent cells. Our study adds to this knowledge by highlighting the role of TE expression in senescent cells, as a cellular adaptation to mediate clearance. Though MAVS is well recognized in the dsRNA recognition pathway, it has also been shown to be activated upon dsDNA recognition [39], hence the contribution of dsDNA in senescence clearance cannot be ruled out and warrants further investigation. A limitation to our analysis was the preparation of cDNA library using poly (A) selection, thus failing to capture most of SINEs and other classes of repeat elements.

Although, TE expression has previously been thought of as a deleterious pro oncogenic event, caused by re-integration and genomic instability [42], recent evidence suggests beneficial effects of TE expression. DNA methyltransferase inhibitors- activate TE and viral recognition pathway, thereby promoting cancer immunogenicity and tumor clearance [12,19,20,32]. Interestingly, expression of TE was recently shown to predict responses to immune checkpoint inhibitor in various cancer [43]. Here we present evidence for evolutionarily conserved non-coding elements, which are active in senescent cells, including in OIS, triggering an interferon-apoptotic cascade. Understanding these pathways and nodes, which regulate the expression of TE in senescence would lead to therapeutic targets to potentiate tissue rejuvenation and arrest oncogenic progression. The orderly execution of these cellular programs in healthy individuals suggests that viral elements in the genome have been co-opted by the mammalian immune system to balance regenerative capacity against tumorigenic risk.

Author contributions

AC performed the bioinformatics analysis of the RNA sequencing and wrote the manuscript. HKE analyzed the data and wrote the manuscript. GR designed and supervised the study, analyzed the data and wrote the manuscript.

Consent for publication

All authors have read the manuscript and consented for publication. Availability of data and material: The data was generated using RNA sequencing from publicly available data.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This project was funded by grants from Leukemia Lymphoma Society-Quest for Cures #0863-15, STOP Cancer, Tower Cancer Research Foundation, and the Keck School of Medicine of the University of Southern California, Jane Anne Nohl Division of Hematology and Center for the Study of Blood Diseases University of Southern California.