Transcriptome analysis of mouse aortae reveals multiple novel pathways regulated by aging

RNA-seq reveals high divergence of transcriptome profiles between young and old aortae

We performed RNA-seq analysis on 5 biological replicates of male aortae from each of the young (age 14.5 weeks) and old (age 92.4 weeks) groups. These old mice displayed typical aging phenotype evidenced by decreased Doppler aortic flow velocity (Supplementary Figure 1A) and left ventricular ejection fraction (EF) value (data not shown), evident fibrosis in carotid arteries (Supplementary Figure 1B), as well as induced expression of senescence marker proteins (p21 and p53) in aortae (Supplementary Figure 1C). We obtained approximately 25 million reads per sample (FastP). ~90% of all reads were mapped to the reference mouse genome (mg38+gencode M22 Annotation, STAR). Differences related to sequencing features were not found between two groups. The principle component analysis (PCA) showed that these samples were tightly clustered within each group and there was a clear transcriptome difference between groups (Figure 1A). Next, we conducted differential expression analysis (DEA), with an adjusted p-value threshold of 0.05 on each set of raw expression measures. Overall, there were 2,321 significantly differentially expressed transcripts, which represented 8.3% of the total number of detected transcripts. 1,084 (46.7%) transcripts were downregulated while 1,237 (53.3%) transcripts were upregulated in the aged aorta. More than 90% of differentially expressed genes were found to be protein coding genes, whereas differentially expressed non-coding genes, including long noncoding RNAs (lncRNAs), small ncRNA (microRNA, small nuclear RNA, and small nucleolar RNA), and pseudogenes only accounted for ~5% (Figure 1B). Hierarchical clustering analysis based on the regularized log transformation of the normalized count data showed that all samples clustered by age, indicating that aortic cells from young and old groups have distinct gene expression profiles (Figure 1C). In order to visualize the magnitude of the fold-changes revealed between the two groups, we generated a volcano plot. Genes with adjusted p<0.05 were colored either green or red according to the direction of the fold-change. A line was drawn at the unadjusted p-value of 0.05 for reference. We found that the most significantly downregulated genes in aged aortae included Edn1, a potent vasoconstrictor [22], Nfil3, a key regulator of circadian rhythm [23], and Atf3, a stress responsive transcription factor [24]. The most significantly upregulated genes in aged aortae included Tnfrsf11b, a regulator of osteoclast development and risk factor for CVDs [25, 26], and Ccl8, a monocyte chemoattractant protein (MCP) family member driving inflammation [27]. In addition, several tumor suppressor genes such as Prune2 [28] and Efemp1 [29] were also induced in the old aorta (Figure 1D). Taken together, these results illustrate that substantial numbers of genes are regulated in the aorta in response to aging.

Functional characterization of gene expression changes between young and old aortae

To identify molecular functions of the aortic transcriptome altered by age, we performed Gene Ontology (GO) enrichment analysis using biological process terms and KEGG pathway analysis. Genes with adjusted p<0.05 and abs (log2FoldChange) > 0 were submitted to Enrichr to identify significantly enriched pathways and transcription factors. Top 30 upregulated and downregulated GO categories/KEGG pathways were presented in Tables 1 and 2, respectively.

As shown in Figure 2A, the majority of upregulated pathways in aged aortae were associated with immune activity and inflammation, such as cytokine-mediated signaling pathway (GO:0019221), cellular response to cytokine stimulus (GO:0071345), and neutrophil activation involved in immune response (GO:0002283). The major immune response pathway, type I interferon (IFN) signaling pathway was also significantly induced. Examination of important inflammatory mediators (cytokines, chemokines, and signaling proteins) revealed significant increases in gene expression with aging (Figure 2C). In accord with GO enrichment analysis, tuberculosis, which is an infectious disease closely related to type I IFN [30], was among the top enriched upregulated KEGG pathways (Figure 2B). These results are consistent with previous findings that chronic low-grade inflammation occurs in aged arteries, which could trigger vascular cell senescence and contribute to age-related impairment of vascular function. As expected, Cdkn2a (p16^INK4a), one of the main cyclin-dependent kinase inhibitors (CDKIs) driving the cell cycle arrest in senescence and a widely accepted specific marker for senescence [31], displayed more than 7-fold increase in the aged aorta (Figure 2D). This data is in support of previous studies demonstrating that p16^INK4a transcriptional activation is a hallmark of senescence in vivo [32–34].

Additionally, both GO enrichment and KEGG pathway analysis uncovered a significant difference between young and old mouse aortae for genes participating in phagocytosis regulation (Figure 2A, 2B), suggesting an increased capacity to engulf cellular debris or apoptotic cells in the aged aorta. Consistently, C4b, C2, and C3, which are central parts of the complement system enhancing the ability of phagocytic cells to clear microbes and damaged cells [35], were higher in old aortae (Figure 2E). The upregulation of phagocytosis and apoptotic cell clearance may be a feedback response to cellular damage triggered by aging.

Among the top 10 downregulated pathways in the aged aorta, ECM organization (GO:0030198) was the most significantly altered (Figure 3A), wherein 41 out of 230 genes were downregulated. A similar result was obtained from KEGG pathway analysis, which identified 22 out of 83 downregulated genes in the ECM-receptor interaction pathway (Figure 3B). We further assessed the expression of representative genes involved in ECM and collagen fibril organization. Aging caused a reduction in transcriptional levels of multiple collagen genes, including Col1a1, Col1a2, Col3a1, Col4a1, and Col5a2, which are abundantly expressed in the vasculature (Figure 3C). Other ECM components or regulators, including Fbn1, Lox, Lamb1, and Serpinh1, also exhibited downregulation with aging (Figure 3C). The related GO biological process terms and KEGG pathways, as well as the involved differentially expressed genes, are presented in Table 3. This observation is consistent with a previous study in the aorta of aging rats [36] and the microarray analysis in the aorta of aging monkeys [37], suggesting that differences in the expression of key ECM genes can directly account for the impairment of collagen/elastin crosslinking, resulting in age-related large artery stiffening. Moreover, several well-known matrix metalloproteinase (MMP) genes, including Mmp3, Mmp9, and Mmp10, were elevated in aged aortae (Figure 3D). This increment of MMPs together with the decrement of ECM biosynthesis would compromise the structural integrity of the vasculature, thus contributing to adverse age-associated vascular remodeling.

Pathways related to protein folding control and stress response are predominantly down-regulated in old aortae

In addition to the cluster of genes involved in ECM organization, we observed that “response to unfolded protein (GO:0006986)” and “chaperone mediated protein folding requiring cofactor (GO:0051085)” were among the top 10 enriched downregulated GO terms in the aged aorta (Figure 3A). Several other pathways associated with protein folding process such as “chaperone-mediated protein complex assembly (GO:0051131)” and “'de novo' posttranslational protein folding (GO:0051084)” also displayed significant downregulation (Tables 2, 4). These pathways consisted of numerous differentially expressed genes controlling protein folding, such as heat shock protein family members (e.g. Hspa8, Hsp90aa1, Hspa1a), Dnaj heat shock protein family (Hsp40) members (e.g. Dnaja1, Dnajb1) and clusterin (Clu) (Figure 4). These findings indicate a functional deterioration in protein folding control in the old aorta, particularly chaperone activities, which has been reported as one of the age-related molecular alterations in the vasculature [1, 38, 39].

Reduced protein folding and endoplasmic reticulum (ER) network efficiency with aging results in the accumulation of unfolded proteins, which provokes ER stress and subsequently cellular ER stress response [40]. On a more general basis, aging vessels are confronted by various endogenous and exogenous stressors such as inflammation, oxidative stress, and DNA damage. Therefore, a notably impaired ability of stress response may be a major determinant for the development of age-related vascular pathologies [1]. Consistent with this notion, we found that the regulation of cellular response to heat (GO:1900034, adjusted p =1.91E-05) and stress (GO:0080135, adjusted p = 4.74E-05) were significantly downregulated in aged aortae (Table 4). Reduced expression of several representative genes implicated in cellular response to stress, including Camk2b, Bag2, Bag3, and Ier5, is shown in Figure 4.

Circadian rhythm genes are differentially expressed in young and old aortae

It has become increasingly appreciated that circadian rhythm disruption precipitates vascular dysfunction [41, 42]. For example, dysregulated expression of circadian clock components results in dysregulation of inflammatory process and ECM dynamics, which are key pathophysiological factors underlying CVDs [42]. To assess the impact of aging on aortic circadian rhythm, we examined the expression of multiple molecular clock genes (Figure 5). Arntl (brain-and-muscle ARNT-like protein, also called BMAL1) and Clock (circadian locomotor output cycle kaput), the two core transcription factors controlling the transcriptional-translational feedback loop (TTFL) of mammalian clock, exhibited a ~50% and ~30% decrease in old aortae, respectively. Other downregulated clock genes were Nfil3 (repressor of Per transcriptional activity) and Npas2 (neuronal PAS domain protein 2, Clock paralogue), concurring with an upregulation of negative feedback loop components, including Dbp (activator of Per1 transcription), Bhlhe41 (Basic Helix-Loop-Helix Family Member E41), and Per2 (repressor of BMAL1/CLOCK complex transcriptional activity) (Figure 5). Altogether, our analysis shows that the transcriptional level of circadian rhythm regulators differs between young and old aortae, supporting the concept that altered circadian rhythm couples with advancing age.

qRT-PCR validation of representative genes for major terms/pathways regulated by aging

To validate age-regulated pathways revealed by aforementioned RNA-seq analysis, we performed qRT-PCR to evaluate the expression of representative genes of several major pathways, including Col1a1, Col3a1, and Lox for ECM organization (Figure 6A), Hsp1a, Hsp1b, and Hsph1 for protein folding control (Figure 6B), as well as Nfil3, Bhlhe41, and Npas2 for circadian rhythm (Figure 6C). qRT-PCR showed that the expression of these representative genes exhibited similar alterations to those derived from RNA-seq analysis. These results confirmed the regulation of aging on ECM organization, protein folding control, and circadian rhythm pathways revealed by RNA-seq.

The differentiated VSMC phenotype is retained in aged aortae

VSMC phenotypic shifting from the quiescent differentiated phenotype to a synthetic dedifferentiated state underlies the pathogenesis of various vascular diseases, including restenosis, atherosclerosis, aneurysm formation, and AVF failure [17–19]. A hallmark of this phenotypic switch is the marked reduction of SMC marker genes. To ascertain whether similar phenotypic switching occurs in arterial aging, we first compared levels of transcripts of SMC markers in young versus old aortae. Bulk RNA-seq revealed a slight and insignificant decrease in SMC marker genes, including Myh11, Cnn1, Acta2, Lmod1, and the master regulator of SMC differentiation program, Myocd in old aortae (Figure 7A). Further, this subtle downregulation of VSMC marker expression was confirmed by western blot and immunostaining of both MYH11 and ACTA2 (Figure 7B–7C). To determine SMC lineage in old aortae, we sought to take advantage of an established SMC reporter mouse line, Myh11-Cre-ER^T2/mTmG system [19]. As expected, the medial layer of old aortae was predominately comprised of GFP labeled SMC-derived cells positive for MYH11 and ACTA2 (Figure 7D). These data suggest that though there is a subtle reduction in the expression of SMC markers, these cells maintain the differentiated SMC phenotype.

Animals

All experiments were performed in compliance with approved protocols of the Institutional Animal Care and Use Committee (IACUC) at Albany Medical College. Both young (2 months) and old (18 months) male C57BL/6 mice were purchased from Jackson laboratory. Both groups were fed with water and standard rodent chow diet for more than one month before sacrifice at the same time for experiments. Young (14.5 weeks) and old (92.4 weeks) mice were sacrificed at the same time followed by tissue isolation for RNA preparation. In order to track VSMCs lineage in old aortae, we bred Myh11-CreER^T2 (JAX, # 007576) to ROSA mT/mG mice (JAX, # 01979) to generate Myh11-Cre/ER^T2-mTmG mice. These lineage tracing mice were injected with Tamoxifen intraperitoneally at an early stage (14 weeks) to activate Cre recombinase in mature VSMCs as previously described to specifically label SMCs with membrane green fluorescence (mG) and other cell types with membrane red fluorescence (mT) [19]. These reporter mice were sacrificed at the age of 54 weeks for the indicated lineage tracing study.

RNA isolation and deep sequencing

Aortae were quickly excised from saline perfused mice and cleaned of periadventitial adipose tissue. Purified aortae (5 aortae/group) were homogenized by a Minilys homogenizer (Bertin Technologies, Rockville, MD, USA) before total RNA extraction using miRNeasy RNA extraction kit (Qiagen, Valenica, CA, USA). RNA integrity was examined with an Experion™ Automated Electrophoresis System by utilizing Experion RNA StdSens analysis kit (Bio-Rad, Hercules, CA, USA). Bulk RNA-seq was conducted by Genomics Research Center at the University of Rochester Medical Center. Detailed information regarding library construction and sequence depth was described previously [20]. Raw data derived from RNA-seq was filtered using FastP-0.20.0 (https://github.com/OpenGene-/fastp). To obtain the high-quality clean reads for subsequent analysis, reads below a length of 35 after read trimming by quality were removed from the sequencing data. Sequence reads were aligned to annotated transcripts on the UCSC Reference Genome (mg38 + gencode M22 Annotation) via STAR-2.7.0f. Subread-1.6.4 feature counts were used to count reads at the gene level. Differentially expressed gene analysis DESeq2-1.22.1 within R-3.5.1 was used to perform data normalization and differential expression analysis with an adjusted p-value threshold of 0.05 on each set of raw expression measures. The ‘lfcShrink’ method was applied to moderate log2 fold-changes for lowly expressed genes. The detailed RNA-seq information of this assay is available in GSE145972 deposited in NIH Gene Expression Omnibus (GEO) database.

Pathway enrichment analysis

Significant up-regulated and down-regulated genes based on p-adj <0.05 and abs (log2FoldChange) > 0 were submitted to Enrichr (https://amp.pharm.mssm.edu/Enrichr/) to identify significantly enriched pathways and transcription factors.

Quantitative reverse transcription-polymerase chain reaction (qRT-PCR)

After RNA samples were extracted from mouse aortic tissues, cDNA synthesis was conducted by using iScript cDNA kit (Bio-Rad). qRT-PCR was performed using Universal SYBR Green Supermix (Bio-Rad) and CFX386 Touch™ Real-Time PCR Detection System (Bio-Rad). mRNA levels were expressed relative to the loading control Gapdh. Technical duplicates of each sample were examined. PCR primers are included in Supplementary Table 1.

Immunofluorescence Staining

Freshly harvested mouse aortae were rinsed in cold 1XPBS and immediately fixed in 4% paraformaldehyde at 4 °C for 24 hours followed by OCT processing. 10 μm frozen sections were prepared for immunofluorescence (IF) staining as previously described [21]. For aortae isolated from mTmG reporter mice, we used trypsin-mediated antigen retrieval (REF TA-015–TR; Lab Vision Corporation, Fremont, CA) to retain the natural mTmG fluorescence prior to the immunofluorescence staining. Primary antibodies used were as follows: anti-MYH11 (1:500, Alfa Aesar, BT-562), anti-ACTA2 (1:2000, Sigma, A2547). Fluorescent signals were captured by a confocal microscope and processed by Photoshop (Adobe, San Jose, CA, USA). All images were captured and processed under equivalent conditions.

Statistical analysis

Statistical analysis was determined by an unpaired two-tailed Student’s t test using GraphPad Prism 6 to compare two groups. Data in graphs were presented as mean ± SD. p<0.05 was considered statistically significant. p<0.05 and p<0.01 were indicated by * and **, respectively, throughout all graphs.