Machine-Learning Analysis of Serum Proteomics in Neuropathic Pain after Nerve Injury in Breast Cancer Surgery Points at Chemokine Signaling via SIRT2 Regulation

Background: Persistent postsurgical neuropathic pain (PPSNP) can occur after intraoperative damage to somatosensory nerves, with a prevalence of 29–57% in breast cancer surgery. Proteomics is an active research field in neuropathic pain and the first results support its utility for establishing diagnoses or finding therapy strategies. Methods: 57 women (30 non-PPSNP/27 PPSNP) who had experienced a surgeon-verified intercostobrachial nerve injury during breast cancer surgery, were examined for patterns in 74 serum proteomic markers that allowed discrimination between subgroups with or without PPSNP. Serum samples were obtained both before and after surgery. Results: Unsupervised data analyses, including principal component analysis and self-organizing maps of artificial neurons, revealed patterns that supported a data structure consistent with pain-related subgroup (non-PPSPN vs. PPSNP) separation. Subsequent supervised machine learning-based analyses revealed 19 proteins (CD244, SIRT2, CCL28, CXCL9, CCL20, CCL3, IL.10RA, MCP.1, TRAIL, CCL25, IL10, uPA, CCL4, DNER, STAMPB, CCL23, CST5, CCL11, FGF.23) that were informative for subgroup separation. In cross-validated training and testing of six different machine-learned algorithms, subgroup assignment was significantly better than chance, whereas this was not possible when training the algorithms with randomly permuted data or with the protein markers not selected. In particular, sirtuin 2 emerged as a key protein, presenting both before and after breast cancer treatments in the PPSNP compared with the non-PPSNP subgroup. Conclusions: The identified proteins play important roles in immune processes such as cell migration, chemotaxis, and cytokine-signaling. They also have considerable overlap with currently known targets of approved or investigational drugs. Taken together, several lines of unsupervised and supervised analyses pointed to structures in serum proteomics data, obtained before and after breast cancer surgery, that relate to neuroinflammatory processes associated with the development of neuropathic pain after an intraoperative nerve lesion.


Introduction
Persistent postsurgical neuropathic pain (PPSNP), defined as pain caused by a lesion of the somatosensory system associated with the surgical procedure [1], poses clinical challenges due to its intensity, relative resistance to current pharmacologic treatments, and 3 p = 0.895) and an almost linear placement of the quantiles in the QQ plot ( Figure 1D). Using the first Bayesian decision limit at x-position 0.18, the resulting two groups of n = 60 and n = 54 samples significantly overlapped with the predefined subgroup structure of non-PPSNP versus PPSNP ( Figure 1E; Fisher's exact test: odds ratio 4.28, 95% confidence interval, CI: 1.384-7.367, p = 0.00468). Consideration of the second Bayesian decision boundary did not yield further significant results and was therefore abandoned.  [26,27], performed during "PC-corr" analysis [33] while attempting group segregation based on the respective PC.  [26,27], performed during "PC-corr" analysis [33] while attempting group segregation based on the respective PC. (B): Bar chart of the loadings of protein markers on PC1, sorted in descending order of magnitude. The proteins are named as in the Proseek panel for consistency. Please refer to Table 1 for standard protein names. (C): Distribution of the patients' individual scores on PC1, described by the Pareto density estimation (PDE) [34], to which a Gaussian mixture model (GMM) with M = 3 modes was fitted. The Bayesian boundaries between the modes are indicated as dashed magenta perpendicular lines. The first boundary at x-position 0.18 provided a suitable GMM based grouping criterion of data set instances as shown in Panel E. (D): Quantile-quantile (QQ) plot of the theoretical and observed quantiles of the data, with line of identity. (E): Heatmap with the original subgroup structure (non-PPSNP versus "PPSNP") and a subgroup structure that resulted from the GMM analysis of the coordinates of the projected samples on PC1 (see panel C). The color scheme green/blue of column 1 repeats that used in panel A for non-PPSNP versus "PPSNP ". The darker red color in columns 2 and 3 indicate data set instances belong to data set instances in Gaussian #1 of panel C, whereas the lighter orange color denotes data belonging to the second and third Gaussians combined in panel C. The GMM-based grouping significantly overlapped with the prior non-PPSNP versus PPSNP group structure (Fisher's exact test [35]: p = 0.00468). The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the library "ggplot2" (https://cran.r-project.org/package=ggplot2 (accessed on 14 March 2022) [37]).

Data Projection-Based Protein Marker Patterns Relevant to Pain-Related Subgroup Separation
The results of PC-Corr analyses indicated that the non-PPSNP versus PPSNP subgroups were best separated when the entire data set recorded before and after surgery was projected onto a lower dimensional plane after probe-level quantile normalization [32] and centering of the data. Significant segregation of the two neuropathic-pain related subgroups was already observed along the first dimension of the PCA projection of the whole data set as mentioned above (PC1; Wilcoxon-Man-Whitney U test p-value < 0.05, AUC-ROC = 0.63, AUC-PR = 0.58), which explained 19.6% of the total variance in the proteomics data ( Figure 1A). The protein with the largest contribution to PC1 was SIRT2 ( Figure 1B). The distribution of the coordinates of the projection of the observations on PC1 was best described by a trimodal Gaussian mixture ( Figure 1C), which showed no significant difference between the fitted and observed distributions (Kolmogorov-Smirnov test: p = 0.895) and an almost linear placement of the quantiles in the QQ plot ( Figure 1D). Using the first Bayesian decision limit at x-position 0.18, the resulting two groups of n = 60 and n = 54 samples significantly overlapped with the predefined subgroup structure of non-PPSNP versus PPSNP ( Figure 1E; Fisher's exact test: odds ratio 4.28, 95% confidence interval, CI: 1.384-7.367, p = 0.00468). Consideration of the second Bayesian decision boundary did not yield further significant results and was therefore abandoned.
A subgroup structure as observed in the PCA-based projection of the proteomics data was further supported by an alternative projection on a trained emergent self-organizing feature map (ESOM). Large U-heights ( Figure 2A) forming a "mountain ridge" separated a small region of 13 data points from the larger region of 101 samples, which indicated the emergence of two main clusters in the data. This agreed with the prior "non-PPSNP" and "PPSNP" group structure (Fisher's exact test: odds ratio 4.26 (95% CI 1.017-25.55, p = 0.03649; Figure 2B). The separated subgroup was smaller than in the analogous PCA-based result; however, all contained data instances also in the smaller subgroups were separated from the majority on the PCA projection ( Figure 2C).  [34], to which a Gaussian mixture model (GMM) with M = 3 modes was fitted. The Bayesian boundaries between the modes are indicated as dashed magenta perpendicular lines. The first boundary at x-position 0.18 provided a suitable GMM based grouping criterion of data set instances as shown in Panel E. (D): Quantile-quantile (QQ) plot of the theoretical and observed quantiles of the data, with line of identity. (E): Heatmap with the original subgroup structure (non-PPSNP versus "PPSNP") and a subgroup structure that resulted from the GMM analysis of the coordinates of the projected samples on PC1 (see panel C). The color scheme green/blue of column 1 repeats that used in panel A for non-PPSNP versus "PPSNP ". The darker red color in columns 2 and 3 indicate data set instances belong to data set instances in Gaussian #1 of panel C, whereas the lighter orange color denotes data belonging to the second and third Gaussians combined in panel C. The GMMbased grouping significantly overlapped with the prior non-PPSNP versus PPSNP group structure (Fisher's exact test [35]: p = 0.00468). The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the library "ggplot2" (https://cran.r-project.org/package=ggplot2 (accessed on 14 March 2022) [37]).
A subgroup structure as observed in the PCA-based projection of the proteomics data was further supported by an alternative projection on a trained emergent self-organizing feature map (ESOM). Large U-heights ( Figure 2A) forming a "mountain ridge" separated a small region of 13 data points from the larger region of 101 samples, which indicated the emergence of two main clusters in the data. This agreed with the prior "non-PPSNP" and "PPSNP" group structure (Fisher's exact test: odds ratio 4.26 (95% CI 1.017-25.55, p = 0.03649; Figure 2B). The separated subgroup was smaller than in the analogous PCAbased result; however, all contained data instances also in the smaller subgroups were separated from the majority on the PCA projection ( Figure 2C). Results of projection of the data, after probe-level quantile normalization and pooled first and second samples, onto an emergent self-organizing map (ESOM; for further details of this artificial neuronal network-based data projection method, see [38,39]). (A): Three-dimensional U-matrix visualization of distance-based structures of the serum concentration of d = 74 proteomic markers following projection of the data points onto a toroid grid of 4000 artificial neurons where opposite Figure 2. Results of projection of the data, after probe-level quantile normalization and pooled first and second samples, onto an emergent self-organizing map (ESOM; for further details of this artificial neuronal network-based data projection method, see [38,39]). (A): Three-dimensional U-matrix visualization of distance-based structures of the serum concentration of d = 74 proteomic markers following projection of the data points onto a toroid grid of 4000 artificial neurons where opposite edges are connected. The dots represent the so-called "best matching units" (BMU), i.e., neurons on the grid that, after ESOM learning, carried a data vector that was most similar to a subjects' data vector.
Please note that one BMU can carry vectors of several cases, i.e., the number of BMUs is not necessarily equal to the number of cases. The U-matrix visualization was colored as a top view of a topographic map with brown (up to snow-covered) heights and green valleys with blue lakes. Watersheds indicate borderlines between different clusters. Two clusters emerged in this way, separated by the white "mountain ridge" at the left of the U-matrix. BMUs belonging to clusters #1 or #2 are colored in green or bluish, respectively. (B): Mosaic plot, visualizing the contingency table between the original group structure and the cluster identified on the U-matrix. The p value of 0.03649 denotes the results of a Fisher's exact test [35]. (C): Heatmap with the original subgroup structure (non-PPSNP versus "PPSNP") and a subgroup structure that resulted from the U-matrix shown in in Panel A. The clusters based on the U-matrix (Panel A) are shown in the 2nd and 3rd column. For comparison, the PCA-based clusters ( Figure 1D) are displayed in the last two columns The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the libraries "ggplot2" (https://cran.r-project.org/package= ggplot2 (accessed on 14 March 2022) [37]), "ggmosaic" (https://cran.r-project.org/package=ggmosaic (accessed on 14 March 2022) [40]) and "Umatrix" (https://cran.r-project.org/package=Umatrix (accessed on 14 March 2022) [41]).

Supervised Machine Learning-Based Identification and Evaluation of Proteomic Markers Informative for Pain-Related Subgroup Segregation
Training the classifiers with all d = 74 proteins included in this analysis was successful in logistic regression, support vector machine, k-nearest neighbors, and random forests, which were able to identify whether an instance of the dataset was acquired from a patient in the non-PPSNP or PPSNP subgroup ( Figure 3A and Table 2). After feature selection using the Boruta method ( Figure 4), d = 19 proteomic markers remained (Table 3). Training the classifiers with these d = 19 markers resulted in better classification performance than with all 74 markers, which is a typical observation in machine learning, where eliminating noise is often rewarded with better results. Now, all classifiers appeared to perform better than change in assigning a sample to the correct neuropathic pain subgroup. In contrast, when using permuted features or the d = 45 proteomic markers of ABC set "C," i.e., the least important items, all classifiers resorted to random class assignment, indicating that (i) the successful classification results were unlikely to be due to overfitting and (ii) the item categorization captured the relevant items ( Figure 3A). edges are connected. The dots represent the so-called "best matching units" (BMU), i.e., neurons on the grid that, after ESOM learning, carried a data vector that was most similar to a subjects' data vector. Please note that one BMU can carry vectors of several cases, i.e., the number of BMUs is not necessarily equal to the number of cases. The U-matrix visualization was colored as a top view of a topographic map with brown (up to snow-covered) heights and green valleys with blue lakes. Watersheds indicate borderlines between different clusters. Two clusters emerged in this way, separated by the white "mountain ridge" at the left of the U-matrix. BMUs belonging to clusters #1 or #2 are colored in green or bluish, respectively. (B): Mosaic plot, visualizing the contingency table between the original group structure and the cluster identified on the U-matrix. The p value of 0.03649 denotes the results of a Fisher's exact test [35]. (C): Heatmap with the original subgroup structure (non-PPSNP versus "PPSNP") and a subgroup structure that resulted from the U-matrix shown in in Panel A. The clusters based on the U-matrix (Panel A) are shown in the 2nd and 3rd column. For comparison, the PCA-based clusters ( Figure 1D) are displayed in the last two columns The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the libraries "ggplot2" (https://cran.r-project.org/pack-age=ggplot2 (accessed on 14 March 2022) [37]), "ggmosaic" (https://cran.r-project.org/pack-age=ggmosaic (accessed on 14 March 2022) [40]) and "Umatrix" (https://cran.r-project.org/pack-age=Umatrix (accessed on 14 March 2022) [41]).

Supervised Machine Learning-Based Identification and Evaluation of Proteomic Markers Informative for Pain-Related Subgroup Segregation
Training the classifiers with all d = 74 proteins included in this analysis was successful in logistic regression, support vector machine, k-nearest neighbors, and random forests, which were able to identify whether an instance of the dataset was acquired from a patient in the non-PPSNP or PPSNP subgroup ( Figure 3A and Table 3). After feature selection using the Boruta method ( Figure 4), d = 19 proteomic markers remained ( Table 2). Training the classifiers with these d = 19 markers resulted in better classification performance than with all 74 markers, which is a typical observation in machine learning, where eliminating noise is often rewarded with better results. Now, all classifiers appeared to perform better than change in assigning a sample to the correct neuropathic pain subgroup. In contrast, when using permuted features or the d = 45 proteomic markers of ABC set "C," i.e., the least important items, all classifiers resorted to random class assignment, indicating that (i) the successful classification results were unlikely to be due to overfitting and (ii) the item categorization captured the relevant items ( Figure 3A).  Results of supervised analyses of the possibility to train machine-learning algorithms with the information of selected proteomic markers to enable them to correctly assign a patient to the subgroup with nerve injury but without neuropathic pain (non-PPSNP) or to the subgroup with nerve injury and neuropathic pain ("PPSNP"). (A): Boxplots of the obtained balanced classification accuracy by different types of machine learning algorithms in assigning sub-jects to the subgroups when training was done with all protein markers or with the markers identified as the most informative in four consecutive item categorization techniques implemented as computed ABC analyses (for the protein markers identified as important, see Table 3). In case the selected proteins carried relevant information for patient subgroup assignment, the classification accuracy should be better than guessing. For comparison, the balanced classification accuracy achieved with permuted characteristics is shown, as well as the balanced classification (balanced) accuracy obtained when using the items placed by the first ABC analysis in subset "C", which captures the least relevant items of a set. The expectations here were that without overfitting the classification (balanced) accuracy should not be better than guessing. The boxes have been constructed using the minimum, quartiles, median (solid line within the box), and maximum. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. (B): Results of the consecutive ABC analysis of the importance of protein markers. In the first ABC analysis, the counts were entered at which each maker occurred among the selected features in 1000 Boruta feature selection analyses on randomly drawn 2/3 of the data sets. In the subsequent ABC analyses, only the counts of occurrence of markers placed in ABC subset A by the previous ABC analysis were entered. The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R packages "ggplot2" (https://cran.r-project.org/package=ggplot2 (accessed on 14 March 2022)) and (https://cran.r-project.org/package=ABCanalysis (accessed on 14 March 2022) [42] accuracy by different types of machine learning algorithms in assigning sub-jects to the subgroups when training was done with all protein markers or with the markers identified as the most informative in four consecutive item categorization techniques implemented as computed ABC analyses (for the protein markers identified as important, see Table 2). In case the selected proteins carried relevant information for patient subgroup assignment, the classification accuracy should be better than guessing. For comparison, the balanced classification accuracy achieved with permuted characteristics is shown, as well as the balanced classification (balanced) accuracy obtained when using the items placed by the first ABC analysis in subset "C", which captures the least relevant items of a set. The expectations here were that without overfitting the classification (balanced) accuracy should not be better than guessing. The boxes have been constructed using the minimum, quartiles, median (solid line within the box), and maximum. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. (B): Results of the consecutive ABC analysis of the importance of protein markers. In the first ABC analysis, the counts were entered at which each maker occurred among the selected features in 1000 Boruta feature selection analyses on randomly drawn 2/3 of the data sets. In the subsequent ABC analyses, only the counts of occurrence of markers placed in ABC subset A by the previous ABC analysis were entered. The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R packages "ggplot2" (https://cran.r-project.org/package=ggplot2 (accessed on 14 March 2022)) and (https://cran.r-project.org/package=ABCanalysis (accessed on 14 March 2022) [42]).

Figure 4.
Example output of the importance analysis of protein markers for the allocation of patient subgroups ("non-PPSNP" versus "PPSNP") according to an analysis based on random forests ("Boruta" [43]). The proteins are named as in the Proseek panel for consistency. Please refer to Table 1 for standard protein names. The importance measure of a feature (here: of the protein markers) results from the decrease in classification accuracy due to the random permutation of feature values. It is calculated separately for all trees in the forest that use the respective feature for classification. Then the mean value and the standard deviation of the loss of accuracy are calculated. The z-score is used in comparison to an external reference, the so-called "shadow" features, which is obtained by permuting the values of the original feature. The boxes were constructed using the minimum, quartiles, median (solid line inside the box) and maximum of these values. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. The black circles indicate outliers from this interval. The green and orange boxes represent "confirmed" or tentatively significant features, respectively, i.e., features that contribute to the classification success. The red boxes are confirmed as non-informative in order to be excluded from further analysis. The empty boxes are the above-mentioned "shadow" features. The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R library "Boruta" (https://cran.r-project.org/package=Boruta (accessed on 14 March 2022) [43]).

Figure 4.
Example output of the importance analysis of protein markers for the allocation of patient subgroups ("non-PPSNP" versus "PPSNP") according to an analysis based on random forests ("Boruta" [43]). The proteins are named as in the Proseek panel for consistency. Please refer to Table 1 for standard protein names. The importance measure of a feature (here: of the protein markers) results from the decrease in classification accuracy due to the random permutation of feature values. It is calculated separately for all trees in the forest that use the respective feature for classification. Then the mean value and the standard deviation of the loss of accuracy are calculated. The z-score is used in comparison to an external reference, the so-called "shadow" features, which is obtained by permuting the values of the original feature. The boxes were constructed using the minimum, quartiles, median (solid line inside the box) and maximum of these values. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. The black circles indicate outliers from this interval. The green and orange boxes represent "confirmed" or tentatively significant features, respectively, i.e., features that contribute to the classification success.
The red boxes are confirmed as non-informative in order to be excluded from further analysis. The empty boxes are the above-mentioned "shadow" features. The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R library "Boruta" (https://cran.r-project.org/package= Boruta (accessed on 14 March 2022) [43]).
Finally, to further narrow the focus on the most relevant proteomic markers, ABC analysis was performed in three further nested steps, whereby the feature set was successively reduced to d = 9, 4, and finally d = 2 protein markers (Table 3). This procedure can be repeated until the ABC curve ( Figure 3B) touches the curve of uniform distribution of feature importance, since this curve marks the condition in which all features had the same chance to contribute to the subgroup separation, from which no particularly important feature can be separated any more. This procedure gradually reduced the classification power, but even with only CD244 and SIRT classification was still better than random assignment for logistic regression, support vector machine and random forests ( Figure 3A). Of note, the observation of SIRT2 as the most prominent marker was consistent with its importance in the PCA projection on the most relevant PC1. Table 2. Performance measures for the correct assignment of patients to the subgroup with nerve injury but without neuropathic pain (non-PPSNP) or to the subgroup with nerve injury and neuropathic pain ("PPSNP"). The performance of machine learning-based random forests classifiers is given; for further algorithms the key data (balanced accuracies) are shown in Figure 3. Classification performance was calculated (i) when training the algorithm with all protein markers or (ii-v) with the markers identified as the most informative in four consecutive item categorization techniques implemented as computed ABC analyses ("reduced data set #2-4; for the protein markers identified as important, see Table 3). For comparison, (vi) the balanced classification accuracy achieved with (permuted characteristics is shown, as well as the balanced classification accuracy obtained when using the items placed by the first ABC analysis in subset "C", which captures the least relevant items of a set. For the protein markers identified as important, see Table 3). For comparison, the balanced classification accuracy achieved with permuted characteristics is shown, as well as the balanced classification accuracy obtained when using the items placed by the first ABC analysis in subset "C", which captures the least relevant items of a set.  Table 3. Details of the d = 19 proteins selected in a first computed ABC analysis that evaluated the counts at which each protein was among the selected features in 1000 Boruta feature selection analyses ( Figure 4) on randomly drawn 2/3 of the data sets, aimed to identify the most relevant proteomic markers for assigning a patient to the subgroup with nerve injury but no neuropathic pain (non-PPSNP) or to the subgroup with nerve injury and neuropathic pain ("PPSNP"). The frequency occurrence in the set of selected features in the Boruta analysis is given in descending order. The p-values of group differences, calculated in the raw untransformed data, are the result of Mann-Whitney U tests [26,27], whereas the effect sizes of the group differences, quantified as Cohen's d [44]. P-values in bold letters indicate significant effects for better visibility. Positive values indicate that the protein marker was observed at higher concentrations in the patients with neuropathic pain "(PPSNP"). The four consecutive ABC analyses reduced the feature set from the initial d = 19 proteins (all table) to finally d = 2 proteins (top two proteomic markers). The proteins are named as in the Proseek panel for consistency. Please refer to Table 1

Discussion
The PPSNP and non-PPSNP subgroups showed different proteomics patterns when classical and machine learning-based feature selection techniques were used to identify the most informative proteins distinguishing these groups. The protein patterns already differed between the groups before nerve injury, whereas there was no clear difference when the proteins were compared before and after nerve injury. Thus, these distinct preinjury protein patterns could reflect protective or predisposing factors associating with the development of PPSNP. The results of these analyses included 19 different serum protein makers from a candidate panel of 74 markers that could eventually be narrowed down to only two proteins with sitruin2 (SIRT2) as a possible predisposing protein for PPSNP. The present analyses were performed in the context of a concerted AI interpretation between data science and biomedical experts, as recently described [45], and conceptually similar to a conversational machine learning approach also recently presented [46], i.e., the results are facilitated by collaboration between different disciplines. Possible biomedical interpretations of the results are outlined below.
The NAD-dependent deacetylase sirtuin 2 (SIRT2) was identified as the most informative protein marker to train machine-learning algorithms to identify samples with neuropathic pain. SIRT2 is a class III histone deacetylase expressed ubiquitously, but more abundantly in the central nervous system than in other tissues [47]. It plays a role in microtubule acetylation and myelination [48], and it is involved in the suppression of NFkB-related inflammatory processes [49][50][51][52]. It is also involved in the regulation of neuroinflammatory processes via activation of microglia [53], which plays an important role in the response to peripheral nerve injury [54] and synaptic plasticity in persistent pain [55]. Another link to persistent pain arises from the role of SIRT2 in learning and memory, which are biological processes in terms of the Gene Ontology (GO) knowledgebase [56] and have emerged as key features of persistent pain in a computational functional genomics analysis [57]. SIRT2 is also involved in cancer where it has been proposed as both a tumor suppressor and tumor promoter [58]. However, its role as a tumor suppressor seems to be more frequently highlighted [59,60], and also in breast cancer [61]. It is also considered as a target for drugs against age-related and/or neurodegenerative disorders [62] and also for cancer [63].
A role of SIRT2 in neuropathic pain has been highlighted in a mouse model of cisplatininduced peripheral neuropathy (CIPN) [64]. In humans, CSF-levels of SIRT2 were also among the protein markers relevant to persistent pain. Painful knee osteoarthritis has been patho-physiologically associated with neuroinflammatory processes and neuroimmune cross-links between the periphery and CNS. The CSF levels of SIRT2 were almost two-fold higher in the knee osteoarthritis patients than in healthy controls (See Table 3 in [65]). In the serum SIRT2 levels, however, there was no difference between the groups. In the present proteomics samples, the serum SIRT2 levels were higher in patients who developed neuropathic pain compared with those who had neuropathy without pain (Table 1). A brief review of what is known about SIRT2 in pain did not provide a clear direction of change. The cited results [65] in humans might be related to a pathology other than nerve lesion after surgery, whereas inflammation in arthritis and neuroinflammation in persistent pain represent a common mechanism. On the other hand, the rodent results are closer to the nerve lesion but were obtained in a laboratory model and in a different species, in contrast to the human origin and the real clinical setting in which both the arthritis study and the present study were performed.
SIRT2 is involved in the dynamics of the microtubule network in peripheral neurons, which forms the basis for axonal transport of proteins, RNA, vesicles, and organelles between the cell body and the axon tip [66]. It has been proposed that the dynamics of this network are maintained at an optimal level by the controlled action of tubulin-acetylating and -deacetylating enzymes [66]. SIRT2 belongs to the latter [67]. Lower tubulin acetylation is associated with lower microtubule stability [68] and lower recruitment of motor proteins to microtubules [69]. Therefore, high levels of SIRT2 in plasma could be a biomarker for lower microtubule acetylation associated with impaired axonal transport in peripheral neurons, and thus be causally involved in neuropathic pain. However, the enzymatic system that maintains the balance may overshoot, as has been shown in Charcot-Marie-Tooth neuropathy [66].
During the present analyses, SIRT2 was accompanied by a second marker, CD244, which remained among the selected features until the selection step (Table 3). CD244 is a cell surface receptor expressed on natural killer cells that activates cytotoxicity [70]. It has also been involved in cancer [71]; however, any direct involvement in pain has not yet been reported, although this is entirely conceivable via its immune modulation. In the present cohort, CD244 was higher in patients with neuropathic pain, which would be consistent with activated immune and inflammatory responses. The patients with painful knee arthrosis also had significantly higher CD244 levels compared with healthy controls in CSF, but not in serum (99).
Because the present analysis focused on reducing the Proseek multiplex inflammation panel [21] to the most relevant proteins associated with PPSNP after breast cancer surgery, it was important to define whether the selection represents, in functional terms, the entire panel or only proteins with specific molecular functions within the whole panel. To this end, an enrichment analysis was implemented as an overrepresentation analysis (ORA [72]) of the annotations to the genes encoding the selected proteins in the Gene Ontology (GO) knowledge base [56], where the current knowledge about genes is formulated using a controlled vocabulary of GO terms (categories) to which the genes [73] are annotated [74]. GO terms are related by "is-a", "part-of", and "regulates" relationships and form a poly-hierarchy represented as a directed acyclic (DAG [75]). The GO database can be searched by three main categories, namely biological processes, cellular components, and molecular functions. The GO category of molecular function, defined as molecular-level activities performed by gene products, such as "catalysis" or "transport" [56], was used as the functional selection of proteins was the main interest in this assessment. Hence, the 19 proteins identified as informative for the presence or absence of neuropathic pain after nerve injury in breast cancer surgery, were submitted to ORA with the whole Proseek multiplex inflammation panel as reference gene set. The analyses were carried out as described previously [76], using our R library "dbtORA" (https://github.com/IME-TMP-FFM/dbtORA (accessed on 14 March 2022) [77]), which in turn uses the data provided with the R packages "org.Hs.eg.db" (https://bioconductor.org/packages/release/data/annotation/html/org.Hs.eg.db.html (accessed on 14 March 2022) [31]) and "GO.db" (https://bioconductor.org/packages/release/ data/annotation/html/GO.db.html (accessed on 14 March 2022) [78]) with the GO base version of 17 March 2021. For comparison, the full Proseek was analyzed against all human genes, using a p-value threshold of 0.05 and false discovery rate correction [79] for multiple testing performed by means of Fisher's exact tests [35]. There, as a basis for selecting the most appropriate terms to describe the functional genomics roles of the genes of interest, so-called "headline terms" were used that to capture the main content of the poly-hierarchy resulting from ORA [80]. This analysis identified the terms GO:0098772 = molecular function regulator, GO:0005515 = protein binding, GO:0005488 = binding, GO:0060089 = molecular transducer activity, GO:0004175 = endopeptidase activity, GO:0008233 = peptidase activity and GO:0008236 = serine-type peptidase activity as the main molecular functions covered by the Proseek panel. Functionally contrasting the 19 genes coding for the 19 selected proteins with the genes coding the proteins of the whole panel was successful only when leaving out a correction; however, then a shift toward chemokines was observed with headline GO terms GO:0048020 = CCR chemokine receptor binding, GO:0001664 = G protein-coupled receptor binding, GO:0008009 = chemokine activity and GO:0042379 = chemokine receptor binding ( Figure 5).
In addition, the signaling pathways involving the currently analyzed proteins were assessed in a reactome pathway-based analysis using the R library "ReactomePA" (http://bioconductor.org/packages/release/bioc/html/ReactomePA.html (accessed on 14 March 2022) [81]) with its default parameter settings. This again pointed at chemokine signaling as also observed in the results of the above ORA, with the pathways involving the finally selected proteins including chemokine receptors bind chemokines, peptide ligand-binding receptors, interleukin-10 signaling, class A/1 (rhodopsin-like receptors), GPCR ligand binding, and G alpha (i) signaling events ( Figure 6).
Further interpretation of the obtained results addressed the therapeutic potential of the present results, and known drugs were screened for an interaction with the d = 19 proteins of particular interest. This was done using the DrugBank database [82] at https://go.drugbank. com (version 5.1.8 dated 3 January 2021, accessed on 16 December 2021). The database was downloaded as an XML file (https://go.drugbank.com/releases/5-1-8/downloads/allfull-database, accessed on 14 March 2022) and processed using the R package "dbparser" (https://cran.r-project.org/package=dbparser (accessed on 14 March 2022) [83]). Cambinol is an experimental inhibitor of SIRT2 and is being investigated for use in cancer treatment. Any of the 19 proteins were listed as human targets for a total of 41 drugs, of which three were classified in the DrugBank as approved (amiloride, danazol, and chondroitin sulfate) and six were investigational drugs (fibrinolysin, ROX-888, CAT-213, CRx-139, LLL-3348, and again chondroitin sulfate), with the latter classified twice in the DrugBank. According to the DrugBank database, danazol is a steroid used to treat endometriosis and severe pain and tenderness associated with benign fibrocystic breasts, and chondroitin sulfate is used for osteoarthritis, which is also consistent with the overlap currently noted in the proteomics of both types of painful conditions. ROX-888 is being developed for severe acute pain and postoperative pain, CAT-213 is an antiallergic agent, and CRx-139 is being developed for the treatment of immune-inflammatory diseases, while LLL-3348 is intended for the treatment of psoriasis. Thus, the identified proteins point to very plausible drugs that clearly have a link to immunity, and the mention of pain among their possible clinical indications is also noteworthy. In addition, the node's text will be colored in blue to indicate that this node is a detail. Yellow: Significant nodes with highest remarkableness in each path from a detail to the root, i.e., the socalled "headlines". The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R library "dbtORA"  [62]). The color coding is as follows: No color: GO terms that are important for the DAG's structure but do not have a significant p-value in Fisher's exact tests. Red: Significantly overrepresented nodes. Green: Significantly underrepresented nodes. Blue: Terms at the end (detail) of a branch of the DAG. In addition, the node's text will be colored in blue to indicate that this node is a detail. Yellow: Significant nodes with highest remarkableness in each path from a detail to the root, i.e., the so-called "headlines". The figure has been created using the R software package (version 4.0.2 for Linux; https://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R library "dbtORA" (https://github.com/IME-TMP-FFM/ dbtORA (accessed on 14 March 2022) [64]) with the DAG creation done with the GraphViz software package (https://graphviz.org (accessed on 14 March 2022) [77]). The observed patterns in proteomics appeared to be present in both samples, i.e., those taken before surgery and chemotherapy and those at 4 to 9 years follow-up, alt-  [36]) and the R library "ReactomePA" (http://bioconductor. org/packages/release/bioc/html/ReactomePA.html (accessed on 14 March 2022) [81]).
The observed patterns in proteomics appeared to be present in both samples, i.e., those taken before surgery and chemotherapy and those at 4 to 9 years follow-up, although in the second sample the patterns associated with neuropathic pain appeared to be more pronounced. This could indicate protective or risk factors that the patients had already before surgery. It strengthens the association of the observed informative proteins with neuropathic pain and not with changes associated with time, different treatments, or cancer progression, which could have, though not specifically, accompanied the development of postoperative neuropathic pain between the two serum samples. However, the difficulties in observing clear differences between the preoperative and postoperative samples may also be related to the ultimately small sample size of the cohort. However, this is outweighed by the plausibility of the results, their partial replication of findings with persistent pain in independent cohorts, and their reflection in contemporary drug development activities. An independent verification of the present set of proteins most relevant to the development of neuropathic pain after intraoperative nerve injury in breast cancer will probably require a similar study, possibly with a narrower hypothesis that can be based on the present results, increasing the power of the study and possibly also enrolling a larger sample for this purpose. The present results are plausible in light of preclinical research, so return to preclinical models in rodents may not seem warranted. On the other hand, potential drugs resulting from the present findings may also need to be tested in patients, giving preference to experimental pain models in healthy subjects. That is, although systematic analyses have shown that experimental human pain models predict the clinical analgesic effects of drug candidates quite well when the right model for the clinical target is selected from a wide range of human experimental pain models [84][85][86], including models that appear to be predictive even for neuropathic pain drugs such as pregabalin [87], the complexity of the current clinical setting, including nerve injury and cancer treatment, may limit the utility of studies in healthy volunteers. However, depending on the particular characteristics and effects of a future new drug, it is difficult to predict the exact steps of drug development.
The present analyses were performed in serum, consistent with the increasing popularity of blood-derived biomarkers over CSF-derived markers as a more convenient and noninvasive approach for biomarker-based individualized prognosis and treatment of pain [88]. However, with the current analytical methods, the CSF samples are still more sensitive to detect differences in proteomics analyses when assessing pain associated with neuropathy (99). Since the present cohort consisted of women treated for breast cancer with drugs that promote peripheral nerve damage [89], the results need to be confirmed with larger cohorts of patients who do not have cancer.

Patients and Study Design
The Coordinating Ethics Committee of the Helsinki and Uusimaa Hospital District had approved the study, which was also registered at ClinicalTrials.gov (NCT02487524). All patients gave informed written consent. The study cohort consisted of a subset of patients from the NeuroPain study [3], which is a follow-up study of the original BrePainGen cohort in which perioperative pain and related psychological and genetic factors were examined in 1000 women undergoing surgery for breast cancer for unilateral, non-metastatic breast cancer at the Helsinki University Hospital between 2006 and 2010 [20]. Breast surgery consisted of either mastectomy or breast-conserving surgery with sentinel lymph node biopsy or axillary lymph node dissection. None of the patients had received neoadjuvant treatment. Postsurgical treatment consisted of chemotherapy, hormonal therapy and radiotherapy, according to the clinical guidelines.
Details of the clinical conditions and patient characteristics have already been described [3]. The NeuroPain cohort was recruited 4-9 years later in 2014-2016 from the BrePainGen cohort to study factors that associate with the development of neuropathic pain in patients who had a surgeon-verified complete or partial resection of the intercostobrachial nerve (ICBN) during surgery. The main inclusion criterion for the current sub-cohort was a surgeon-verified ICBN injury without persistent postsurgical neuropathic pain (non-PPSNP group) or with definite PPSNP and clinically meaningful pain intensity on a numerical rating scale (NRS, 0-10) ≥4, and no active cancer.

Acquisition of Pain-Related Information
At the preoperative visit, patients rated their pain during the past week in the area to be operated on and elsewhere, separately, on an 11-point numerical scale (NRS) (0 = no pain, 10 = worst pain imaginable). At the follow-up visit, sensory examination was performed to establish a diagnosis for PPSNP according to the latest NP grading criteria [1]. Other pain-related information collected from the patients included rating pain intensity on an 11-point numeric scale (0 = no pain, 10 = worst pain imaginable) by completing the Brief Pain Inventory (BPI) [90] for the worst pain experienced in the surgical area and elsewhere during the past week.

Blood Samples and Quantification of Serum Concentrations of Inflammatory Proteins
At the follow up visit, blood samples were collected for standard laboratory analysis of high-sensitivity C-reactive protein (hs-CRP) and oroso-mucoid (ORM), lipids (total cholesterol, high-density lipoproteins, low-density lipoproteins, and triglycerides), and 25-hydroxyvitamin-D). The results of these assessments have been reported previously [3].
For the proteomics analyses (Olink Analysis Service Uppsala, Uppsala, Sweden), blood samples were collected both before surgery in the BrePainGen study and at the follow up visit in ethylenediaminetetraacetic acid (EDTA) tubes and centrifuged at 3000 min −1 for 10 min. Serum was then transferred to cryotubes and the samples were immediately frozen and stored at −80 • C. The samples were collected and prepared by the same research nurse both preoperatively and at follow-up. The frozen samples were shipped on dry ice to Olink Proteomics, Uppsala, Sweden, for assay. The details of the assay have been described in detail by Wiberg et al. [91]. In brief, 92 proteins from the Proseek multiplex inflammation panel (https://bio-protocol.org/bio101/r9741259 (accessed on 14 March 2022) [21]) were quantified using a proximity extension assay (PEA) that involves two separate antibodies that bind to the same protein in a sample. Each antibody is coupled to a cDNA strand that is ligated on approach, extended by a polymerase, and finally detected using a Biomark HD 96 real-time dynamic PCR array (Fluidigm, South San Francisco, CA). Two incubation controls comprising green fluorescent protein and phycoerythrin were included in the assay to determine the lower limit of detection and to normalize the measurements. A normalized protein expression value (NPX) was calculated for each protein in the sample by normalizing the Ct values by subtracting the values for the extension control and an inter-plate control. The scale was shifted by a correction factor (normal background noise) [91]. Further details about initial laboratory data processing can be obtained at https://www.olink.com (accessed on 16 December 2021).

Summary of the Concept of Data Analysis
The goal of the study was to identify proteins from the Proseek multiplex inflammation panel that are most informative in discriminating patients with and without PPSNP after similar nerve injury during breast cancer surgery. The goal was translated into the task of "feature selection", i.e., reducing data dimensionality by filtering out uninformative or redundant variables to simplify models for easier interpretation by field researchers [92]. Feature selection prior to training computational algorithms is a standard practice for improving classifier performance and reducing the computational burden of training and applying the algorithms. However, in addition to its main application of automatically assigning cases to classes or subgroups, supervised machine learning can also be used to discover structures in the data in order to obtain a description that provides better insights about the dataset. This knowledge discovery approach assumes that, if a classifier can be trained to identify whether a patient belongs to the PPSNP or non-PPSNP subgroup better than by guessing, then the features, i.e., the proteins in the dataset needed by the classifier to accomplish this task, contain relevant information about the addressed patient subgroup structure. In this way, the most informative proteins can be identified. In this use of feature selection, creating a powerful classifier is not the final goal, but feature selection takes precedence over classifier performance. This means that the analysis is considered as successful when the class assignment is just better than guessing and the variables needed for this assignment have been identified, and not necessarily that the classifier is further tuned.
Examples of feature selection methods [92] established in biomedical research [93] include classical approaches such as principal component analysis (PCA [94,95]), regression-based methods such as Least Absolute Shrinkage and Selection Operators (LASSO [96]), and methods based on generally well-performing machine learning methods such as the "Boruta" method [43] or an item categorization-based selection of the most important features for a classifier's performance [97], both of which use the commonly used random forests machine learning classification algorithm as their basis [98,99]. For the present analysis, PCA and the "Boruta" method were used as a representation of a classical statistical approach and an established supervised machine learning approach. To evaluate whether the selected features actually contain information relevant to the subgroup structure in the present patient cohort, the identified features were then used to train a set of classifiers of different types, so as not to rely on the specifics of a single method, but to use a range of methods to internally validate the obtained results. The task here was to achieve better classification than random assignment to PPSNP or non-PPSNP subgroups, and this should not be similarly possible with the other proteins that were not selected as informative for this subgroup structure, nor should it be achieved when the classifiers were trained with permuted proteomics information, i.e., when the internal relationships of the protein levels to the pain-related subgroups were intentionally broken.
In its main components, the data analysis follows the previously proposed workflow for omics data from chronic patients [100] and is shown in a schematic drawing in Figure 7. The necessary programming work was performed in the R language [101] using the R software package [36], version 4.0.2 for Linux, which is available free of charge in the Comprehensive R Archive Network (CRAN) at https://CRAN.R-project.org/ (accessed on 14 March 2022). Analyses were performed on an Intel Core i7-10510U (Intel Corporation, Santa Clara, CA, USA) notebook computer running Ubuntu Linux 20.04.1 LTS 64-bit (Canonical, London, UK)). The detailed descriptions of the data analysis are provided in the following sections.

Quantitative Information Analyzed
Pain-related information consisted of the presence or absence of PPSNP, scaled as [0, 1]. The proteomic panel included initially d = 92 different proteins [91]; however, d = 74 variables could be included in the analyses as the remaining proteins were below the detection level. The proteomic variables consisted of normalized serum protein expression value (NPX) [91], acquired before and 6.6 ± 1.2 (mean ± standard deviation) years after surgery. Thus the proteomic information provided a 74 × 114 (d × 2n) sized data space D = (x i , y i ) x i ∈ R X , y i ∈ Y{1, 2}, i = 1 . . . n , which contained the information, x i on d = 74 proteomic markers acquired at two time points from n = 57 patients, and an output data space, y i , that included the criteria for the grouping into two classes, i.e., the two patient groups comprising "nerve injury and no NP" (non-PPSNP)" and "nerve injury and NP" (PPSNP). The proteomics data set was complete and did not require imputations. Raw data, separated by subgroups and time of sampling, are shown in Figure 8.

Quantitative Information Analyzed
Pain-related information consisted of the presence or absence of PPSNP, scaled as [0,1]. The proteomic panel included initially d = 92 different proteins [91]; however, d = 74 variables could be included in the analyses as the remaining proteins were below the detection level. The proteomic variables consisted of normalized serum protein expression value (NPX) [91], acquired before and 6.6 ± 1.2 (mean ± standard deviation) years after surgery. Thus the proteomic information provided a 74 × 114 (d × 2n) sized data space , | ∈ ℝ , ∈ Y 1,2 , 1 … , which contained the information, xi on d = 74 proteomic markers acquired at two time points from n = 57 patients, and an output data space, yi, that included the criteria for the grouping into two classes, i.e., the two patient

Data Projection-Based Assessment of Proteomics Data Structures Relevant to Pain-Related Subgroup Separation
PCA was performed using the recently proposed "PC-corr" approach [33]. This is an algorithm that facilitates PCA to find a data transformation that optimizes subgroup segregation by retrieving the correlations of the features that produce the segregation of the subgroups along a principal component (PC). It calculates different quality measures for each combination of PC, normalization and centering, and uses different transformations of the data. If its results consist of non-significant separations that are evaluated by quantitative analyses (expressed as p-value, AUC and AUPR) using any type of normalization and dimension, then a nonlinear dimensional reduction is required, since the data is difficult to linearize by different types of normalization. If it turns out that the significant separations, assessed by means of a Mann-Whitney U test [26,27], correspond to certain types of normalization and in dimensions that are not within the first three dimensions of the embedding, then the data has nonlinearities that can be treated by normalizing the data. Therefore, significant group separations in PC1-3 were sought in the PC-corr results as a basis for deciding on the most appropriate data transformation. This analysis was performed using an R script provided with the description of the PC-corr analysis (pccorrv2.R, https: //github.com/biomedical-cybernetics/PC-corr_net (accessed on 14 March 2022) [33]). The results of this analysis indicated that the data set should be probe-level quantile normalized [32] for further analysis. This was performed using the R library "preprocessCore" (https://www.bioconductor.org/packages/release/bioc/html/preprocessCore.html (accessed on 14 March 2022) [102] groups comprising "nerve injury and no NP" (non-PPSNP)" and "nerve injury and NP" (PPSNP). The proteomics data set was complete and did not require imputations. Raw data, separated by subgroups and time of sampling, are shown in Figure 8.  Table 1 for standard protein names. The box plots show the raw values of proteomic marker levels in the plasma of the patients, separately for the first (before surgery) and second (at follow-up 4-9 years later) plasma sample and for the patients with nerve injury but no neuropathic pain ("non-PPSNP") and patients with nerve injury in whom neuropathic pain developed "PPSNP". The boxes were constructed using minimum, quartiles, median (solid line inside the box) and maximum. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. The presentation of the data has been arbitrarily split into two panels to enhance visibility. SIRT2 as a major result of the analysis is highlighted in red; for statistical details, see Table 1 Table 1 for standard protein names. The box plots show the raw values of proteomic marker levels in the plasma of the patients, separately for the first (before surgery) and second (at follow-up 4-9 years later) plasma sample and for the patients with nerve injury but no neuropathic pain ("non-PPSNP") and patients with nerve injury in whom neuropathic pain developed "PPSNP". The boxes were constructed using minimum, quartiles, median (solid line inside the box) and maximum. The whiskers add 1.5 times the interquartile range (IQR) to the 75th percentile or subtract 1.5 times the IQR from the 25th percentile. The presentation of the data has been arbitrarily split into two panels to enhance visibility. SIRT2 as a major result of the analysis is highlighted in red; for statistical details, see Table 1) The figure has been created using the R software package (version 4.0.2 for Linux; http://CRAN.R-project.org/ (accessed on 14 March 2022) [36]) and the R library "ggplot2" (https://cran.r-project.org/package=ggplot2 (accessed on 14 March 2022) [37]).
In the relevant PCs resulting from the PCA described above, subgroup structures consistent with the prior classification (before versus after surgery, PPSNP versus non-PPSNP) were sought by means of Gaussian mixture modeling. Specifically, the distribution of the coordinates of the data set instances (observations) on the principal component space was described by the Pareto density estimation (PDE), which is a kernel estimator of the probability density function (PDF) that has been designed for group discovery [34]. Modal structures were analyzed by fitting Gaussian mixture models (GMM) to the PDE, using our interactive R tool "AdaptGauss" (https://cran.r-project.org/package=AdaptGauss (accessed on 14 March 2022) [103]). The quality of the fit was monitored using the root mean squares, and finally assessed using a Kolmogorov-Smirnov test [104] of the distribution of fitted versus observed data and visual inspection of the quantile-quantile plots of quantiles of the observed data versus the theoretical quantiles according to the fitted model. The assignment of subjects to the identified subgroups was determined using the Bayesian Theorem [105], which provides the decision limits for assigning a single observation to mode M i based on the calculation of posterior probabilities. The correspondence of the group assignment based on the Gaussian modes in the relevant PCs with the a priori subgroup distribution was statistically evaluated using Fisher's exact tests [35].
As an alternative data projection method, self-organizing maps of artificial neurons were used [106] in a modification where the network consisted of a two-dimensional toroid grid with 50 rows and 80 columns [107] that has been shown to be well suited to subgroup detection in biomedical data [38]). Each neuron holds, in addition to a position vector on the two-dimensional grid, a further vector carrying "weights" of the same dimensions as the input dimensions. The weights were initially drawn randomly from the sets of data variables and subsequently adapted to the data during the learning phase with 20 epochs. Following training of the neural network, an ESOM was obtained that represented the subjects on a two-dimensional toroid map as the localizations of their respective "best matching units" (BMU). On the top of the obtained grid of trained neurons, the distances between the data points were calculated using the so-called U-matrix [39,108]. Each value (height) in the U-Matrix represents the average high-dimensional distance of one prototype in relation to all immediately adjacent prototypes in terms of grid position. The corresponding visualization technique uses a topographic map including the coloring, which facilitates the recognition of distance-and density-based structures. Large "heights" in brown and white colors represent large distances between the data. These calculations were performed using the R package "Umatrix" (https://cran.r-project.org/package= Umatrix (accessed on 14 March 2022) [41]).

Supervised Machine-Learning Based Assessment of Proteomics Data Structures Relevant to Pain-Related Subgroup Separation
As an established method of feature selection in machine learning that precedes training of various different types of classifiers in different research environments, the random forest-based Boruta approach [43] was used to identify the most informative protein makers for partitioning the patient cohort into PPSNP and non-PPSNP subgroups. "Boruta" provides a decision on whether a variable is important or not for the classification task, which is derived from a 100-fold cross-validation approach followed by statistical evaluation of the variables importance with p-values defaulting to 0.01 [43]. These calculations were performed with the R package "Boruta" (https://cran.r-project.org/package=Boruta (accessed on 14 March 2022) [43]) with the default hyperparameter settings.
To further enhance the validity of the feature selection, the Boruta approach was nested into a 1000 cross-validation scenario using each time 2/3 of the data set randomly drawn class-proportionally from the original data set by means of using Monte Carlo resampling [109] implemented in the R library "sampling" (https://cran.r-project.org/ package=sampling (accessed on 14 March 2022) [42]). The features selected by the Boruta algorithm during each run were collected, and the final set of proteins was assembled in descending order of the frequency with which they were among the selected features in the 1000 cross-validation Boruta runs. The cutoff value for the selection was set using the computed ABC analysis [110]. This item categorization method divides each set of positive numbers into three non-overlapping subsets "A", "B", and "C" [111], of which category "A" contains the "important few" that have been retained in the present analyses. The exact computations of the set limits "A/B" and "B/C" have been described elsewhere [110]; the calculations were performed using our R package "ABCanalysis" (https://cran.r-project. org/package=ABCanalysis (accessed on 14 March 2022) [110]).

Supervised Machine Learning-Based Evaluation of Identified Proteomic Markers to Distinguish Pain-Related Patient Subgroups
The final step of the data analysis consisted in an evaluation of the identified proteomic markers to provide, in a variety of classification algorithms, suitable information about the segregation of the patient cohort into PPSNP or non-PPSNP subgroups. Classifier training and testing was performed in a 100-fold cross-validation design using disjoint training (2/3 of the data) and test (1/3 of the data) data subsets obtained by means of Monte-Carlo random resampling. Classification performance was evaluated primarily on the basis of balanced accuracy [112]. Further performance criteria included the area under the receiver operator curve (AUC-ROC [113]), sensitivity, specificity, precision, recall, positive and negative predictive value [114,115] and the F1 measure [116,117]. These calculations were performed with the R libraries "caret" (https://cran.r-project.org/package=caret (accessed on 14 March 2022) [118]) and "pROC" (https://cran.r-project.org/package= pROC (accessed on 14 March 2022) [119]).
The classifiers were trained with the selected proteomic markers, as these were of most interest in this evaluation of the results obtained in the previous steps of data analysis. If these markers enabled the algorithms to assign patients to pain subgroups better than by guessing, the selected proteins could be considered informative for this clinical subgrouping. To control for possible overfitting, all machine learning algorithms were additionally trained with randomly permuted proteomic markers, with the expectation that a classifier trained with these data should not perform better than guessing, i.e., give a balanced accuracy or an AUC-ROC around 50 %. Furthermore, classifiers were trained with all protein markers, and again with the protein markers that were not selected during feature selection, in order to ensure that the selection had indeed identified the most informative markers.

Conclusions
Present analyses pointed in particular to sirtuin 2, with its role in neuroinflammatory processes and in learning and memory, as a key marker in the development of PPSNP. Results extended to 18 other proteins that were informative in distinguishing between samples from patients with neuropathic pain and those without neuropathic pain, without a clear distinction between samples before or after surgery. This suggests that the proteomic patterns were not simply a consequence of the development of neuropathic pain or other influences after surgery but reflected risk or protective factors that were already present before surgery. The identified informative proteins had a remarkable number of target proteins for approved or investigational drugs that have pain, including postoperative pain or chest pain, as a clinical target, providing remarkable support for the relevance of the present results.