The data for this study come from two parts: the immune repertoire data is obtained through the GSE158055 dataset in the GEO database. The single-cell sequencing data is sourced from the website http://covid19.cancer-pku.cn/#/dimensional-reduction, which provides a complete h5 file containing B cell, CD4T cell and CD8T cell related data. Following data collection, the two datasets were matched, and we specially screened for B cells (BCRH_cdr3aa) and T cells (TCRA_cdr3aa) that include the cdr3 amino acid sequences for further analysis.

We used the Scanpy library in Python for single-cell sequencing data analysis (13). The analysis proceed as follows: firstly, we performed K-Nearest Neighbor(KNN) clustering separately on the selected B cells and T cells to identify distinct cell subpopulations. After clustering, we identified the marker genes of each cell subpopulation based on the significant expression level. Using the identified marker genes, we employed manual annotation methods, supported by relevant literature, to classify the cell subpopulations. Finally, we visualized the cell annotation results using UMAP diagrams, which effectively illustrate the distribution and characteristics of the annotated cell subpopulations.

The kmers algorithm segments each specific cdr3 amino acid sequence based on selected lengths. In this study, we chose k values of 2, 3 and 4. Based on these lengths, we combined the 20 common amino acid types to generate 400 amino acid concatemers for k=2, 8000 amino acid concatemers for k=3, and 160,000 amino acid concatemers for k=4. Next, we constructed a count matrix for CDR3 amino acid concatemers according to the different states of the same cell type, assembling the dataset based on the various k values. This process resulted in seven datasets containing CDR3 amino acid concatemers. During the construction of the count matrix, any concatemers that did not appear were filled with 0 to ensure accurate counting.

Given the large number of features in the dataset, we first counted the number of filled 0 values in each column and set a threshold to eliminate features containing outliers.

The threshold filtering mainly considered the expression of feature in the sample, which determined whether to retain the features. Next, we employed analysis of variance (ANOVA) for feature dimensionality reduction, retaining only those features with a significance level lower than 0.05. After dimensionality reduction, we divided each dataset into a training dataset and a test dataset in a 4:1, using seven different classification algorithms to classify and compare each cell type based on structural similarity. The algorithms included: LR, SVM, MLP, KNN, DC, RF and XGB. In the final presentation, we selected the classification algorithm that exhibited the best performance. Furthermore, the grid search algorithm was applied to introduce hyperparameters to filter out the best classification model, with a p-value range of 0.001-0.049. Comparing multiple classification models through changes in p-value to find the best combination and improve the generalization ability. To address the imbalance in data distribution, we utilized the ROC curve and the average AUC value as evaluation indicators for the classification results. The advantage of ROC curve was that it was intensive to data imbalance and effectively demonstrate the classification performance across the entire threshold range. The average macro AUC was calculated independently for each category, avoiding the influence of large categories and not relying on the number of categories (14, 15).

To explore the connections between cell subpopulations, we measured the interaction between them based on the membrane protein pairs present on the surface of the cells in different states (16). This assessment provided insights into the interactions among cell subpopulations. For calculating the metabolic pathway scores of different cell subpopulations, we accessed the KEGG database (KEGG, https://www.kegg.jp/) and utilized the KEGGREST package in R to screen 84 human-related metabolic pathways (17). We then combined the expression levels of metabolic genes associated with each pathway to compute a score for each pathway in various states. Additionally, we conducted significance testing for each metabolic pathway to highlight the differences among cell subpopulations under different states. The t test method was used for pairwise comparison, while the Kruskal Wallis test method was used for multiple comparison. This approach allowed us to identify key metabolic pathways that may be relevant to the functional characteristics of the cell subpopulations in response to different conditions.

We propose a novel framework that integrates k-mer cutting, dataset screening, dimensionality reduction and classification algorithms for a joint analysis of immune repertoire data and single-cell sequencing data. This framework achieves heterogeneous classification of immune cells based on the structural similarity of CDR3 amino acid sequences, providing new insights for personalized disease treatment. The framework consists of five main parts: the first part realizes cell clustering and cell annotation through single-cell sequencing data, and then uses the immune repertoire data to extract the cdr3 amino acid sequences of the same cell type under each cell type and construct the corresponding dataset. The second part is to perform kmers processing on the cdr3 amino acid sequences in the dataset to cut out amino acid concatemers of different kmer lengths. The third part is to generate the corresponding count matrix for the amino acid concatemers according to different cutting lengths. The fourth part is to reduce the dimension of the count matrix of the amino acid concatemers using feature selection methods. The fifth part is to classify and compare the count matrix after dimensionality reduction using seven different classification algorithms to achieve heterogeneous classification of immune cell subsets. A schematic diagram of the framework is provided in Figure 1 , illustrating the workflow and interconnections among these components.

We investigated the shared cells between the immune repertoire data and single-cell sequencing data, ultimately identifying 278,298 B cells and 220,968 T cells, as illustrated in Supplementary Figure 1 . Using single-cell sequencing data, we generated B cell and T cell profiles of COVID-19 patients, as shown in Figure 2A and Figure 2B . For B cells, we identified three cell subpopulations, with the marker genes of each cell subpopulation displayed in Supplementary Figures 2A–C . We found that the proportion of B_c1_TCL1A cell subpopulation increased significantly after COVID-19 infection, and the proportion of B_c3_TNFRSF1B subpopulation cells decreased significantly, as shown in Figure 2C . For T cells, we identified eight cell subpopulations, and the marker genes for each cell subpopulation are shown in Supplementary Figures 2D–K . We found that the proportion of CD4T_c1_FOS and CD8T_c4_LDHA cell subpopulation increased after COVID-19 infection. The proportion of CD4T_c4_CCR7 and CD8T_c1_CCL5 cell subpopulation shown a notable decrease, as shown in Figure 2D . We examined the inflammatory factors related to immune cells and discovered that the expression of MT-CO2 and MT-ND3 genes was significantly elevated in both B cell and T cell subsets following COVID-19 infection. The B_c3_TNFRSF1B populations shown a significant increase in the expression of NFKB1, NFKB2, ID3, SKIL and TNF genes. Conversely, B_c1_TCL1A, CD4T_c1FOS and CD8T_c3_GZMK subpopulations exhibited overexpression of DUSP1 and JUN genes. This result displayed in Figures 2E, F .

To characterize the structure of cdr3 amino acid sequences among different immune cell subsets, we applied the KDDC framework to achieve heterogeneous classification of immune cells. We constructed a count matrix using cdr3 sequence amino acid concatemers of varying lengths and compared the classification results to reflect the heterogeneity of different immune cell subsets. Using variance analysis, we performed dimensionality reduction on the count matrix and evaluated several classification algorithms for comparison. ultimately, we selected the most suitable classification algorithm for our study. The classification results of each algorithm are shown in Table 1 . Our findings indicate that the classification performance for the three B cell subsets was notably strong. Specifically, the average AUC value for the B_c1_TCL1A cell subset using kmer3 reached 96.5%, as shown in Figure 3A . The average AUC value for the B_c2_CD27 cell subset across kmer2, kmer3 and kmer4 reached an impressive 99%, depicted in Figure 3B . Additionally, the average AUC value for the B_c3_TNFRSF1B cell subset with kmer2, kmer3 and kmer4 reached 93.2%, illustrated in Figure 3C . The relevant parameters of the best model can refer to Supplementary Table 1 . Afterwards, we conducted an interactive analysis on the three kmers datasets that demonstrated the best classification results. We discovered share amino acid concatemers among the kmers dataset, which reflect the structural similiarities with cdr3 sequences. The structure of DLW and FDL were high importance in B_c1_TCL1A and B_c2_CD27. Meanwhile, the structure of YDY and RGGF were great importance in B_c2_CD27 and B_c3_TNFRSF1B, as displayed in Figure 3D .

The cell subsets within T cells were classified and compared using the KDDC framework, with the results of different classification algorithms summarized in Tables 2 and 3 . We discovered that the CD4T_c2_AQP3 subset had a strong classification performance in kmer4, achieving an average AUC value of 93.5%, as illustrated in Figure 4A . The CD4T_c3_LEF1 subset demonstrated effective classification with kmer4, yielding an average AUC value of 92.7%, as shown in Figure 4B . The classification performance for the CD8T_c2_CD8B subset was optimized using kmer2, kmer3 and kmer4, attaining an impressive average AUC value of 99%, represented in Figure 4C . The CD8T_c3_GZMK subset exhibited a better classification performance in kmer2 and kmer3, resulted in an average AUC value of 91.2%, depicted in Figure 4D .

We investigated the interactions among 11 immune cell subsets by selecting immune cell surface protein pairs. This approach allowed us to evaluate the connectivity of the immune cell subsets in COVID-19 patients under different infection states, based on the affinity of each pair of proteins and the expression of the encoding genes. When comparing the post-COVID-19 infection state to the normal state, we observed a weakened binding affinity between the CD4T_c4_CCR7 and CD8T_c2_CD8B subset. Conversely, the binding strength of CD8T_c1_CCL5 to both CD8T_c4_LDHA and CD4T_c1_FOS was found to be enhanced. In comparisons between the mid-term and severe states following COVID-19 infection, we noted a significant reduction in the binding affinity between the B_c1_TCL1A and CD4T_c1_FOS subset, as illustrated in Figures 5A–C .

Furthermore, we explored the specific metabolic pathways associated with the 11 immune cell subpopulations. By comparing the post-COVID-19 infection state with the normal state, we observed that these cell subpopulations exhibited active engagement in various metabolic pathways, revealing metabolic heterogeneity across different infection states, as shown in Figure 6A . Meanwhile, we also conducted significant hypothesis tests on the metabolic pathways involved, as shown in the Supplementary Figure 3 . The Glycerophospholipid metabolism, Oxidative phosphorylation, and Purine metabolism metabolic pathways were significantly expressed in B_c1_TCL1A, B_c2_CD27, and CD8T_c1_CCL5 cell subpopulations, particularly notable during the mid-stage of infection, as illustrated in Figure 6B . The CD4T_c4_CCR7 subpopulation shown significant changes in the scores of Pyrimidine metabolism, Nicotinate and nicotinamide metabolism metabolic pathways post-COVID-19 infection, as followed in Figure 6C . Additionally, the metabolic scores for the B_c3_TNFRSF1B, CD4T_c2_AQP3, CD4T_c3_LEF1, CD8T_c2_CD8B, and CD8T_c3_GZMK subpopulations displayed a trend of increasing followed by a decrease, as summarized in Figure 6D . The CD8T_C4_LDHA subpopulation exhibited notable changes in the scores for the Starch and sucrose metabolism metabolic pathways following COVID-19 infection, as summarized in Figure 6E .

Finally, we employed the KDDC algorithm to classify COVID-19 patients in different states based on BCR and TCR sequencing data. We first analyzed the distribution of cdr3aa sequence lengths in both BCR and TCR across different infection states. The cdr3aa length distribution exhibited significant differences post-COIVD-19 infection, as depicted in Figures 7A , B . The KDDC framework was applied to compare datasets with different kmer lengths. We identified the datasets and classification algorithms yielded the best classification performance, as summarized in Table 4 . The classification performance of COVID-19 patients based on BCR sequencing data was exceptional, achieving an average AUC value of 98.8% in kmer3, as shown in Figure 7C . In contrast, the classification performance of COVID-19 patients based on TCR sequencing data was good, with an average AUC value of 83.3% in kmer4, as followed in Figure 7D .