Using Deep Learning to Detect Cell Aging from Protein Sequences

Senescence‐associated genes were first curated from the CellAge database. The database can be accessed at https://genomics.senescence.info/cells/index.html↗. To ensure functional specificity, only genes that were experimentally validated to induce cellular senescence in vitro (i.e., annotated as "Senescence Effect = induce") were retained, while genes labeled as unclear or with no reported senescence‐inducing effect were excluded. This process yielded 163 high‐confidence senescence‐inducing genes, as listed in Table S1.

To capture senescence states observed in physiological contexts, we further incorporated three in vivo senescence resources: SASP Atlas, SenMayo genes (Saul et al. 2022), and SenPy (Sanborn et al. 2025). Genes consistently identified in at least two of these datasets were retained, yielding 47 additional in vivo‐supported senescence signature genes. These genes represent conserved and tissue‐relevant senescence states across diverse biological systems. The full gene list and source annotations are provided in Table S2. In total, 210 senescence‐associated genes were included in the positive set.

To construct the negative gene set, we selected genes associated with cellular states that are functionally distinct from senescence. These included genes involved in cytoskeletal maintenance, housekeeping and ribosomal functions, cell cycle and proliferative signaling, programmed cell death (apoptosis), core metabolic processes, and receptor‐mediated signaling. These pathways collectively represent non‐senescent physiological or proliferative cellular conditions, providing a biologically meaningful contrast to senescence. In total, 148 non‐senescent genes were included in the negative set. The full list of negative genes and their functional annotations is provided in Table . S3

Protein sequences corresponding to all positive and negative gene sets were retrieved from the Universal Protein Knowledgebase (the UniProtKB database) (UniProt Consortium 2021). Because the predictive model operates on the protein sequence level, all annotated protein isoforms for each gene were collected. To reduce redundancy and avoid model bias caused by overrepresented sequence families, the initial dataset of 292,759 positive and 373,098 negative protein sequences was clustered using MMSeqs2 (Steinegger and Söding 2017) with an identity threshold of 0.3. This clustering procedure yielded a final representative dataset of 43,381 positive and 32,673 negative sequences. Only these representative sequences were used for subsequent model training and validation to avoid data leakage (Bernett et al. 2024). The dataset was divided into a training set (80%) and a validation set (20%) for model development and performance evaluation. The code and benchmark dataset are publicly available at https://github.com/HanliJiang13/SenSeqNet↗, ensuring transparency and reproducibility of the study.

To assess the generalizability of SenSeqNet, we performed external validation using two in vitro datasets (GSE13330 and GSE63577) and one in vivo dataset (SenCID; Tao et al. 2024). GSE13330 and GSE63577 profile human fibroblasts undergoing replicative or DNA damage–induced senescence, while SenCID provides experimentally verified senescence signatures from human tissues. In total, 26 senescence‐related genes were collected from these datasets, and their corresponding protein sequences were processed using the same UniProtKB retrieval. The validated gene names and their dataset sources are provided in Table S4.

As illustrated in Figure 2, SenSeqNet consists of two main components: ESM‐2 and the hybrid LSTM‐CNN architecture. Beginning with the first component, effective representation of sequence data is an essential step in building a robust prediction model, especially when dealing with protein sequences (LeCun et al. 2015). Thus, we employed the ESM‐2 (Lin et al. 2023), known for its exceptional ability to capture deep evolutionary relationships in protein sequences (Figure 2A). ESM‐2 is trained on approximately 65 million unique protein sequences using a masked language modeling objective, which involves predicting the identity of randomly selected amino acids within a sequence based on their surrounding context (Du et al. 2023). This training enables ESM‐2 to learn the dependencies between amino acids, allowing it to effectively capture both local and global contextual information.

The unsupervised nature of the ESM‐2 model's training allows it to internalize sequence patterns across diverse proteins, effectively capturing structural information linked to these patterns (Nguyen and Hy 2024). ESM‐2 is highly effective at processing and embedding from large‐scale protein datasets, making it particularly suitable for sequence representation tasks. The model's performance across various bioinformatics applications has been well‐documented (Mall et al. 2024), demonstrating its ability to discern subtle patterns within sequences, which are important for accurate prediction tasks. There are multiple ESM‐2 models, each differing in architectural configuration, where an increase in the number of transformer layers corresponds to a proportional increase in the number of model parameters. To identify the most effective variant for our dataset, we evaluated each variant's output and utilized the results to train two distinct models: a CNN to assess the spatial features and an LSTM to analyze the sequential dependencies within the results. By balancing the accuracy, computational efficiency, and scalability of these models, we determined the most effective ESM‐2 variant for our dataset. The selected model variant was then used for embedding generation throughout our study. The detailed results of this evaluation will be discussed in the Results section.

Followed by ESM‐2, the embedded sequences are input into the second main component, a hybrid LSTM‐CNN architecture (Figure 2B). In this setup, a Bidirectional LSTM network is first employed to process the sequence data, leveraging its capacity to model temporal dependencies in both directions, thus capturing comprehensive long‐term relationships between amino acids in the protein sequences. These sequential features are then passed to the CNN, which focuses on identifying spatial patterns within the LSTM‐extracted features. This combined approach enables the model to capture both temporal and spatial characteristics. This two‐part architecture, combining the evolutionary context captured by ESM‐2 with the temporal and spatial analysis capabilities of the LSTM‐CNN, enables SenSeqNet to effectively identify crucial features in protein sequences. This comprehensive approach significantly enhances the model's ability to distinguish senescence‐related proteins from non‐senescent proteins, as demonstrated by the results.

In this study, we explored several ML and DL architectures to develop a reliable model for detecting cell aging using protein sequences (LeCun et al. 2015). Given the lack of prior research in this specific domain of cell aging, our approach involved comprehensive experimentation with a range of models commonly employed in protein sequence prediction to understand their capabilities and limitations. The architectures evaluated include CNNs (Zhang et al. 2023), Recurrent Neural Networks (RNNs) (Müller et al. 2018), LSTMs (Tsukiyama et al. 2021), and bidirectional LSTMs. Recognizing the strengths of both CNNs in extracting spatial features and bidirectional LSTMs in capturing sequential dependencies (Khademi et al. 2022), we also experimented with hybrid models that combined these two architectures. The combined CNN and bidirectional LSTM model was designed to leverage the advantages of both, aiming to improve the accuracy and robustness of cell aging detection from protein sequences. Additionally, we constructed several ML models, including Bagging (Breiman 1996), Support Vector Machines (SVM) (Hearst et al. 1998), Random Forest (RF) (Breiman 2001), and XGBoost (Chen and Guestrin 2016), to provide a comparative understanding of their performance against DL approaches. Below, we present a comprehensive description of each DL model. Hyperparameters such as learning rate, dropout rate, batch size, number of layers, and hidden dimension were optimized for each model using Optuna (Akiba et al. 2019), a framework for hyperparameter optimization. To illustrate parameter dependence, we provide a representative sensitivity analysis for the LSTM variant of SenSeqNet (Figure S1), showing how variation in these hyperparameters influences model performance across the Optuna search.

The following commonly used evaluation metrics were employed to assess the effectiveness of the models (Ao et al. 2022): sensitivity (Sn), specificity (Sp), Matthews correlation coefficient (MCC) (Chicco and Jurman 2020), accuracy (Acc), F1 score, and AUC. The formulas for these metrics are provided below:Sn=TPTP+FNSp=TNTN+FPMCC=TP×TN−FP×FNTP+FPTP+FNTN+FPTN+FNAcc=TP+TNTP+FP+FN+TNF1score=2×TP2×TP+FP+FNwhere: True positive (TP): The number of correctly identified cases of cell aging; True negative (TN): The number of correctly identified non‐aging cases; False positive (FP): The number of non‐aging cases incorrectly identified as aging; False Negative (FN): The number of aging cases incorrectly identified as non‐aging; AUC represents the area under the receiver operating characteristic (ROC) curve and is used to evaluate the overall performance of the model across different thresholds.

We evaluated our model's performance at the gene level and conducted pathway analysis. The trained model generated predictions for all protein sequences, classifying each as senescence‐associated or not. To determine gene‐level performance, we aggregated predictions across all protein sequences belonging to each gene, calculating the proportion of correctly classified sequences as the per‐gene accuracy score. This approach accounted for genes with multiple protein products and provided a robust measure of classification performance. To identify a high‐confidence set of senescence‐associated genes, we applied stringent filtering criteria: (1) genes must be supported by at least 30 test sequences to ensure statistical reliability, and (2) genes must achieve a classification accuracy of at least 0.9. These thresholds were selected to balance between maintaining high confidence and retaining sufficient genes for downstream analysis. The resulting high‐confidence gene set was subjected to pathway enrichment analysis using g:Profiler web server (Reimand et al. 2007, 2016). Enrichment significance was assessed using the hypergeometric test with false discovery rate (FDR) correction for multiple testing. Pathways with FDR < 0.05 were considered significantly enriched. The gene ratio for each pathway was calculated as the proportion of high‐confidence genes present in the pathway relative to the total pathway size.

As detailed in the Methods section, we evaluated several variants of the Evolutionary Scale Modeling (ESM) for feature extraction from our protein sequence dataset. The goal was to identify the model that would provide the most accurate features while maintaining computational efficiency. Here, we present the results of this evaluation and the rationale for selecting the optimal ESM variant. The ESM variants evaluated included ESM‐2_t6, ESM‐2_t12, ESM‐2_t30, and ESM‐2_t33, each varying in complexity and embedding dimensions. The performance of these models was assessed based on their accuracy from our LSTM‐CNN architecture (Table 1).

The ESM‐2_t33 model achieved the highest accuracy at 86.43%, outperforming the other variants. This accuracy advantage makes ESM‐2_t33 the best choice for feature extraction. Therefore, we have selected ESM‐2_t33 for use throughout our experiments. To provide a more intuitive illustration of how ESM‐2_t33 embeddings differentiate the classes, we performed t‐SNE on the extracted features (Figure 3A). In contrast to the manual feature extraction methods—ACC (Figure 3B), CKSAAP (Figure 3C), and DDE (Figure 3D)—which exhibit heavily overlapping clusters, the ESM‐2_t33 embeddings form distinct groupings that closely align with the original class labels. This clear separation indicates that ESM‐2_t33 captures critical sequence‐level signals pertinent to senescence, rather than arbitrary noise. Moreover, the t‐SNE plots offer a complementary perspective to numerical metrics such as accuracy, specificity, and MCC, thereby reinforcing our confidence that ESM‐2_t33 provides a robust foundation for downstream classification tasks.

To assess the generalizability of SenSeqNet, we extracted features from multiple PLMs, including ProtTrans (T5), ProteinBERT, ESM‐1b, ESM‐1v, and ESM‐2, and subsequently fed these embeddings into our LSTM‐CNN architecture. Table 2 presents the performance metrics across all tested configurations.

Our results reveal a stark contrast between the ESM series and the other PLMs evaluated. ESM‐2 achieved the highest accuracy (86.43%), with strong sensitivity (90.51%), specificity (81.01%), and an F1‐score of 88.38%. ESM‐1b and ESM‐1v also exhibited robust performance, attaining accuracies of 84.57% and 83.48%, respectively. These findings suggest that ESM‐based embeddings capture key senescence‐related patterns in protein sequences, even when considering longer sequences exceeding 1024 residues.

ESM‐2 was selected for several key reasons. First, scalability and efficiency: ESM‐2 processes longer sequences (over 1024 residues) without requiring extensive chunking and merging, but feature extraction with ESM‐1b and ESM‐1v required approximately twice as much computational time compared to ESM‐2. Second, active development: As the most recent iteration of the ESM family, ESM‐2 benefits from ongoing updates, improving scalability and offering potential for further accuracy gains. Third, comparable high‐accuracy: ESM‐2 achieved slightly higher accuracy (86.43% vs. 84.57%), superior sensitivity (90.51%), and a robust Matthews correlation coefficient (MCC) of 0.7221, ensuring reliable predictions.

In contrast, ProtTrans (T5) and ProteinBERT, even when applying standard chunking and averaging strategies to handle sequences exceeding 1024 residues, performed significantly worse than ESM‐based models. ProtTrans (T5) achieved only 51.54% accuracy, 31.70% sensitivity, 68.71% specificity, an F1‐score of 37.77%, and an MCC of 0.0044, reflecting near‐random classification. ProteinBERT performed moderately better, with 60.94% accuracy, but still fell far short of the performance achieved by ESM‐based embeddings.

These results underscore the challenges faced by more generalized PLMs in distinguishing senescence‐associated proteins from non‐senescent proteins. In contrast, ESM‐2 provided robust performance, computational efficiency, and scalability, making it the optimal choice for SenSeqNet and well‐suited for large‐scale analyses of aging‐related proteins.

We compared the performance of SenSeqNet with several DL models, including CNN, RNN, LSTM, Bidirectional LSTMs, and various hybrid models. Each model was tested to determine its effectiveness in detecting cell aging from protein sequences, focusing on key performance metrics such as accuracy, AUC, Sn, Sp, and MCC using an independent test set. The comparative results of these models are summarized in Table 3.

The CNN, known for its ability to capture spatial features, achieved an accuracy of 80.38% and an MCC of 0.6002. While effective in identifying local patterns within the protein sequences, the CNN struggled to capture temporal dependencies critical for detecting cellular senescence. To address this limitation, RNNs were introduced to model sequential data. The RNN showed a small improvement, with an accuracy of 80.66% and an MCC of 0.6055, indicating its capability to better capture sequential patterns in the data. Next, more advanced temporal models, such as LSTMs, were explored. The LSTM, designed to handle long‐term dependencies in sequential data, further improved accuracy to 81.30%, with an MCC of 0.6165. Building on this, the Bidirectional LSTM, which processes information in both forward and backward directions, enhanced the ability to capture complex dependencies, achieving an accuracy of 81.89% and an MCC of 0.6343. This model excelled at processing temporal patterns but lacked focus on spatial feature extraction. We also evaluated several hybrid architectures. Among these, the CNN‐LSTM configuration demonstrated the lowest performance, with an accuracy of 78.02% and an MCC of 0.5570, suggesting that processing spatial features before sequential features may not fully leverage the strengths of both architectures. The pLSTM‐CNN model performed better, achieving an accuracy of 80.94% and an MCC of 0.6164, though it had higher computational demands due to the simultaneous processing of spatial and sequential features.

The best‐performing architecture was the SenSeqNet, which effectively combined the LSTM's ability to capture sequence dependencies with the CNN's proficiency in refining spatial features. SenSeqNet achieved the highest accuracy of 86.43%, an AUC of 0.94, and an MCC of 0.7221. These results demonstrate that the combination of LSTM and CNN architectures, when sequenced in this order, leads to a more comprehensive understanding of protein sequences, resulting in superior predictive accuracy. Overall, the results highlight that while single models such as CNNs and LSTMs are effective in their respective domains of spatial or temporal feature extraction, the hybrid approach of SenSeqNet provides the best of both worlds, outperforming both individual models and other hybrid configurations. The confusion matrix and ROC curve (Figure 3E,F) further confirm the stability and generalizability of SenSeqNet in detecting cellular senescence. To better understand the model's misclassifications, we conducted detailed gene‐level analysis (Table S5), which revealed that false positives predominantly occurred in specific gene families including ribosomal proteins, while true senescence markers showed consistently low error rates, providing biological insight into the model's decision boundaries.

To further evaluate the performance of SenSeqNet, we also compared its performance with various commonly used ML models including RF, XGBoost, SVM, and Bagging. The highest accuracy achieved among these machine models was 79.67% by RF, while XGBoost showed a notably high Sn of 94.19% but compromised Sp at 57.95%. Overall, none of the traditional models surpassed 80% in accuracy, with moderate MCC ranging between 0.5521 and 0.5844. The results are summarized in Table 4.

In comparison, SenSeqNet exhibited consistently superior performance, achieving the highest accuracy, specificity, F1‐score, and MCC, highlighting its effectiveness in extracting complex, nonlinear features from large‐scale protein sequence data. This disparity in performance underscores the capability of SenSeqNet to better handle the extensive and intricate relationships inherent in biological datasets, especially when dealing with large volumes of data, where traditional ML models may struggle to fully capture contextual dependencies.

To further validate SenSeqNet's predictive capability and simulate real‐world usage scenarios, we tested the model on an external dataset comprising 26 senescence‐related genes (Table ) that were not included in any phase of model training or development. This dataset yielded 52,227 protein sequences for validation. SenSeqNet classified 40,502 of these 52,227 sequences as senescence‐associated, achieving a true positive rate of 77.55% (Figure ). Detailed per‐sequence and per‐gene validation results are provided in Tables and, respectively, which include both binary predictions (positive/negative) and their associated confidence scores. The strong performance on these independent sequences confirms SenSeqNet's ability to generalize beyond its training data. This validation approach directly mirrors real‐world applications where researchers encounter newly discovered or poorly characterized proteins that require classification without prior experimental validation, demonstrating SenSeqNet's utility as a practical tool for identifying senescence‐associated sequences in exploratory research contexts. S4 S2 S6 S7

To evaluate the model's ability to identify senescence‐associated genes, we assessed per‐gene classification performance across our test set (Figure 4). The relationship between per‐gene accuracy and the number of supporting test sequences revealed that senescence‐related genes (positive class) consistently achieved high‐accuracy scores, with most clustering above 0.9 accuracy when supported by sufficient sequence data (Figure 4A). Notably, the model demonstrated exceptional performance on well‐validated senescence markers such as CKB, SYK, and DAO, all of which achieved near‐perfect accuracy despite having large numbers of test sequences. In contrast, certain ribosomal protein genes including RPL21, RPL26, and RPL28 showed markedly poor performance despite adequate sequence support, suggesting these genes may possess sequence features that confound the model's senescence‐specific patterns. This pattern was further confirmed by the distribution analysis, which showed that positive genes predominantly concentrated in the high‐accuracy range (> 0.9), while negative control genes exhibited a broader accuracy distribution, particularly when sequence support was limited (Figure 4B). To identify robust biological signals, we selected high‐confidence positive genes meeting stringent criteria (> 30 supporting sequences and > 0.9 accuracy), which yielded a set of reliably classified senescence‐associated genes. Pathway enrichment analysis of these high‐confidence genes revealed strong enrichment for established senescence‐related pathways, with the most significant associations found in fluid shear stress and atherosclerosis, senescence and autophagy in cancer, and negative regulation of cell population proliferation (Figure 4C). Notably, pathways related to positive regulation of cell population proliferation, cell migration, and response to oxygen levels all showed significant enrichment, consistent with known senescence biology where these processes are dysregulated.