An Erythroid Gene Expression Signature in Preoperative Whole-Blood RNA-Seq Associated with Rebound Pain After Brachial Plexus Block: A Prospective Observational Study

Background

Peripheral nerve blocks are widely used in upper extremity surgery to provide reliable regional analgesia and reduce early postoperative opioid exposure.^1–3 Although these techniques confer well-documented perioperative analgesic benefits, an increasing body of evidence has described rebound pain, defined as a marked increase in pain as the numbness from the block resolves.^4,5 Rebound pain has been associated with sleep disruption, delayed functional recovery, and increased use of rescue analgesics.^6,7

To date, research has focused on perioperative and procedural factors, including local anesthetic type, adjuvant use, block duration, psychological factors, and surgical characteristics.^8–10 However, these factors may explain only part of the substantial interindividual variability in rebound pain severity. Notably, postoperative pain trajectories often diverge after block resolution, even among patients with similar demographics, comparable surgical procedures, and similar regional techniques.¹¹

The mechanisms underlying rebound pain after peripheral nerve block remain incompletely understood.¹ Current evidence suggests that rebound pain may involve spontaneous hyperexcitability of peripheral nerve fibers, inflammation secondary to tissue injury, neurotoxic effects of local anesthetics, and local nerve damage.^12–16 Small fiber neuropathy, which affects thinly myelinated Aδ fibers and unmyelinated C fibers and is difficult to diagnose and treat, may also influence postoperative pain trajectories and rebound pain assessment.¹⁷

Although neural and local tissue mechanisms are directly relevant to rebound pain, sampling these tissues before surgery is generally impractical and may expose patients to unnecessary risk. Whole-blood transcriptomic profiling may instead provide a feasible way to assess systemic gene expression patterns relevant to pain biology.¹⁸ Erythroid differentiation is closely linked to heme synthesis and hemoglobin production, while heme metabolism has been implicated in nociceptive biology.^19,20 For preoperative systemic profiling, whole blood preserves transcriptomic information from multiple circulating cell populations, whereas isolated cell populations provide cell-type-specific profiles but require additional processing that may introduce processing-related variability.²¹ However, bulk whole-blood RNA-seq may be influenced by blood cell composition, so erythroid-related signals should be interpreted with this potential confounding in mind. Although transcriptomic studies have been conducted in several acute and chronic pain models,^22,23 gene expression patterns specific to rebound pain have not been characterized.

In this study, we compared preoperative blood gene expression profiles between patients who developed rebound pain after nerve block and those who did not to delineate transcriptomic features associated with differential pain responses. Protein-protein interaction analysis and machine learning analyses were used to identify coordinated gene clusters and candidate transcripts, and key findings were validated in an independent cohort using qRT-PCR. With further validation, preoperative biomarkers associated with rebound pain may help identify patients at increased risk before surgery and support individualized analgesic planning, perioperative counseling, and closer postoperative follow-up after block resolution.

Methods

This was a single-center, prospective observational study. The study was carried out in accordance with the Declaration of Helsinki, and the research protocol was approved by the Ethics Committee of Xiangyang Central Hospital, Affiliated Hospital of Hubei University of Arts and Science, Xiangyang, China. All participants provided written informed consent. Reporting followed the Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) statement.²⁴

Study Setting and Population

Patients aged 18–75 yr with American Society of Anesthesiologists (ASA) physical status I–III who were scheduled to undergo open reduction and internal fixation (ORIF) of upper-extremity fractures under brachial plexus block at Xiangyang Central Hospital (Xiangyang, China) were eligible.

The exclusion criteria were as follows: concurrent surgery at other anatomical sites; emergency surgery; peripheral neuropathy or pre-existing nerve injury of the operative upper limb; prior axillary surgery on the operative side; more than one previous surgery on the operative upper limb; allergy to local anesthetics or contraindications to postoperative analgesics (eg. nonsteroidal anti-inflammatory drugs, ibuprofen); pregnancy; diabetes mellitus or vascular disease; chronic use of opioids; antipsychotics; or antidepressants; infection at the injection site; acute porphyria; or renal disease requiring dialysis; contraindications to regional anesthesia; requirement for supplemental general anesthesia because of an inadequate block or inadequate anesthesia; or conversion to general anesthesia; chronic pain syndrome; inability to complete perioperative questionnaires because of language barriers or cognitive impairment; or participation in similar studies within the preceding 3 months. Participants were recruited from April 2024 until the required sample size was achieved in September 2024.

Clinical Management

Midazolam (2 mg) and dexamethasone (8 mg) were given intravenously preoperatively. Brachial plexus block was performed by trained anesthesiologists under real-time ultrasound guidance, with the approach selected according to the patient’s surgical site: 20 mL of 0.4% ropivacaine for interscalene block, 15 mL of 0.4% ropivacaine for axillary block, or a combination of both approaches. To assess sensory block adequacy prior to skin incision, cold testing (using a 0°C glass ampoule) was conducted across multiple anatomical territories. Patients with incomplete blocks were converted to general anesthesia and excluded from the analysis.

Throughout the surgery, sedation was provided with dexmedetomidine: a loading dose of 1 μg/kg infused over 10 min, followed by a continuous infusion of 0.5 μg·kg⁻¹·h⁻¹ until the end of surgery. For postoperative pain management, patients were transferred to the ward and received a standardized regimen consisting of oral ibuprofen 400 mg every 8 h and patient-controlled analgesia (PCA) with sufentanil (background infusion 5 μg/h, bolus 3 μg, lockout 15 min).

The primary outcome was the incidence of postoperative rebound pain within 48 h after brachial plexus block. Pain intensity was assessed using an 11-point numeric rating scale (NRS), where 0 indicated no pain and 10 indicated the worst imaginable pain. NRS scores were recorded preoperatively and at 08:00 and 20:00 on postoperative days 1, 2, and 3 by trained research staff. Rebound pain was defined as a transition from well-controlled pain, defined as NRS ≤ 3, to severe pain, defined as NRS ≥ 7, within the first 48 h after block performance. This definition was adapted from the criteria proposed by Tim et al.²⁵ Patients who met this criterion were assigned to the rebound pain (RP) group, whereas those who did not were assigned to the no rebound pain (NRP) group.

Transcriptome Analysis

RNA Extraction

RNA extraction, library preparation, and sequencing were performed by Wuhan Metware Biotechnology Co. Ltd. (Wuhan, China) in accordance with the manufacturer’s protocols (Illumina, San Diego, CA, USA). Before induction of anesthesia, 2 mL of peripheral venous blood was collected from each participant into EDTA tubes. TRIzol^® (6 mL) was added, the sample was mixed gently by inversion, incubated for 5 min, and then stored at −80°C until processing. RNA concentration was measured using a fluorometer (Thermo Fisher Scientific, Waltham, MA, USA), and RNA integrity was assessed with a Qsep400 Bioanalyzer (Bioptic, Changzhou, China). Libraries were constructed only from samples that met quality requirements.

mRNA Library Construction

Poly(A)+ mRNA was enriched using oligo(dT) magnetic beads. Purified mRNA was fragmented under controlled conditions using fragmentation buffer. First-strand cDNA was synthesized by reverse transcription with random hexamer primers. For strand-specific library preparation, dUTP was substituted for dTTP during second-strand synthesis to prevent amplification of the uracil-containing strand by the high-fidelity DNA polymerase, thereby preserving strand specificity. End repair and dA-tailing were performed. Adapter ligation was performed. Libraries were purified and size selected using magnetic beads to obtain insert sizes of 250–350 bp, then PCR-amplified and purified again. Final libraries were eluted in nuclease-free water. Library concentration was quantified by fluorescence, and fragment size distribution was evaluated using the Qsep400 system. Qualified libraries were pooled based on effective concentration and target sequencing depth and sequenced on an Illumina^® platform to generate 150-bp paired-end reads.

RNA-Seq Quality Control, Alignment, and Quantification

Raw sequencing reads were subjected to quality control using fastp v0.23.2. Adapter-containing reads were removed. Paired reads were discarded if the proportion of N bases in either read exceeded 10% of the read length or if bases with Q ≤ 20 accounted for more than 50% of either read. Clean reads were used for subsequent analyses. Sequencing quality was evaluated using raw reads, clean reads, clean bases, overall sequencing error rate, Q20, Q30, and GC content for each sample. Per-base sequencing error-rate and GC-content distributions were also inspected as part of quality assessment. Clean reads were aligned to the human reference genome from the Ensembl database, annotation version Homo sapiensGRCh38.109, using HISAT2 v2.2.1 with default parameters. Gene-level raw read counts were generated using featureCounts v2.0.3 with default parameters. FPKM and TPM values were calculated for descriptive expression-level visualization, whereas raw count matrices were used as input for differential expression analysis. All RNA-seq samples were processed and sequenced together in a single experimental batch by the same sequencing service provider using the same platform and workflow. Therefore, no predefined technical batch covariate, such as library-preparation batch or sequencing run, was available for batch-effect adjustment, and batch-effect correction was not performed. PCA, multidimensional scaling, and sample-distance heatmap analyses were performed as unsupervised sample-level quality assessment and are shown in Supplementary Figure 1–3.

Gene Set Enrichment Analysis

Gene set enrichment analysis (GSEA) was conducted to identify biological processes associated with susceptibility to rebound pain using a ranked list of all expressed genes.²⁶ Genes were ranked according to the differential expression test statistic derived from the RP versus NRP comparison. Enrichment was assessed against Gene Ontology gene sets using the fgsea package. Enrichment significance was evaluated using the normalized enrichment score (NES), and multiple testing across gene sets was controlled using the false discovery rate (FDR). Gene sets with FDR-adjusted P values < 0.05 were considered significantly enriched. Among the significantly enriched gene sets, positively and negatively enriched terms were ranked separately according to the absolute value of NES in descending order. The top 10 positively enriched and top 10 negatively enriched terms were first selected. Within each direction, the five terms with the lowest FDR-adjusted P values were further prioritized as core enriched terms for downstream interpretation and visualization.

Differential RNA Expression Analysis and GO Enrichment Analysis

Differential gene expression (DGE) analysis was performed using DESeq2 v.1.38.3²⁷ by comparing the RP group with the NRP group in the age-, and sex-matched RNA-seq cohort. Raw count matrices were normalized using the median-of-ratios method implemented in DESeq2 to account for differences in sequencing depth and RNA composition across samples. P values were adjusted for multiple testing across all tested genes using the Benjamini-Hochberg procedure to control the false discovery rate.²⁸ Fold changes were moderated using adaptive shrinkage. Genes were considered differentially expressed if they met both criteria: |log2FC| ≥ 0.585, corresponding to fold change ≥ 1.5, and FDR-adjusted P value ≤ 0.05.

Protein–Protein Interaction Network and Gene Module Identification

Protein‒protein interaction (PPI) analysis of the DEGs was conducted through the STRING database (https://string-db.org).²⁹ We uploaded the list of DEGs to the STRING database with the following settings: network type: full STRING network, meaning of network edges: confidence, minimum required interaction score: highest confidence (0.4). We subsequently imported the results into Cytoscape (version: v3.10.0) for visualization and subsequent analysis.³⁰ To identify core hub genes, we leveraged the CytoHubba plugin in Cytoscape and evaluated gene connectivity using six canonical scoring approaches: MNC, MCC, EPC, Degree, Stress, and Radiality (full names: Maximum Neighborhood Component, Maximal Clique Centrality, Edge Percolated Component, Degree centrality, Stress centrality, and Radiality centrality).³¹ Hub genes that appeared in the top 20 rankings of all six methods were ultimately selected as the prioritized hub genes.

Selection of Hub Genes via Machine Learning

Prioritized hub genes initially screened by the Cytohubba plugin were further refined using two machine learning approaches: least absolute shrinkage and selection operator (LASSO) regression and random forest. LASSO regression was implemented using the glmnet R package, and 5-fold cross-validation was performed to determine the optimal regularization parameter (λ). The data was split into five subsets for training and validation, ensuring the model was tested on different portions of the data, which helps reduce overfitting. The optimal λ was selected using λ.1se (one standard error above the minimum λ), which results in a simpler model by limiting the number of selected features, thus reducing the risk of overfitting.³² For Random Forest analysis, the randomForest R package was used, training the model with 800 trees and evaluating gene importance using the MeanDecreaseGini index. Only genes with a MeanDecreaseGini score greater than 1 were retained, ensuring that only the most important genes were included in the final model. Additionally, 10-fold cross-validation was applied to assess the model’s stability and generalizability, further minimizing the possibility of overfitting. By combining results from both methods, three genes were consistently selected.

Predictive Performance Analysis of GYPB

Receiver operating characteristic (ROC) analysis was performed to evaluate the ability of GYPB expression to discriminate between patients with and without rebound pain in the RNA-seq cohort (n = 24; 12 with rebound pain, 12 without rebound pain). ROC curves were constructed using GYPB expression values, and the area under the ROC curve (AUC) was calculated as a measure of overall discriminative performance. To account for the small sample size, 95% confidence intervals for the AUC were estimated using a stratified bootstrap procedure with 5000 resamples, in which samples were resampled with replacement within each outcome group; the 2.5th and 97.5th percentiles of the bootstrap AUC distribution were used as the 95% confidence interval. Statistical significance of the observed AUC was further assessed using a permutation test: group labels were randomly permuted 5000 times, and the AUC was recalculated for each permutation to generate a null distribution. An empirical P value was calculated as the proportion of permuted AUCs greater than or equal to the observed AUC.

Validation of GYPB Expression by qRT-PCR

To independently validate the differential expression of GYPB identified in the discovery cohort, an additional independent validation cohort was recruited for qRT-PCR analysis. This cohort was registered at the Chinese Clinical Trial Registry under the registration number ChiCTR2500099214. The validation cohort included 45 patients, comprising 27 patients in the RP group and 18 patients in the NRP group. Patients were recruited using the same inclusion and exclusion criteria as those applied to the discovery cohort, and peripheral venous blood samples were collected specifically for qRT-PCR confirmation. Total RNA was extracted from peripheral venous blood using TRIzol^® reagent (Takara) according to the manufacturer’s protocol. Complementary DNA (cDNA) was synthesized using the PrimeScript^TM RT reagent Kit with gDNA Eraser (Takara). qRT-PCR was performed using SYBR Premix Ex Taq^TM II (Takara) on ABI 7500 system. The primers used for qRT‒PCR were as follows:

GYPB Forward: 5′-TGGCAATGCACACTTCAACC-3′,

GYPB Reverse: 5′-GGAGCTGGTACAGTGAAACGA-3′.

Relative mRNA expression was calculated using the 2^−ΔΔCt method, with β-actin as the internal reference gene.

Statistical Analysis

The sample size for the transcriptomic analysis was estimated a priori using a negative binomial model. Assuming a coefficient of variation of 0.3, a minimum fold change of 1.5, and a two-sided significance level of 0.05, 12 patients per group would provide an estimated statistical power of 0.88. Given the exploratory nature of RNA-seq analysis and the limitations of model-based power estimation, these results were interpreted as hypothesis-generating. On the basis of previous studies reporting an incidence of rebound pain after brachial plexus block of 45.2%–56.3%,⁸ the lower estimate of 45.2% was used to determine the total enrollment requirement. A binomial probability model showed that 37 evaluable patients would provide a >95% probability of identifying at least 12 patients with rebound pain. After allowing for a 10% dropout rate, the planned enrollment sample size was 41 patients. In total, 47 patients were ultimately enrolled.

Statistical analyses were performed with R software (version: 4.5.1) and IBM SPSS Statistics (version: 26.0). Prior to analysis, the Shapiro–Wilk test was applied to verify the normality of continuous variables. Continuous variables with a normal distribution were summarized as the mean (± standard deviation, SD) and were compared between the two groups via independent-samples t tests; those with a non-normal distribution were reported as the median (25th–75th percentile, IQR) and were analyzed via the Wilcoxon rank-sum test (Mann–Whitney U-test). Statistical significance was defined as a two-tailed P value < 0.05 for all analyses. Given the limited sample size, multivariable analysis was not performed because the number of events per variable would be limited. Potential confounding was therefore addressed by baseline comparisons and propensity score matching for the transcriptomic analysis.

Results

Characteristics of the Participants

From April 2024 to September 2024, 51 consecutive patients met the inclusion criteria. Among these, 3 patients had their anesthesia method changed, and 1 patient refused postoperative blood sampling. These 4 patients were not included in the final analysis, resulting in 47 patients being analyzed. The 47 patients were stratified into the RP group and NRP group according to the definition of rebound pain. For transcriptome analysis, a matched subset was generated using 1:1 nearest-neighbor propensity score matching based on age and sex. This matching procedure yielded 12 matched pairs, including 12 patients in the RP group and 12 patients in the NRP group. No additional patients were excluded after matching.

Baseline characteristics were compared between the RP group (n=28) and the NRP group (n=19) (Table 1). No significant differences were observed in demographic variables (age, BMI), hematological parameters (WBC, NEUT, and LYM counts), procedural variables (tourniquet application time, ropivacaine dose), sex distribution, comorbidities (cerebrovascular disease, hypertension, diabetes mellitus), education level, ASA classification (II/III), nerve block site, or surgical type (P >0.05). However, two procedural indicators differed significantly between the two cohorts: surgical time (median [interquartile range]: 85.00 [60.00, 98.25] min in the RP group vs. 70.00 [49.00, 74.50] min in the NRP group; P<0.05) and PCA sufentanil dose (mean ± SD, 275.14±12.36 μg in the RP group vs. 207.42±11.46 μg in the NRP group; P<0.05).

Table 1 Demographic and Clinical Data

To evaluate whether whole-blood transcriptomic differences might be influenced by blood-cell composition, we further compared CIBERSORT-estimated leukocyte fractions and routinely measured preoperative hematological indices in the RNA-seq cohort. The top ten estimated leukocyte subsets showed no systematic group-wide separation between RP and NRP samples. Similarly, no significant between-group differences were observed in RBC, HGB, HCT, MCV, MCH, MCHC, RDW-SD, or RDW-CV. These exploratory analyses are shown in Supplementary Figures 4–5.

GSEA Enrichment Results Against the GO Database

GSEA of GO terms revealed distinct transcriptional programs associated with susceptibility to rebound pain. The gene sets related to ribosome biogenesis, rRNA maturation, ncRNA processing, tRNA processing, and catalytic activity acting on tRNA were significantly enriched at the positive end of the ranked gene list (Figure 1A). These enrichment patterns suggest enhanced translational capacity and RNA-processing activity in patients with rebound pain.

Figure 1 GSEA of GO pathways between rebound pain and non-rebound pain patients. (A) GSEA plot of top 5 most significantly upregulated GO pathways with smallest adjusted P values in RP patients, enriched in positive-value gene regions. (B) GSEA plot of top 5 most significantly downregulated GO pathways with smallest adjusted P values in RP patients, enriched in negative-value gene regions.

In contrast, gene sets associated with erythrocyte development, gas transport, hydrogen peroxide catabolic process, innate immune response in mucosa, and organ- or tissue-specific immune response were preferentially enriched at the negative end of the ranked gene list (Figure 1B). These enrichment patterns may reflect peripheral blood transcriptomic alterations, although their mechanistic relevance to rebound pain remains to be further established. Detailed GSEA data for GO terms are provided in Supplementary Table 1.

Analysis of Differential Genes and Enriched Pathways

Among the 27,551 genes detected in both RP and NRP groups, 110 were identified as differentially expressed genes (DEGs), including 108 upregulated and 2 downregulated DEGs (Figure 2A). The distribution of upregulated and downregulated DEGs is visualized in a volcano plot, and the 10 DEGs with the largest fold changes are annotated (Figure 2B). Heatmap analysis of the DEGs revealed that samples from the RP and NRP groups clustered distinctly: the RP samples were primarily distributed on the left side with relatively high expression, whereas the NRP samples were mainly concentrated on the right side with relatively low expression (Figure 2C). Detailed DEG data was provided in Supplementary Table 2. The raw sequencing data have been deposited in the NCBI SRA (BioProject: PRJNA1392155).

Figure 2 Identification of DEGs between rebound pain patients and non-rebound pain patients. (A) Number of upregulated and downregulated genes. (B) Volcano plot of DEGs. The upregulated genes are marked in red; the downregulated genes are marked in blue. (C) Heatmap showing the DEGs between rebound pain patients and non-rebound pain patients.

GO analysis revealed that differentially expressed genes were predominantly enriched in erythrocyte related biological processes, including erythrocyte development, erythrocyte differentiation, erythrocyte homeostasis, and hemoglobin metabolic process. In addition, pathways involved in myeloid cell development and differentiation, gas transport, and iron ion homeostasis were also significantly overrepresented. Genes such as ALAS2, AHSP, EPB42, SLC4A1, RHAG, TAL1, KLF1, and BPGM contributed to these enriched terms, indicating the presence of erythrocyte- and hemoglobin-related transcriptional signals in the differential expression profile (Figure 3). Detailed GO enrichment data was provided in Supplementary Table 3.

Figure 3 GO enrichment analysis of differentially expressed genes between patients with rebound pain and non-rebound pain. Top 10 GO biological processes enriched among differentially expressed genes between rebound pain and non-rebound pain patients, together with representative shared genes contributing to each pathway, ranked by adjusted P value.

PPI Network Analysis of DEGs

DEGs were entered into the STRING database for PPI network construction and 65 nodes and 500 edges were ultimately identified from the resulting network (Figure 4A). In addition, we identified the hub cluster with the genes shared by the top 20 hub genes ranked by each of the six different CytHubba algorithms, which consisted of 13 genes (Figure 4B). The hub gene cluster containing 13 genes (SLC4A1, GYPB, EPB42, AHSP, ALAS2, ANK1, BPGM, CA1, HBM, KLF1, OSBP2, RHAG, SNCA) exhibited prominent clustering characteristics in the PPI network (Figure 4C).

Figure 4 Protein–protein interactions (PPIs) between DEGs and selection of hub genes. (A) PPI network (110 DEGs filtered into the PPI network that contained 162 nodes and 756 edges). (B) Distribution of common genes among the top 20 hub genes identified by different CytHubba algorithms. (C) Hub gene cluster containing 13 genes identified by different CytHubba algorithms.

Identification of the Hub Gene GYPB

Based on the 13 candidate hub genes screened from the hub cluster using CytoHubba, we applied LASSO regression and random forest analyses to identify genes associated with the occurrence of rebound pain. Following the 5-fold cross-validation procedure, LASSO regression identified three hub genes associated with rebound pain in the model (Figure 5A). Random forest analysis identified three hub genes associated with rebound pain (Figure 5B). Intersecting the results from the two methods yielded one shared hub gene (Figure 5C). Accordingly, GYPB was selected as a candidate gene for subsequent predictive analysis and experimental validation.

Figure 5 Identification of hub genes via two machine learning strategies. (A) Coefficient profiles corresponding to the log(λ) sequence derived from the LASSO regression model. The optimal λ value was determined using 5-fold cross-validation and the λ.1se criterion. (B) The random forest algorithm was utilized to screen 3 hub genes, all of which presented MeanDecreaseGini values greater than 1. (C) Venn diagram depicting the overlapping gene set identified by the two machine learning approaches.

Validation of the Predictive Performance of GYPB for Rebound Pain

ROC analysis was performed to assess the predictive performance of GYPB expression for distinguishing RP from NRP patients in 24 samples, including 12 patients in each group. The AUC was 0.792, with a stratified bootstrap 95% CI of 0.585–0.999, suggesting preliminary discriminative ability with limited precision due to the small sample size (Figure 6A).

Figure 6 ROC analysis and permutation test evaluating the discriminative performance of GYPB for rebound pain. (A) ROC curve illustrating the ability of GYPB expression to distinguish RP from NRP patients (AUC = 0.792, 95% CI: 0.585–0.999). (B) Distribution of AUC values derived from 5000 permutation tests; the red dashed line denotes the observed AUC based on GYPB expression, with a permutation P value of 0.013.

To further evaluate whether the observed discriminative performance could be attributed to random chance, a permutation test was performed using 5000 random label permutations. The observed AUC of 0.792 was significantly higher than the majority of permuted AUCs (P = 0.013) (Figure 6B). Nevertheless, given the small sample size, this result should be regarded as exploratory and requires validation in a larger independent cohort.

Validation of GYPB Through qRT-PCR

RNA-seq differential expression analysis revealed a statistically significant upregulation of GYPB in the RP group (P < 0.05) (Figure 7A). To further validate this expression signature, qRT-PCR was performed using preoperative peripheral blood samples from an independent validation cohort of 45 patients, including 27 patients in the RP group and 18 patients in the NRP group. The detailed demographic and clinical characteristics of this validation cohort are provided in Supplementary Table 4. To maintain consistency with the transcriptome analysis, the same 1:1 nearest-neighbor propensity score matching method based on age and sex was applied to the validation cohort using a maximum matching strategy, yielding 15 matched pairs for qRT-PCR analysis. In this matched validation subset, qRT-PCR analysis showed that GYPB mRNA expression was significantly higher in RP patients than in NRP patients, corresponding to a 3.390-fold increase (RP: n = 15; NRP: n = 15; 95% CI: 1.258–3.523; P < 0.05) (Figure 7B).

Figure 7 Expression profiles of GYPB in RNA sequencing and qPCR assays. (A) The expression level of GYPB in the RP group was significantly greater than that in the NRP group. (B) The qPCR assays further confirmed that the relative mRNA abundance of GYPB in the RP group was markedly greater than that in the NRP group. Statistical significance is indicated as follows: *P < 0.05; ***P < 0.001.

Discussion

This prospective transcriptomic study identified a distinctive preoperative gene signature in patients who developed rebound pain after brachial plexus block. Among the differentially expressed transcripts, GYPB showed the most consistent association, with higher expression confirmed in an independent cohort. These findings raise the possibility that preoperative systemic hematologic-related transcriptomic patterns, particularly erythrocyte homeostasis and erythrocyte-mediated gas transport signatures, may be associated with interindividual differences in pain escalation after sensory block resolution.

The preoperative alterations in the erythroid signature observed in this study may help generate exploratory biomarker hypotheses for rebound pain susceptibility. Among these signals, GYPB emerged as a key representative gene of interest. GYPB, which encodes glycophorin B, stabilizes the erythrocyte membrane and mediates interactions between red blood cells, vascular endothelium and immune components.^33–35 These findings raise the hypothesis that erythroid-related systemic physiology may be linked to rebound pain susceptibility, but they do not establish a direct mechanistic relationship between erythrocyte biology, tissue oxygenation, and pain sensitivity.

Recent research focusing on surgical and neuropathic pain has focused on the role of systemic physiology in modulating pain vulnerability.^36–38 When systemic stressors such as operative anemia occur, they may impair perioperative rehabilitation and negatively influence clinical outcomes.^39,40 Although direct evidence linking erythrocyte biology to rebound pain after regional anesthesia is limited, studies in other pain conditions provide indirect biological plausibility for a relationship between oxygenation-related physiology and pain phenotypes. For example, in complex regional pain syndrome, near-infrared spectroscopy studies have shown reduced deep tissue oxygen saturation in the affected limb, findings that have been correlated with pain intensity and functional impairment.⁴¹ In parallel, leukocyte transcriptome profiling after trauma has highlighted heme-handling programs that track systemic inflammatory severity and septic complications,⁴² and blood transcriptomics has been used to explore susceptibility to persistent pain phenotypes in other clinical contexts.^22,43 In our study, GO and GSEA analyses suggested erythrocyte function, gas transport, ribosomal, and RNA-processing-related transcriptomic differences between RP and NRP groups. Collectively, these findings suggest that erythroid-related transcriptomic variation may be associated with rebound pain susceptibility, although the underlying mechanistic relationship remains to be clarified. Although blood-cell composition, systemic inflammation, and perioperative stress-related states remain possible alternative explanations, exploratory CIBERSORT and routine hematological analyses did not reveal systematic group-wide differences between RP and NRP samples; however, residual confounding cannot be fully excluded.

These preoperative transcriptomic signatures may represent preliminary biomarker candidates for rebound pain susceptibility following brachial plexus block. These findings suggest that individual hematologic-related characteristics before surgery could be associated with rebound pain susceptibility after neural blockade resolves. From a clinical perspective, they may encourage future studies evaluating whether preoperative hemoglobin trends, red-cell indices, oxygenation status, and transcriptomic information could improve risk stratification in patients scheduled for regional anesthesia. This perspective is consistent with existing patient blood management and enhanced recovery frameworks,⁴⁴ and may conceptually complement them by highlighting physiological variation not fully captured by routine laboratory tests.⁴⁵ However, clinical application of these transcriptomic markers remains premature and will require validation in larger independent cohorts.

This study has several limitations that should be acknowledged. Its modest sample size and single-center design limit generalizability. Additionally, the observational nature of the data means that causal relationships cannot be established. Although exploratory CIBERSORT and routine hematological analyses did not reveal systematic group-wide differences between RP and NRP samples, residual confounding from unmeasured blood-cell composition, systemic inflammation, or perioperative stress-related physiology cannot be fully excluded. In addition, the modest RNA-seq sample size limited the feasibility of comprehensive adjusted multivariable modeling. The group differences in operative time and sufentanil consumption could introduce residual confounding, although we attempted to mitigate these influences in our analyses. Although ROC curve analysis combined with permutation testing was employed to evaluate the discriminative performance of the model in this study, the relatively small sample size may limit the stability of the AUC estimates and the statistical power of the permutation test. Consequently, these findings should be interpreted as preliminary and hypothesis-generating, and further validation in larger, independent cohorts with appropriate adjustment for clinical and biological confounders is warranted.

Conclusion

These findings suggest that preoperative GYPB upregulation and erythroid-related whole-blood transcriptomic signatures are associated with rebound pain following brachial plexus block. This association may reflect baseline hematologic-related physiological variation, but should be considered preliminary and requires validation in larger independent cohorts with appropriate adjustment for clinical and biological confounders.