A Computational Framework to Characterize the Cancer Drug Induced Effect on Aging Using Transcriptomic Data

Cancer treatments such as chemotherapies may change or accelerate aging trajectories in cancer patients. Emerging evidence has shown that “omics” data can be used to study molecular changes of the aging process. Here, we integrated the drug-induced and normal aging transcriptomic data to computationally characterize the potential cancer drug-induced aging process in patients. Our analyses demonstrated that the aging-associated gene expression in the GTEx dataset can recapitulate the well-established aging hallmarks. We next characterized the drug-induced transcriptomic changes of 28 FDA approved cancer drugs in brain, kidney, muscle, and adipose tissues. Further drug-aging interaction analysis identified 34 potential drug regulated aging events. Those events include aging accelerating effects of vandetanib (Caprelsa®) and dasatinib (Sprycel®) in brain and muscle, respectively. Our result also demonstrated aging protective effect of vorinostat (Zolinza®), everolimus (Afinitor®), and bosutinib (Bosulif®) in brain.


INTRODUCTION
Cancer survival has been significantly improved because of better screen/diagnosis strategies and more effective treatments. The number of cancer survivors is projected to increase to 22.2 million by 2030 in the United States. However, the cancer survivors are at increased risk for accelerated aging and related health conditions (Guida et al., 2019). The normal aging process is described as a gradual accumulation of molecular and cellular damage, which eventually leads to systematic dysregulation. Cancer treatments such as chemotherapies can also cause genotoxic and cytotoxic damage to normal tissues in the cancer patients. This may change or accelerate the aging process in multiple organs of patients in long term. Emerging clinical studies have shown that cancer treatments can lead to a broad spectrum of aging-related health conditions at a younger age for the cancer survivor (Guida et al., 2021).
Clinical follow-up studies indicate that cancer treatment contributes to the early onset of agingrelated symptoms in the young cancer survivors, such as frailty, incident comorbidities, functional loss, and cognitive decline (Guida et al., 2019;Guida et al., 2021). Furthermore, observational studies have shown that survivors of adult-onset cancers have a higher burden of mobility limitations (Keating et al., 2005), comorbid conditions (Alfano et al., 2017), pain (Alfano et al., 2017), and a greater risk of functional and cognitive impairments compared with healthy, age-matched controls (Hewitt et al., 2003). These preclinical and clinical findings suggest that cancer treatments may lead to a change or acceleration of the aging trajectory. Several cancer drugs have been shown to cause cell damage through mechanisms similar to those mediating the aging process (Guida et al., 2019). However, for most of the FDA approved cancer drugs, the molecular and cellular changes underlying the interaction of cancer treatment and altered aging trajectory are unknown. This limitation has constrained the efforts to identify, predict, and mitigate the agingrelated consequences for cancer survivors (Guida et al., 2021).
Emerging evidence has shown that transcriptome and other types of "omics" data can be used to study molecular changes and trajectory of the aging process (Edwards et al., 2007;Peters et al., 2015). For example, by studying the epigenome and transcriptome landscapes of mice in different age groups, Benayoun and collogues have revealed widespread induction of inflammatory responses during the aging process (Benayoun et al., 2019). Omics data analyses have shown that under-expression of mitochondrial gene in tissues from aged donors (Yang et al., 2016). Study in transcriptomes across multiple species with varied lifespans demonstrated feasibility of using gene expression analysis to characterize the molecular signatures of longevity (Ma et al., 2018). DNA-methylation aging markers have been identified using the epigenetic profile as a "clock" of aging (Horvath and Raj, 2018). Despite recent development in transcriptomic and epigenetic research in the normal aging process, limited work has been done to characterize the cancer drug-induced aging process in cancer survivors.
This gap of knowledge is in part because of the paucity of samples and data that can be obtained after cancer treatment. The LINCS L1000 (L1000) transcriptome database is a comprehensive gene expression knowledgebase of pre-and post-treatment cell lines (Subramanian et al., 2017). This collection of post-treated expression profiles is an important resource for finding the connections between drugs, therapeutic effects, and disease states. In addition to the L1000 data, the availability of the multiple normal tissue transcriptomic data from the Genotype-Tissue Expression (GTEx) database (Consortium et al., 2017) allows us to robustly characterize an aging-associated signature in each tissue. In this regard, we hypothesize that we can integrate the L1000 and GTEx databases to investigate the scope and molecular changes of cancer drug-induced aging processes by comparing the transcriptional profiles of normal aging tissues and post-treatment normal cell lines. To test this hypothesis, we have developed a computational framework to identify aging-associated signatures from GTEx, the druginduced expression signatures from L1000 database, and further compared these two signatures to study the drugaging interaction and underlying molecular changes during this process. Among the significant drug-aging interactions are aging-accelerating effect of vandetanib in brain and dasatinib in muscle. Meanwhile, the protective effect of vorinostat and bosutinib on brain aging, and everolimus' different roles in muscle and brain aging processes are characterized.
Our integrative study successfully characterized the drug-aging interaction maps of 28 FDA approved cancer drugs in four normal tissues.

Expression Data
To date, the L1000 database (Subramanian et al., 2017) has included more than 1.3 million post-treated expression profiles from both normal and cancer cell lines. The majority of these gene expression profiles comprise transcriptional responses of human cancer cells to chemical and genetic perturbations. Phase II L1000 data of the best inferred gene space were obtained from https://lincsproject. org/LINCS/.
The Genotype-Tissue Expression (GTEx) database (Consortium et al., 2017) contains gene expression profiles from multiple tissue sites across nearly 1,000 individuals. Version 8 (V8) data for 46 tissue types with at least 84 samples for following analysis were downloaded from GTEx Portal (https://www.gtexportal.org/).

Statistical Analysis
Spearman correlation coefficient (Schober et al., 2018) between each gene's expression and the age of the corresponding donors in each tissue type were calculated. The p-value <0.05 was used to select significant correlations. Gene Set Enrichment Analysis (GSEA) (Subramanian et al., 2005) were performed on the genes ranked based on their correlation with donor ages. Gene Ontology (Ashburner et al., 2000;Wu et al., 2020), KEGG (Kanehisa and Goto, 2000), and curated gene sets from the Molecular Signatures Database (MSigDB) (Liberzon et al., 2015) were used to identify the pathways that enriched with aging/treatment-associated genes. All the computational and statistical analyses were performed using R (version 3.6.2) and Python (version 3.8.0).

Identification of Aging-Associated and Drug-Induced Signatures
In the GTEx dataset, Spearman correlation analysis was performed to identify the gene expression that was significantly correlated with the ages of the donors in different tissue types respectively. In each tissue, genes with p-values <0.05 and Spearman correlation coefficient rho> 0.2/ < −0.2 were selected as aging-associated signature genes. In some tissue types, we found very few aging-associated genes using the above criterion which might be caused by the small sample size. Then we chose 500 genes with top positive correlation and 500 genes with top negative correlation with donor ages as the signature genes for the following study. The signature r contains signature genes with positive correlation assigned with "+1," signature genes with negative correlation assigned with "−1," and the rest assigned with "0".
In the L1000 dataset, drug-induced signatures were built on genes that significantly changed after treatment using experiments with high reproducibility among replicates. For the transcriptomes treated with the same drug, genes that have significant changes with z score >1 or < −1 in 20% or more experiments were selected as signature genes similar to previously described in (Stathias et al., 2018). After drug treatment in the normal cells, each drug-induced signature z includes significantly upregulated genes assigned with "+1," significantly downregulated genes assigned with "−1," and the rest of the genes marked as "0".

Integrating the Drug-Induced Signatures and Aging-Associated Signatures to Identify
Cancer Drug-Aging Interactions 9,792 overlapped genes between the LINCS L1000 database and GTEx data were used to investigate the interaction between drug-induced signatures and aging-associated signatures as shown in Figure 1. The specificity (S) score and concordance ratio (CR) similar to (Stathias et al., 2018) were used to quantify the signature interaction and evaluate the druginduced aging effects. The specificity score was calculated as the percentage of shared genes between the drug-induced signature in a specific cell line and the aging-associated signature in a matched tissue type comparing to the total number of genes in the aging-associated signature in the matched tissue type (Eq. 1). The range of the specificity score is between 0 and 1. Higher specificity score indicates higher overlap between two signatures.
The concordance ratio was defined as the ratio of the number of shared genes in both signatures that are upregulated/downregulated after drug treatment and positively/negatively correlated with increased donor ages versus the number of shared genes in both signatures that are upregulated/downregulated after drug treatment but negatively/positively correlated with increased donor ages (Eq. 2). Concordance ratio where z and r are the vectors of drug-induced signature and aging-associated signature respectively. An adjusted concordance ratio is further obtained by adding a pseudocount of 1 to both of numerator and denominator of the original ratio, which can be used to characterize the cases with denominator at zero when calculating the concordance ratio. For example, when no overlapped genes or no genes with opposite alteration direction were found between aging-associated signature and drug-induced signature, adjusted concordance ratio will provide a proper calculation. Moreover, a permutation test will be performed for each adjusted concordance ratio to show its significance. In each permutation, an adjusted concordance ratio will be calculated for the randomly selected aging-associated signature and drug-induced signature keeping the number of positive and negative genes remains the same with the original signatures. This process will be repeated for 10,000 permutations and a normal distribution is approximated for permutated adjusted concordance ratios, which will be used to calculate the p-value for the observed adjusted concordance ratio. A significant adjusted concordance ratio <1 suggests a potential protective effect from aging, and a significant concordance ratio >1 indicates more likely the drug can accelerate aging-related process in this specific tissue. Further enrichment analysis was performed on genes that were shared between aging-associated signature and drug-induced signature, using pathway annotations from Gene Ontology (GO) (Gene Ontology, 2004) by Enrichr (Chen et al., 2013) with Fisher's exact test.

Aging-Associated Gene Expression in the GTEx can Recapitulate the Well-Established Aging Hallmarks
The GTEx database contains gene expression profiles from multiple tissue sites across nearly 1,000 individuals aged from 20 to 71. We first sought to determine if gene expression data can be used to characterize the molecular changes during the normal aging trajectory. For each of 46 tissue types that has at least 84 samples in the GTEx database (V8), we performed GSEA analyses on the genes ranked based on their expression correlation with donor ages. This analysis revealed that the aging-associated transcriptomic changes are highly consistent with established aging hallmarks (Lopez-Otin et al., 2013). As shown in Figure  tissue types) ( Figure 2). These results suggest the aging process in multiple tissue types can be characterized by the gene expression profiles.
In the amygdala and hippocampus regions of the brain, the aging-associated genes are highly overlapped with characteristics associated with neurodegeneration phenotypes such as "IMPAIRED SOCIAL INTERACTIONS" ( Figure 3A). Specifically, PTEN induced kinase 1 (PINK1) is one of the core genes that are significantly decreased in aged donors in the IMPAIRED SOCIAL INTERACTIONS function ( Figure 3B). Previous studies FIGURE 3 | The aging-associated genes in normal brain and muscle tissues are highly overlapped with brain and muscle aging phenotypes. (A) GSEA analysis of aging-associated genes using gene sets derived from the Human Phenotype Ontology. (B) Decreased gene expression in normal brain and muscle tissues of aged donors.
Frontiers in Pharmacology | www.frontiersin.org June 2022 | Volume 13 | Article 906429 5 have shown that PINK1 plays a vital role in mitochondrial maintenance (Wilhelmus et al., 2011), and PINK1 dysregulations may contribute to neurodegeneration and aging (Kitagishi et al., 2017), which is consistent with our results here. Dysfunction of PINK1 related pathway has been reported in patients with Parkinson's and Alzheimer's diseases Martin-Maestro et al., 2016), which makes it a candidate therapeutic target for those diseases (Du et al., 2017).
Meanwhile, neurodegeneration associated function "CEREBRAL WHITE MATTER ATROPHY" (Patel and Barkovich, 2002;Poretti and Boltshauser, 2015) demonstrated enrichment with genes that significantly decreased in aged donors' cerebellum and cerebellar hemisphere of brain (FDR <0.05, Figure 3A). Among 29 core genes in this functional module, Hydroxysteroid 17-Beta Dehydrogenase 4 (HSD17B4) showed significantly low expression as the donor's age increases (Spearman p < 10 −4 , Figure 3B). HSD17B4 has been reported to play an important role in the catalyzing β oxidation of very long chain fatty acids (VLCFA) (Violante et al., 2019), and mutations in this gene can cause progressive cerebral atrophy in clinical cases (Amor et al., 2016;Landau et al., 2020).
In the muscle skeleton samples, the significantly downregulated genes in aged donors were enriched in functions of "AXIAL MUSCLE WEAKNESS" and "DISTAL LOWER LIMB MUSCLE WEAKNESS" (FDR <0.15, Figure 3A). For example, striated preferentially expressed gene (SPEG) and Filamin C (FLNC), involved in the two pathways respectively, showed significantly decreased expression in muscle samples of aged donors (Spearman p < 10 −7 and p < 10 −6 , Figure 3B). Consistent with our results, SPEG deficiency in skeletal muscle was found to cause defective calcium handling and excitation-contraction coupling, further lead to congenital myopathies (Huntoon et al., 2018). Meanwhile, previously research revealed a clinical case with FLNC mutations experienced weakness in the lower limbs and proximal muscles (Chen et al., 2019). Moreover, previous study has shown that genetic variants in FLNC is one of the prevalent causes of myopathies and cardiomyopathies (Verdonschot et al., 2020). Our results suggest muscle weakness, the main phenotype of aging process in skeletal muscle, can be characterized by the transcriptome profiles from donors of different ages.

Cancer Drug Induced Transcriptomic Changes Are Highly Resembled to Normal Aging Processes in Multiple Tissue Types
In addition to the molecular and cellular processes, we also investigated if aging-associated genes in multiple tissue types can enrich previously identified cancer drug-induced pathways. The GSEA analysis revealed that cisplatin induced genes are highly enriched in the overexpressed genes in 27 out of 46 tissue types of the aged donors, including multiple regions of brain, adipose, muscle and kidney (FDR <0.1, Figure 4A). For example, genes with positive expression correlation with increased donor ages in the brain cortex region are highly overlapped with the core genes that are induced after cisplatin treatment in cancer cells (KERLEY RESPONSE TO CISPLATIN UP (Kerley-Hamilton et al., 2005), FDR <0.1, Figure 4B). In particular, cyclin dependent kinase inhibitor 1A (CDKN1A) is one of shared genes exhibited significant over-expression as the donor's age increase in kidney and brain cortex region ( Figure 4C). This observation is consistent with previous report of CDKN1A's important role in the kidney injury (Zaidan et al., 2020), brain aging (Belsky et al., 2015) and cisplatin cell killing effect (Zamagni et al., 2020). Clinically, it has been shown that cisplatin treatment can lead to kidney injury (Ozkok and Edelstein, 2014) and cognitive impairment (John et al., 2017) in cancer patients.
In addition to cisplatin, we have also found the association between irradiation therapy and aging in 19 out of 46 sample types (GHANDHI BYSTANDER IRRADIATION UP and WARTERS IR RESPONSE 5GY, FDR <0.1, Figure 4A) (Kim et al., 2014). These observations suggested the impact of cisplatin and radiation treatments on aging process in multiple normal tissues.

Drug-Induced Transcriptomic Alterations Recapitulate the Mechanism of Action of the Treatment in L1000 Data
According to the record of National Cancer Institute, there are 28 FDA approved cancer drugs having enough treatment data points (e.g., dosage and duration) in at least one of four normal/ immortalized cell lines in the L1000 data ( Figure 5A). Those four cell lines are originated from brain (neural progenitor cells We first sought to determine if drug-induced transcriptomic alterations can recapitulate the cancer drugs' known effect on tissue. In this regard, we performed GSEA analyses on druginduced transcriptomic data of the 28 FDA approved drugs (Liberzon et al., 2015). These analyses revealed that the druginduced transcriptomic alterations are consistent with the known effect in 22 out 28 FDA drugs (Supplementary Table S1). For example, steroid biosynthesis pathway (METABOLISM OF STEROIDS) was enriched with upregulated genes after treated with selective estrogen receptor modulators (SERMs) such as tamoxifen and raloxifene in HA1E cells (FDR< 0.0001, Figures  5B,C). SERMs are a group of non-steroidal drugs with the ability to bind to estrogen receptors and can upregulate steroid metabolismrelated genes by interacting with sterol regulatory element-binding protein (SREBP-2) (Fernández-Suárez et al., 2021). Moreover, important enzymes involved in the "METABOLISM OF STEROIDS" pathway such as farnesyl diphosphate synthetase (FDPS), 7-dehydrocholesterol reductase (DHCR7), methylsterol monooxygenase1 (MSMO1) were upregulated by tamoxifen treatment in a dose-dependent manner ( Figure 5F).
For chemotherapy drugs such as mitoxantrone and doxorubicin, cell cycle related pathway "CELL CYCLE CHECKPOINTS" was significantly downregulated after these two treatments (FDR<0.0001, Figures 5D,E). Mitoxantrone are known to work against cancer by killing fast-growing cells. Consistent with its mechanism of, cell cycle regulatory genes such as neuroepithelial cell transforming 1 (NET1), TGFB induced factor homeobox 1 (TGIF1), and checkpoint kinase 2 (CHEK2) were downregulated by mitoxantrone treatment in a dose-dependent manner ( Figure 5G). These results suggested the drug-induced transcriptomic alterations and functions are consistent with the mechanisms of action for the treatment.

Interaction Between Cancer Drug-Induced Signatures and Aging-Associated Signatures in Multiple Normal Tissue Types
To determine the drug-aging interaction, we first identified agingassociated signatures in each tissue type in the GTEx dataset and drug-induced signatures in each normal cell line in the L1000 dataset as described in the Section 2. In GTEx database, there are 17 tissue types/regions with sufficient samples (ranging from 85 to 803) that match the 4 normal cell lines tissue origins in the L1000 data (Supplementary Figure S1A). The distribution of the correlation between gene expression and the donor's age was demonstrated in each tissue type ( Figure 6A). Then the number of selected genes varies at different rho thresholds (Supplementary Figure S1B). A reasonable number of genes were selected at rho>0.2 in the aging-associated signatures in different tissues. Interestingly, the total number of signature genes and the percentage of positively/negatively correlated genes are different in each tissue type (Supplementary Figure S1C). For example, a large proportion of signature genes are downregulated in some of the brain tissues of the aged donors while more genes are upregulated in adipose and muscle tissues of the aged donors. This difference might be affected by different sample sizes, but also indicate tissue-specificity of aging-associated signature. Meanwhile, we have observed that the different number of genes whose expression are altered by different drugs, which may be caused by the distinct mechanism of actions (MOA) of the drugs (Supplementary Table S3). Chemotherapy drugs such as doxorubicin and mitoxantrone tend to induce an extensive transcriptomic change in the cell, while targeted therapy drugs may induce less changes in the gene expression.
For each drug-induced signature and aging-associated signature, we calculated a specificity (S) score and an adjusted concordance ratio (CR) to quantify the signature interaction and evaluate the drug-induced aging effects (Figure 1 and see details in Section 2). The drug-aging interaction with a higher specificity score indicates that a drug-induced signature has more overlapped genes with an aging-associated signatures in matched tissue type. On the other hand, the significant adjusted concordance ratio indicates the drug's accelerating (i.e., CR > 1) or protective (i.e., CR < 1) effect on the aging trajectory of corresponding tissue. There are 326 pairs of drug-induced signature and aginginduced signature with matched cell tissue origins between GTEx and L1000 datasets. 34% of the drug-aging interactions have a specificity score larger than 0.45 ( Figure 6B). For example, chemotherapy drug mitoxantrone and doxorubicin showed high aging specificity scores at 0.749 and 0.756 in kidney. Meanwhile, 36% of the drug-aging interactions have a significant adjusted concordance ratio with p < 0.01, including highly significant concordance ratios between vandetanib and aging of multiple brain regions ( Figure 6B).
Using a threshold of specificity>0.45, CR < 0.9 or >1.1 and concordance p < 0.01, we characterized 34 potential drugregulated aging events in 15 tissue types/regions ( Figure 6C; Supplementary Table S2). Among them, ten drugs were predicted to have 22 events of aging-accelerating effects in 14 tissue types/regions. In particular, doxorubicin was predicted to have aging-accelerating effect in kidney aging with specificity at 0.756 and significantly high adjusted concordance ratio at 1.4 (p < 10 −4 ), which is consistent with doxorubicin's established side effect in kidney damage and aging (Su et al., 2015;Yemm et al., 2018). In addition to the aging accelerating effect, our analyses also identified five drugs showing 12 aging protective events in 9 tissue types/regions.

Drug-Aging Interaction Analysis Revealed That Vandetanib (Caprelsa ® ) Treatment May Accelerate Brain Aging
Vandetanib (Caprelsa ® ) is a tyrosine-kinase inhibitor (TKI) of vascular endothelial growth factor receptor (VEGFR), epidermal growth factor receptor (EGFR), and RET tyrosine kinase (LoRusso and Eder, 2008). It has been approved by FDA to treat medullary thyroid cancer. In our study, vandetanib Frontiers in Pharmacology | www.frontiersin.org June 2022 | Volume 13 | Article 906429 9 treatment in NPC cells induced gene expression changes demonstrating specificity scores higher than 0.45 and significantly high adjusted concordance ratios with agingassociated signatures in multiple brain regions including amygdala, anterior cingulate cortex, cerebellar hemisphere, cerebellum, cortex, frontal cortex, hypothalamus, nucleus accumbens, and substantia nigra (CR > 1.49, p < 0.01, Figure 6C; Supplementary Table S2). Pathways such as "Mitotic cell cycle phase transition" (Enrichr p-value <0.05), "Negative regulation of cellular senescence" (Enrichr p-value <0.05) and "Regulation of cellular amide metabolic process" (Enrichr p-value <0.001) were significantly enriched with the genes that are downregulated in both vandetanib-treated NPC cells and aged substantia nigra tissues ( Figure 6D). As a drug that can penetrate the blood-brain barrier (Subbiah et al., 2015), vandetanib treatment may cause dysregulation of cell senescence, mitotic cell cycle, amide metabolism in brain tissue and accelerate aging process of brain. Consistent with our discovery, a recent study has demonstrated that vandetanib exert a deleterious effect on the dopaminergic system in a Parkinson's disease model (Requejo et al., 2018).

Drug-Aging Interaction Analysis Revealed That Dasatinib (Sprycel ® ) Treatment May Accelerate Muscle Aging
Our analysis has identified another tyrosine-kinase inhibitor, dasatinib, tend to trigger aging process in muscle tissue.
Dasatinib (Sprycel ® ) was approved to treat Philadelphia chromosome-positive (Ph+) chronic myeloid leukemia and acute lymphoblastic leukemia for both children and adult patients (Keam, 2008). The dasatinib-induced signature in SKL cells has a high specificity score (S = 0.596) and significantly high adjusted concordance ratio (CR = 2.21, p < 10 −4 ) with agingassociated signature in muscle tissues. Moreover, the overlapped genes between two signatures were significantly enriched in the interferon gamma mediated signaling pathway (Enrichr p < 0.001) and alpha/beta T cell activation (Enrichr p < 0.05). These results suggest dasatinib may trigger inflammatory response which is observed over-expressed in muscle tissue of aged donors ( Figure 6E). Indeed, Naif I AlJohani and co-authors (AlJohani et al., 2015) reported a case of inflammatory myopathy in a patient with chronic myeloid leukemia after treated with dasatinib. There was also a clinical study demonstrating that chronic myeloid leukemia patients on TKI therapy showed significantly more muscle fatigue than control groups (Janssen et al., 2019).

Drug-Aging Interaction Analysis Revealed That Vorinostat (Zolinza ® ) and Bosutinib (Bosulif ® ) Treatment May Exert Protective Effect of Brain Aging
The histone deacetylase (HDAC) inhibitor, vorinostat (Zolinza ® ), was identified to provide protective effect of brain aging. Vorinostat is a FDA approved drug for cutaneous T cell lymphoma treatment (Iwamoto et al., 2013). Vorinostat-induced signature in NPC cells had overall high specificity scores, but significantly low adjusted concordance ratios with agingassociated signatures in multiple brain tissues (hippocampus: S = 0.629, CR = 0.744, p < 10 −5 ; anterior cingulate cortex: S = 0.639, CR = 0.768, p < 10 −5 ; amygdala: S = 0.654, CR = 0.83, p < 10 −9 ). Overlapped genes between aging-associated signature and vorinostat-induced signature were significantly enriched in mitochondrial metabolism pathways such as aerobic respiration (Enrichr p < 0.05) ( Figure 6F). Our results suggest that vorinostat treatment may reverse mitochondrial dysfunction in brain tissue, which is among the main factors involved in neurodegeneration. Interestingly, in vitro study by Surabhi Shukla et al. (2016) showed that vorinostat can independently induce neuritogenesis, and vorinostat treatment can confer memory reinstatement in several cognitive decline mouse models (Kilgore et al., 2010;Benito et al., 2015). Moreover, a phase I clinical trial (NCT03056495) is underway to evaluate the neuroprotective effect of vorinostat in patients with mild Alzheimer disease (VostatAD01, 2017), which supports our discovery of the interactions between vorinostat and brain aging.
Another drug showed protective effect of brain aging is bosutinib (Bosulif ® ). It is a second generation tyrosine-kinase inhibitor for Bcr-Abl and Src family kinases, which is used for the treatment of chronic myeloid leukemia (Isfort and Brümmendorf, 2018). Bosutinib-induced signature in NPC cells had an overall high specificity, but significantly low concordance ratio with aging-associated signatures in multiple brain regions (hippocampus: S = 0.607, CR = 0.647, p < 10 −3 ; cortex: S = 0.583, CR = 0.589, p < 10 −7 ; substantia nigra: S = 0.579, CR = 0.632, p < 10 −2 ). Bosutinib's predicted aging protective effect of brain is supported by previous and ongoing studies. In preclinical studies, bosutinib was reported to facilitate the clearance of α-synuclein in Parkinson's disease (Hebron et al., 2014). Bosutinib treatment also enhanced the clearance of neurotoxic proteins α-amyloid and tau, leading to cognitive improvement in Alzheimer's disease mouse models (Lonskaya et al., 2013a;Lonskaya et al., 2013b;Hebron et al., 2018). In 2016, a phase I clinical trial (NCT02921477) began to evaluate bosutinib's effect on patients with mild cognitive impairment and the results support an overall positive outcome after 1 year treatment of bosutinib (National Library of Medicine, 2016; Mahdavi et al., 2021). In 2019, a phase II trial (NCT03888222) was started to test bosutinib as a possible treatment in dementia with Lewy bodies (National Library of Medicine, 2019).

Drug-Aging Interaction Analysis Revealed That Everolimus (Afinitor ® ) Treatment Demonstrated Aging Protective Effect in Brain but Aging Accelerating Effect in Muscle
The mTOR inhibitor everolimus (Afinitor ® ) was used to treat multiple cancers (Motzer et al., 2010;Baselga et al., 2011;Yao et al., 2011) and to suppress immunity to prevent rejection in patients having organ transplantation (Tedesco-Silva et al., 2022). Our interaction results showed different effects of the same drug on aging of brain and muscle tissues. Everolimusinduced signature in SKL cells showed aging-accelerating effect in muscle tissue (S = 0.610, CR = 2.765, p < 10 −8 ).
The overlapped positive genes between the two signatures were enriched in functions "DNA damage response by p53," "interferon-gamma mediated signaling pathway" and Frontiers in Pharmacology | www.frontiersin.org June 2022 | Volume 13 | Article 906429 "positive regulation of chemotaxis" (Enrich p < 0.05, Figure 6G). Consistent with our prediction, previous studies reported that muscle wasting can be induced by everolimus treatment in cancer patients (Albiges et al., 2011;Gyawali et al., 2016). Interestingly, the interaction between everolimus-induced signature and aging-associated signature in brain hippocampus demonstrated agingprotective effect (S = 0.650, CR = 0.424, p < 10 −5 ). In vivo studies from multiple research groups showed that infusion of everolimus restored cognitive function in Alzheimer's disease models (Fanoudi et al., 2018;Cassano et al., 2019). As a cancer drug that is able to penetrate the blood-brain barrier (Subbiah et al., 2015), everolimus's potential mechanism in cognitive protection is the protection of intact blood-brain barrier and limited entry of proinflammatory and neurotoxic factors into the brain tissue Van Skike and Galvan, 2018).

DISCUSSION
Although laboratory and clinical evidence has shown that several cancer drugs can change the cancer patient aging trajectory, for most of the FDA approved cancer drugs, whether and how they influence the aging process of cancer patient is unknown. More importantly, the molecular and cellular changes underlying the interaction of cancer treatment and altered aging trajectory remain elusive. In this study, we have developed a computational framework to integrate the drug post-treatment transcriptomic data and human tissue transcriptomics data to investigate the druginduced aging processes in multiple human organs. This strategy allows us to first identify aging-associated signatures in normal aging tissues. We then link the drug treatment effect with aging by comparing the aging-associated signatures to the drug-induced expression signatures using specificity scores and adjusted concordance ratios. By applying this computational framework to L1000 and GTEx databases, we have successfully identified experimentally validated druginduced effects on aging such as aging-accelerating effects of vandetanib and dasatinib on brain and muscle, respectively. Our result also identified aging-protective effect of vorinostat, everolimus, and bosutinib on brain. With the availability of new data from L1000, GTEx and other databases, our computational strategies can be extended to future analysis of more drugs and tissues. Our study is an important complement to the clinical and preclinical study. The results from our study can serve as an initial screen of cancer treatments that are likely to induce aging-related consequences in patients. This will provide potential candidate drugs for long term clinical and preclinical study. Moreover, our study can also identify drugs that provide protective effects to tissue aging, which could be considered as candidate drugs that could treat drug-induced side effect of aging or co-treatments to prevent possible aging process. In the future, we will be interested to further utilize this framework to characterize the drugs' synergistic effect on aging for combination therapy.
Another advantage of our study is that our analysis can help characterize the underlying molecular changes during the drug-aging interaction. For example, Requejo and collogues have used rat model to demonstrate that vandetanib treatment significantly increased behavioral impairment in 6-hydroxydopamine induced preclinical model of Parkinson's Disease (Requejo et al., 2018). They also observed morphological changes in brain after vandetanib treatment including the decrease of TH-immunopositive striatal volume and the decreased axodendritic network in the substantia nigra. However, the mechanism of how this tyrosine kinase inhibitor influences the dopaminergic system is unknown. In this study, our analyses not only recapitulate the deleterious effect of vandetanib on Parkinson's Disease, but also reveal that vandetanib treatment causes dysregulation of cell senescence, mitotic cell cycle, amide metabolism in NPC cells. These cellular processes were significantly downregulated in aged substantia nigra human tissues. These observations, combined with the published preclinical animal model result, suggest that vandetanib-mediated cell senescence and metabolism disruption in substantia nigra region may be the mechanism of vandetanib's deleterious effect on Parkinson's Disease. By providing the underlying molecular and cellular changes of the drug-induced aging process, our study can support future mechanistic studies and the development of therapeutic strategies to mitigate the drug-induced aging process.
Our study has limitations. One limitation of our study is that it is hard to evaluate the dosage effect of each cancer drug on the aging process in patients. Our drug-induced expression signatures are learned from in vitro cell line assays of the L1000 database. The dosages used to treat the normal cell lines in our analysis not necessarily represent the drug concentration in the patient's normal tissues (e.g., brain) during cancer treatment. In the future study, we will consider the pharmacokinetics and drug tissue distribution (e.g., blood-brain barrier) data to adjust for the practical drug concentration that a cancer patient's uptake in different tissues. Another limitation of our study is that we established the aging signature in different tissues assuming the tissues are from "healthy" people undergoing natural aging process. However, in reality, aging patients usually have comorbidities such as diabetes or dementia which may influence drugs' effect on aging. In the future study, we will construct aging signatures with different disease conditions if we can obtain the comorbidities information from the patients. We will build a multivariable regression model that includes the comorbidities as confounding factors to learn the interaction among drugs, aging, and comorbidities.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
MZ conceived the study. YZ, YW, and MZ contributed to the method development. YZ, YW, and MZ performed bioinformatics analysis. YZ, KS, DY, and MZ performed interpretation. YZ and MZ wrote the manuscript while all authors assisted in manuscript review, editing and preparation. All authors read and approved the final manuscript.

FUNDING
This study was supported by the Shear Family Foundation, the American Cancer Society Research Scholar Award (132632-RSG-18-179-01-RMC), and National Cancer Institute (1R01CA222274 and R01CA255196).