Memory B cells and their transcriptomic profiles associated with belimumab resistance in systemic lupus erythematosus in the maintenance phase

The factors contributing to the treatment efficacy of belimumab in patients with systemic lupus erythematosus (SLE) in the maintenance phase are unknown. Here, we collected blood samples from patients with SLE (n=44) treated with belimumab before and three and six months after treatment. RNA-Seq of whole blood was performed, and gene expression was quantified. Immune cell type enrichment analysis estimated immune cell subtype proportions and gene expression in each subtype. The Systemic Lupus Erythematosus Disease Activity Index 2000 (SLEDAI-2K) < 4 at six months was set as the primary efficacy criterion. Non-responders exhibited upregulated B cell proliferation signals before treatment, associated with an increased number of memory B cells. A higher proportion of memory B cells before treatment predicted poor response (p=5.1×10^-4). This was also associated with changes in disease activity and glucocorticoid dose at six months compared with baseline. Belimumab did not affect memory B cell proportion during the treatment time course, in contrast to naïve B cells. Higher memory B cell proportion was associated with higher type-I interferon (IFN) scores and lower white blood cell and complement C4 levels. Transcriptomic analysis of memory B cells in non-responders revealed significant upregulation of immunoglobulin genes (Ig). Memory B cells and high Ig expression in them were identified as a treatment-resistant factor of belimumab in SLE patients. Lower C4 and white blood cell counts may serve as clinical markers of higher memory B cells.

1 Introduction

Systemic lupus erythematosus (SLE) is a chronic autoimmune disorder that affects a variety of organ systems and markedly impairs health-related quality of life (1). Patients with SLE often have increased blood levels of the B lymphocyte stimulator (BLyS; also termed BAFF) (2). It is a member of the tumor necrosis factor ligand family that induces B cell proliferation and immunoglobulin secretion (3). Belimumab is a recombinant, human immunoglobulin G1 lambda monoclonal antibody that binds soluble BAFF protein and inhibits its biological activity (4, 5). In terms of cell type specificity, it has a different mode of action from anifrolumab, which blocks the type I interferon receptor signaling across multiple immune cell types (6).

This drug is currently recommended to be used to achieve disease remission as well as to decrease the dose of oral glucocorticoid (GC) in the maintenance phase (7). However, as with other biological agents, some patients are known to exhibit poor responses to this drug. The recent application of new treatment options, such as anifrolumab, underscores the necessity of developing predictive models for treatment response to belimumab. To date, lower disease activity (5), higher complement (8), lower ds-DNA antibody (8), established organ damage (9), lower BlyS (10, 11), a higher number of B cells (12), a smaller number of differentially expressed genes induced by belimumab treatment (13) are reported to be associated with poor response to this drug.

However, there are several limitations of the previous studies. The first one is that most of those studies focus on high disease activity status. Since patients with high disease activity were reported to show better responses to belimumab (8), non-responders should be enriched in moderate to low disease activity patients. This indicates that focusing on such patients has the potential power to elucidate factors associated with treatment resistance. The second limitation is that few of the previous studies focus on the genome-wide transcription difference between responders and non-responders. A comprehensive approach, such as comparing genome-wide expression profiles, has the potential power to detect unprecedented findings.

In this study, we assessed the treatment resistance factor of belimumab through clinical and molecular level profiling of moderate to low disease activity patients.

2 Methods

2.1 Study patients

The study design is shown in Figure 1A. We enrolled 46 patients diagnosed with SLE according to the European League Against Rheumatism/American College of Rheumatology classification criteria (The 2019 EULAR/ACR classification criteria) (14). Among the 46 patients, 38 were from Kyoto University Hospital, Kyoto, Japan, and the remaining eight were from Kobe City Medical Center General Hospital, Kobe, Japan. We excluded samples derived from two patients in which belimumab treatment was discontinued before six months due to a change of hospital (n=1) and withdrawal of treatment due to the patient’s request (n=1). We also enrolled healthy individuals (n=17) matched for age and sex (Supplementary Table 1).

Figure 1

Figure 1. (A) Schematic of study design. BLM: belimumab. Created with BioRender.com. Box plots comparing type I IFN score (B), the expression level of TNFSF13B (BAFF) (C), and TNFRSF13C (BAFF-R) (D) between SLE patients and healthy controls. CPM, Count per million (E–I) Box plots comparing each clinical index between responders and non-responders and its time course. "M0" stands for before treatment, "M3" for three, and "M6" for six months after treatment. "P(R vs. NR)" represents P-value comparing responders and non-responders, "P(ALL)" for the P-value comparing paired samples between 0M vs. 3M, or 3M vs. 6M. "P(R)" for P-value comparing 0M vs. 3M, 3M vs. 6M within responders, and "P(NR)" stands for the P-value obtained by comparison within non-responders. (J) Volcano plots showing differentially expressed genes between before treatment and three months after treatment. In (E–I), each dot represents each specimen.

2.2 Clinical evaluation and definition

At every clinic visit, disease activity was evaluated using the systemic lupus erythematosus disease activity index (SLEDAI-2K) (15). Clinical characteristics, including age, sex, treatment profiles (the use of prednisolone (PSL) and immunosuppressive agents), and clinical measurements (white blood cells, hemoglobin, complement) were obtained from medical records. The information on SLICC/ACR damage index (SDI) (15) and information on disease onset was obtained from Kyoto Lupus Cohort (16–18) for 38 samples in Kyoto University Hospital and from Kobe Lupus Cohort for 8 samples in Kobe City Medical Center General Hospital. Patients with a SLEDAI-2K score < 4 after six months of treatment with belimumab were defined as responders; those with a score of 4 or higher were classified as non-responders.

2.3 Sample collection, RNA Sequencing, transcriptome analysis

We collected blood samples from patients before and approximately 3 and 6 months after treatment, respectively. There were no missing samples from the before-treatment specimens for 44 patients; however, five samples were missing at the three-month time, and four were missing at the six-month time after treatment. Whole blood samples were stored in PAXgene tubes (QUIAGEN). RNA was extracted by PAXgene Blood RNA Kit (QUIAGEN). Library preparation was performed using TruSeq stranded Total RNA Library PrepKit with Ribo-Zero Globin Human. Sequencing was conducted by NovaSeq6000 in the 100-bp paired-end mode. Sequencing reads were trimmed using Trimmomatic ver. 0.36 (19) (leading: 20, trailing: 20, slidingwindow: 4:15, minlen: 36) and aligned to hg38 reference genome using STAR (ver. 2.7.3a) (20). Gene counts were generated by RSEM (ver. 1.3.1) (21) using Homo_sapiens.GRCh38.95.gtf from the Ensembl database (22). Gene counts were normalized by size factor implemented in DESeq2 (23) and converted to count per million (CPM), and log₂(CPM+1) was calculated. For the genome-wide differential gene expression test, we filtered low-count genes (0 count in ≥ 30% of 140 specimens, which include specimens from patients and healthy controls) and performed the Wald test using DESeq2. Type I IFN score was calculated by the mean expression of the four genes (IFI27, IFI44, IFI44L, and RSAD2) (6). Gene set enrichment analysis (24) was performed for 7,658 GO pathways in MSigDB (ver 7.5) (25) with 10,000 permutations. Overrepresentation analysis of gene sets was performed using Metascape (26) with default settings.

2.4 Immune cell type enrichment analysis

We performed an immune cell type enrichment analysis using CIBERSORTx (27) to estimate the proportion of each of the 22 cell types and the gene expression in each cell type in each specimen, utilizing the leukocyte signature matrix (LM22) (28) as a reference panel. Differential cellular proportion tests between responders and non-responders were performed on 13 cell subtypes (Table 1), excluding nine cell types (Mast cells activated, Dendritic cells resting, T cells follicular helper, T cells regulatory (Tregs), Eosinophils, Macrophages M0, Macrophages M1, Macrophages M2, and T cells gamma delta) whose estimated proportion was zero in over 95% of specimens (>133/140 specimens). We imputed expression of 57,189 genes (all genes quantified by RSEM) in memory B cells for each specimen using high resolution (‘hires’) function of CIBERSORTx, utilizing the same reference panel (LM22). Differential gene expression tests between responders and non-responders in memory B cells were performed on 4,154 genes, excluding 53,035 genes whose expression was estimated to be the same value or “NA” in all specimens. The Mann-Whitney test was used, and false discovery rate (FDR) was calculated by the Benjamini-Hochberg method.

Table 1

Table 1. Comparison between responders and non-responders.

2.5 Single-cell data analysis

To compare the expression levels of Ig within memory B cells across cell subpopulations, publicly available single-cell data of human peripheral blood mononuclear cells (PBMC) (29) were downloaded, comprising both CITE-Seq and scRNA-Seq data. Specifically, data from cells annotated as memory B cells (n=3,285) were utilized. Each immune subset was defined as follows: unswitched memory B cell (CD27⁺IgD⁺), switched memory B cell (CD27⁺IgD^-), non-activated memory B cell (CD27⁺CD38^lo), and activated memory B cell (CD27⁺CD38^hi). The definitions of CD27^+/- and IgD^+/- were based on whether cellular protein levels were nonzero or zero. Regarding CD38, data from two antibodies (“CD38-1” and “CD38-2”) were available, and we analyzed “CD38-1” data, as the data of CD38-1 and CD38-2 were well-correlated within memory B cells (data not shown), and no superiority was observed between them. The cutoff line for CD38-high (CD38^hi) or CD38-low (CD38^lo) cells was determined based on the median value (15) of protein levels in memory B cells. Ig expression levels were calculated as the average expression levels of IGHM, IGHA2, IGHG2, IGHA1, and IGHG1.

2.6 Statistical test

Unless otherwise specified, significant differences were tested using the Wilcoxon signed rank test for the two groups with paired samples and the Mann-Whitney test for the two groups with unpaired samples. The analysis of the usefulness of identified variables in distinguishing responders from non-responders employed receiver-operating-characteristic (ROC) curve techniques. Basically, we calculated the true/false positive rate and cutoff value to distinguish non-responders from responders; however, as for white blood cell count and C4, because the values were lower in non-responders, we calculated cutoff values to detect responders from non-responders. An area under the curve (AUC) below 0.6 was considered worthless, between 0.6 and 0.7 as poor, between 0.7 and 0.8 as fair, between 0.8 and 0.9 as good, and above 0.9 as excellent predictive accuracy. The statistical analysis and test were performed by the Python package sklearn (ver 0.23.2) or scipy (ver. 1.5.2).

3 Results

3.1 Patient characteristics and the effect of belimumab on the clinical and transcriptome feature

A total of 44 belimumab-naïve SLE patients were analyzed in this study (see Methods). Approximately 90% were female, with a median disease duration of 13.5 years (interquartile range (IQR) 6.0-23.5) (Supplementary Table 1). All patients were receiving GC treatment (PSL dose 7.5 (5.0-10.0) mg/day), with SLEDAI-2K indicating low to moderate disease activity (median 4.0, IQR 2.0-6.0). Additionally, to examine the transcriptomic features before belimumab administration, we compared the type I IFN score and the expression of TNFSF13B (BAFF) and TNFRSF13C (BAFF-R) with healthy individuals. The results showed higher type I IFN score and higher expression and TNFSF13B and lower expression of TNFRSF13C compared to healthy controls, consistent with previous reports (30, 31) and reflecting the characteristics of SLE (Figures 1B–D). Subsequently, we observed changes in clinical traits after belimumab administration in these patient groups (Figures 1E–I and Supplementary Figure 1). Initially, at three months post-administration, SLEDAI-2K decreased, while CH50 and C4 increased, and anti-DNA decreased within all the patients (Figures 1E, G, H; Supplementary Figure 1). However, from three to six months post-administration, there were no significant changes in these parameters. At six months, whether SLEDAI-2K < 4 or not divided the patient group into approximately halves (SLEDAI-2K < 4: 24/44, 54.5%) (Figure 1E). Therefore, we used SLEDAI-2K<4 at six months post-belimumab administration as the primary efficacy criterion (Responders: n=24, Non-responders: n=20). Non-responders had higher baseline SLEDAI-2K and lower SDI, and higher IgG and lower C3, C4, and CH50 levels among laboratory findings compared to responders (Table 1). We also tested post-administration changes unique to responders or non-responders. The responders showed a decrease of SLEDAI-2K and an increase of C4 at three months (p=0.009, 0.003, respectively, Figures 1E, G), while non-responders showed a reduction of anti-DNA antibody titers at three months (p=0.033, Figure 1H). However, the PSL dosages were reduced in all groups at all time points except for the non-responders at three months. We then observed changes in the transcriptome because of belimumab administration utilizing both responders and non-responders data, finding significant reductions in the expression levels of four genes (COL19A1, DEFA1B, TCL1A, and PCDH9) three months post-treatment compared to pre-treatment (FDR<0.05) (Figure 1J; Supplementary Table 2). We investigated the expression of the four genes among immune cell subtypes using the Database of Immune Expression (32) (Supplementary Figure 2). We found that COL19A1, TCL1A, and PCDH9 were genes with high expression levels in B cells. On the other hand, DEFA1B is poorly expressed in immune cells in physiological states, but it has been reported to be highly expressed in SLE patients (33). Similary, we observed changes between three months and six months utilizing both responders and non-responders, but there were no significant gene expression changes (Supplementary Table 3).

3.2 B cell proliferation pathway was upregulated in non-responders before treatment

Next, we conducted differential expression gene analysis utilizing pre-treatment data between responders and non-responders (Figure 2A; Supplementary Table 4). Significant differential expression was not observed in any genes, including TNFSF13B and TNFRSF13C (FDR>0.05, Supplementary Figure 3). However, when we tested the potential difference at the pathway level by gene set enrichment analysis (GSEA) using pre-treatment expression data, there was the most enrichment in the pathways related to B cell proliferation (Figure 2B). One of the most significant pathways, the “B cell proliferation” pathway, contained 100 genes, of which 7 had p-values below 0.05, all showing higher expression in non-responders (Figure 2C). The average expression of these 100 genes was significantly lower in both responders and non-responders compared to healthy controls (p < 0.05); however between responders and non-responders, responders showed higher expression with borderline significance (p=0.049, Figure 2D). This gene module expression was still higher in the non-responders than in the responders after six months (p=0.03). The seven genes shown in Figure 2C also showed similar dynamics during treatment (Supplementary Figure 4). We also found that the average expression of these seven genes could predict the treatment response with fair accuracy (AUC = 0.76, Figure 2E). We tested the association between the “B cell proliferation” pathway and various phenotypes and found negative correlation with PSL dose, white blood cell count, C4 and positive correlation with type I IFN pathway (Supplementary Table 5).

Figure 2

Figure 2. (A) Volcano plots showing different expressed genes between responders (R) and non-responders (NR) before treatment. (B) Results of gene set enrichment analysis before treatment. (C) Heatmap showing the expression of the seven genes in each sample. Expression was standardized across all samples, and the Z scores were shown. (D) Box plots comparing the average expression of the 100 genes included in the gene set "GOBP_B_CELL_PROLIFERATION" and their time course. Each dot represents each specimen. Red and blue boxes represent NR and R groups. As for the P values, please refer to the legend of Figure 1. (E) ROC curve for non-response to belimumab. Red plots indicate the cut-off point at the highest accuracy for predicting non-responders.

3.3 Memory B cells were abundant in non-responders before treatment

To investigate the immune cell subtypes associated with resistance to belimumab, we employed immune cell type enrichment analysis to estimate the proportions of each immune cell subset in each specimen. We have confirmed these estimates showed high correlations with clinical measurements (Supplementary Figure 5). Correlation analysis between B cell subtypes (memory B cell, naïve B cell, plasma cell) and the mean expression of the 100 genes in the B cell proliferation pathway utilizing SLE patient-derived specimens (n=123) showed the highest correlation with memory B cell (Supplementary Table 6). Therefore, we first compared the proportion of memory B cell types of pre-treatment data between responders and non-responders, finding a higher proportion in non-responders (p=5.1×10^-4) (Figure 3A). This association was also confirmed when a responder was defined as having a SLEDAI-2K score of 4 or less at six months, instead of the definition of less than 4 described in the Methods (Supplementary Figure 6A). We further divided the patient groups based on the decrease of SLEDAI-2K and PSL reduction at six months post-treatment, finding higher memory B cell proportions in patients with [(SLEDAI-2K at six months – SLEDAI-2K at baseline)] ≥ 0 and those unable to reduce the dosage of PSL (Figures 3B, C). Furthermore, the proportion of memory B cells before treatment exhibited a positive linear correlation with SLEDAI-2K increase (Figure 3D), and a positive linear correlation with [(PSL dose at 6 month) – (PSL dose at baseline)] (Supplementary Figure 6B), and [(PSL dose at 6 months)/(PSL dose at baseline)] (Supplementary Figure 6C), with suggestive significance. Overall, a higher proportion of memory B cells was identified as a treatment resistance factor, showing consistency across multiple response criteria.

Figure 3

Figure 3. (A–C). Box plots comparing the proportion of memory B cells to total white blood cells between responders and non-responders as defined in Methods (A), between patients with a decreased SLEDAI-2K [(SLEDAI-2K at six months) – (SLEDAI-2K before treatment) < 0] and those with an increased SLEDAI-2K (B), and between patients with a decreased prednisolone (PSL) dose ((PSL dose at six months) / (PSL dose before treatment) < 1) and those with an increased PSL dose (C). (D) Relation between memory B cell (%) before treatment and [(SLEDAI-2K at six months) – (SLEDAI-2K before treatment)] (E) ROC curve for non-response to belimumab. Red plots indicate the cut-off point at the highest accuracy for predicting non-responders. (F, G) Box plots comparing the proportion of memory B cells (E) and naïve B cells (F) between responders and non-responders and their time course. Each dot represents each specimen. Red and blue boxes represent NR and R groups. As for the P values, please refer to the legend of Figure 1.

To compare the results of memory B cells with other immune cell subtypes, we tested the difference between non-responders and responders for the remaining cell subtypes (Table 1). Although naïve CD4⁺ T cells showed significance (p =0.04) with higher expression in non-responders, this significance was lower than that observed for memory B cells. Therefore, it was suggested that an increase in memory B cells is a more prominent hallmark of treatment resistance in SLE patients. We found the proportion of memory B cells before treatment predicted the response with good accuracy (AUC=0.8) (Figure 3E). In order to search for alternative markers of memory B cells, correlations of the ratio of memory B cells with other traits, such as type I IFN score and clinical traits, were calculated using samples from enrolled patients. As a result, a positive correlation with type I IFN score and negative correlations with the number of white blood cells and complement C4 were observed (Supplementary Table 7). Among those indices, the number of white blood cells exhibited poor predictive accuracy; however, type I IFN score and C4 had fair predictive accuracy (Supplementary Figure 7).

3.4 Belimumab had no effect on the abundance of memory B cells in contrast to naïve B cells

Next, changes in the B cell subtype due to treatment with belimumab were examined. It was found that there were no significant changes in the proportion of memory B cells due to belimumab administration (Figure 3F). Non-responders exhibited the tendency of a higher proportion of memory B cells compared to responders across the treatment time course. Additionally, the impact of belimumab administration on naïve B cells and plasma cells was investigated. A significant decrease was observed in naïve B cells at three months of treatment, with no change observed from 3 to 6 months of treatment (Figure 3G). On the other hand, no significant changes were observed in plasma cells (Supplementary Figure 8). These findings indicate potential bias in the effects of belimumab treatment across B cell subtypes.

3.5 Upregulated genes in memory B cells of non-responders

To characterize the transcriptomic feature in memory B cells of non-responders, we estimated each gene’s expression levels in each specimen’s memory B cells using the immune cell type enrichment analysis and compared between non-responders and responders. Although none of the genes showed significance with FDR < 0.05, 225 genes showed nominal differences (p < 0.05), with 164 of them being higher in non-responders (Figure 4A; Supplementary Table 8). Enrichment analysis of these 164 genes revealed the strongest enrichment in the complement activation pathway (q-value = 0.02) (Figure 4B; Supplementary Table 9). It was found that immunoglobulin genes (IGHM, IGHA2, IGHG2, IGHA1, and IGHG1) contributed to the enrichment (Supplementary Table 9, Figure 4C), and all the five genes were upregulated in non-responders than in responders (Figure 4A, Supplementary Table 8). Hereafter, we call these five genes as Ig. We calculated the average expression levels of Ig in memory B cells for each sample and found that non-responders had higher levels of Ig compared to responders (Figure 4D). Although the expression levels decreased following belimumab administration, non-responders showed higher levels than responders throughout the treatment course, up to six months post-administration.

Figure 4

Figure 4. (A) Volcano plots showing differentially expressed genes between responders (R) and non-responders (NR) before treatment in memory B cells. Red dots represent p<0.05. (B) The results of enrichment analysis of genes upregulated in non-responders before treatment in memory B cells. (C) Heatmap showing the expression of the five immunoglobulin genes (IGHM, IGHA2, IGHG2, IGHA1, and IGHG1) in memory B cells in each sample. Expression was standardized across all samples, and the Z scores were shown. (D) Box plots comparing the average expression of the five immunoglobulin genes and their time course. Each dot represents each specimen. Red and blue boxes represent NR and R groups. As for the P values, please refer to the legend of Figure 1.

Finally, we aimed to identify cell populations that express Ig in memory B cells and their related cell types. We first attempted to perform analysis utilizing single-cell data derived from patients with SLE (34), but there weren’t a sufficient number of cells to perform analysis (data not shown). Alternatively, we used publicly available single-cell data of human PBMCs (Figure 5A) (29). When comparing populations between unswitched memory B cells (CD27⁺IgD⁺) and switched memory B cells (CD27⁺IgD^-), there was no difference in Ig expression (Figure 5B). However, Ig expression was higher in activated memory B cells (CD27⁺CD38^hi) than in non-activated memory B cells (CD27⁺CD38^lo) (Figure 5C). CD38 is a surface glycoprotein that serves as an activation marker on B cells (35). Furthermore, Ig expression levels showed a positive correlation with the expression of CD38 mRNA as well as CD38 protein (Figures 5D, E). Therefore, these results indicate that activated memory B cells are the primary source of Ig expression.

Figure 5

Figure 5. Analysis of publicly available single-cell RNA and CITE-Seq data. (A) Schematic of this analysis. (B, C) Comparison of Ig expression between switched and unswitched memory B cells (B) or between activated and non-activated memory B cells (C). (D, E). Association between the expression of Ig and CD38 mRNA levels (D) and CD38 protein levels (E) in memory B cells. The horizontal line represents the unique molecular identifiers (UMI) count of CD38 mRNA (D) and protein (E) derived from scRNA-seq and CITE-Seq data. Red lines indicate the line regressed by linear regression analysis.

4 Discussion

We identified transcriptomic and cellular characteristics related to the response of belimumab. The B cell proliferation pathway was upregulated in non-responders, which was associated with an increase in memory B cells. The proportion of those cells was found to be higher in non-responders before treatment. Belimumab did not induce significant changes in this transcriptomic and cellular feature during the treatment time course. More precisely, only 4 genes decreased their expressions after three months of treatment (Figure 1J) and no differentially expressed genes were detected between three months after and six months after treatment (Supplementary Table 3). Consequently, the difference in memory B cell proportion between responders and non-responders persisted after six months of treatment. In contrast, the number of naïve B cells decreased rapidly with belimumab administration in all patients (responders + non-responders, Figure 3G), consistent with previous reports (12, 36). These results suggest a potential bias in the effect of belimumab among B cell subtypes.

What drives memory B cell proliferation? We found memory B cells, B cell proliferative signals and type I IFN were all found to be positively correlated with each other (Supplementary Tables 5–7). This suggests that an increase in type I IFN may contribute to the activation of proliferative signaling and subsequent increase in memory B cells. Indeed, a recent report has shown that type I IFNs are associated with memory B cell proliferation in chronic infectious conditions (37).

We found the expression of Ig was higher in memory B cells in non-responders than in responders, which persisted after six months. Our data-driven approach defined “Ig” as a gene set that contains not only IGHG/A, markers for class switching, but also IGHM, a marker for non-class switching. This explains why the expression of such defined Ig did not differ between unswitched and switched memory B cells (Figure 5B). Our single-cell analysis using publicly available data indicated Ig is predominantly expressed by activated memory B cells, defined as CD27⁺CD38^hi (Figures 5C–E). These cells can also be referred to as plasmablasts. Together, our analysis suggests that belimumab may have limited effects on these activated memory B cells. Further studies are warranted to confirm the involvement of activated memory B cells and their Ig expression during belimumab treatment.

Based on the results of this study, it is crucial to establish criteria for predicting responsiveness to belimumab in clinical settings. We have discovered that a higher proportion of memory B cells predicts non-responsiveness with high accuracy (Figure 3E). However, promptly assessing the proportion of memory B cells in clinical settings is challenging. We explored the substitute markers and found a significant correlation between higher levels of memory B cells and lower levels of C4 and white blood cells, and lower levels of C4 predict the treatment response with fair accuracy (Supplementary Figure 7).

Then, what alternative treatments could be considered for patients who were predicted to be poor responders to belimumab? Rituximab may help decrease the number of memory B cells. Additionally, since the number of memory B cells showed a positive correlation with the type I IFN score, anifrolumab (anti-type I IFN receptor monoclonal antibody) could be one of the candidates. Furthermore, although not currently approved for SLE, Janus kinase (JAK) inhibitors are also potential candidates since JAK is one of the essential components of IFN signaling (38, 39). They would modulate B cell activation, differentiation, and functions by reducing STAT3 signaling, leading to decreased switched memory B cell formation and decreased immunoglobulin production, as shown in the previous in-vitro and clinical study (40). Lastly, continuous use of belimumab, irrespective of initial treatment response, may aid in reducing memory B cells, as reported in a study demonstrating that seven-year treatment depleted the number of memory B cells (41). Further clinical research is warranted to establish a treatment strategy for SLE with higher number of memory B cells.

There are several conflicting results compared to previous reports. Firstly, it was found that patients with higher disease activity before treatment showed poor responsiveness to belimumab (Table 1), which contradicts previous findings (8). Differences in the disease activity of the analyzed patients could partly explain this. While the previous report focused on patients with high disease activity, our study targeted patients with low to moderate activity. Another potential explanation is the difference in treatment response definition. While the previous study employed the SLE Responder index, which includes the changes in SLEDAI, the present study employed the criterion of responders when SLEDAI-2K < 4 after six months. We will discuss this response criteria issue further in the next paragraph. Secondly, patients with high SDI showed good responses, contrary to previous findings (9). This could be explained by the confounding effect of disease duration; responders tended to show more extended disease duration periods (Table 1). Third, in contrast to several reports indicating a transient increase in memory B cells after Belimumab treatment initiation (42, 43), there was no significant change in the proportion of memory B cells (Figure 3F). However, this phenomenon could not be replicated in another study (12), and there are suggestive increases, especially in responders, after 3 months of treatment initiation.

Based on the primary definition of response criteria in our study (SLEDAI-2K <4 at six months), responders showed low disease activity before treatment with minimal changes in clinical indices during the six months (Figures 1E–I and Table 1). Additionally, some non-responders showed improvements in certain parameters, such as increased C3 levels and decreased anti-dsDNA antibody titers (Figures 1F, H and Table 1). These observations reflect an inherent challenge in SLE clinical research: defining treatment response in patients with low to moderate disease activity is particularly challenging due to the limited margin for improvement in conventional metrics. Despite these limitations in response definition, our key finding regarding memory B cell proportions remains consistent. The association between higher memory B cell proportions and non-response was confirmed across multiple response criteria, including definitions based on SLEDAI-2K changes (Figure 3B) or PSL dosage reductions (Figure 3C). Furthermore, we demonstrated a linear correlation between memory B cell proportions and SLEDAI-2K changes (Figure 3D), supporting the biological relevance of this association regardless of response definition. In Figure 2D, the average expression levels of the 100 genes in the B cell proliferation pathway were higher in healthy controls compared to both non-responders and responders (SLE patients). This could be attributed to the effects of treatments, such as steroids or immunosuppressants. Although non-responders showed higher expression levels in the B cell proliferation pathway compared to responders, their levels remained lower than those of healthy controls.

There are several limitations to this study. Firstly, the proportion of immune cell subset and gene expression in each cell subset were not experimentally determined, such as through flow cytometry or RNA-Seq in each cell subtype but were based on estimation by cell type enrichment analysis. Nevertheless, this estimation showed a good correlation with clinical measurements (Supplementary Figure 5). Secondly, our cohort followed patients for six months, but longer-term observation may yield additional perspectives, as discussed above (41). Thirdly, none of the Ig components showed significance at the conservative threshold (FDR <0.05). Therefore, the contribution of each gene contribution should be interpreted with caution. However, the enrichment analysis, which evaluates the genes as a ranked collection rather than individual components, still yields statistically significant results (q-value =0.02, Figure 4B; Supplementary Table S9). Fourth, our analysis using single-cell RNA-seq and CITE-seq data was derived from non-SLE subjects, which warrants further study to identify the cellular source of Ig expression. Another limitation is the relatively small sample size compared with previous clinical studies. Considering the heterogeneity of SLE, confirmation of the obtained results in more extensive and independent cohorts is warranted.

Data availability statement

The gene expression dataset supporting the conclusions of this article is available in the Zenodo repository, https://doi.org/10.5281/zenodo.14557188. These data can be downloaded without restriction.

Ethics statement

The studies involving humans were approved by Kyoto University Graduate School and Faculty of Medicine Ethics Committee (approval number: R2591), and Kobe City Medical Center General Hospital (approval number: 20272). The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.

Author contributions

TI: Writing – original draft, Writing – review & editing, Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Supervision, Validation, Visualization. HY: Conceptualization, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Visualization, Writing – original draft, Writing – review & editing. KK: Conceptualization, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Visualization, Writing – original draft, Writing – review & editing. SS: Writing – original draft, Writing – review & editing. SA: Writing – original draft, Writing – review & editing. RN: Writing – original draft, Writing – review & editing. HT: Writing – original draft, Writing – review & editing. RH: Writing – original draft, Writing – review & editing. MS: Writing – original draft, Writing – review & editing. KM: Writing – original draft, Writing – review & editing. AO: Writing – original draft, Writing – review & editing. HO: Writing – original draft, Writing – review & editing. MT: Writing – original draft, Writing – review & editing. FM: Writing – original draft, Writing – review & editing. AM: Writing – original draft, Writing – review & editing. KO: Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Supervision, Validation, Writing – original draft, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This study received funding from GlaxoSmithKline (GSK). However, GSK had no additional roles in the study design, data collection, data analysis, the decision to publish, or the preparation of the manuscript.

Acknowledgments

We thank all patients and physicians involved in the care of patients. We extend our gratitude to Ms. Chikako Endo for her technical assistance. We also thank Dr. Hiroshi Doi for helping with the establishment of the sample collection system and Ms. Masayo Motokawa and Ms. Aki Shimizu for maintaining the system.

Conflict of interest

KO (until March 2021) and HY (from April 2021) received research funding from GSK as principal investigators for a research project, which includes the present study. HT received research grants from the GSK, outside the work. RH received research grants and speaker fees from GSK, outside the work. HO, AO, and MT belong to the department financially supported by Asahi Kasei, outside the work. AM received research grants and speaker fees from Chugai, Astellas Pharma, and Asahi Kasei, outside the work. KO received speaker fees from AstraZeneca and GSK, outside the work.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2025.1506298/full#supplementary-material

References

1. Kaul A, Gordon C, Crow MK, Touma Z, Urowitz MB, van Vollenhoven R, et al. Systemic lupus erythematosus. Nat Rev Dis Primers. (2016) 2:1–21. doi: 10.1038/nrdp.2016.39

PubMed Abstract | Crossref Full Text | Google Scholar

2. Cheema GS, Roschke V, Hilbert DM, Stohl W. Elevated serum B lymphocyte stimulator levels in patients with systemic immune-based rheumatic diseases. Arthritis Rheum. (2001) 44:1313–9. doi: 10.1002/1529-0131(200106)44:6<1313::AID-ART223>3.0.CO;2-S

PubMed Abstract | Crossref Full Text | Google Scholar

3. Moore PA, Belvedere O, Orr A, Pieri K, LaFleur DW, Feng P, et al. BLyS: member of the tumor necrosis factor family and B lymphocyte stimulator. Science. (1999) 285:260–3. doi: 10.1126/science.285.5425.260

PubMed Abstract | Crossref Full Text | Google Scholar

4. Baker KP, Edwards BM, Main SH, Choi GH, Wager RE, Halpern WG, et al. Generation and characterization of LymphoStat-B, a human monoclonal antibody that antagonizes the bioactivities of B lymphocyte stimulator. Arthritis Rheum. (2003) 48:3253–65. doi: 10.1002/art.11299

PubMed Abstract | Crossref Full Text | Google Scholar

5. Halpern WG, Lappin P, Zanardi T, Cai W, Corcoran M, Zhong J, et al. Chronic administration of belimumab, a BLyS antagonist, decreases tissue and peripheral blood B-lymphocyte populations in cynomolgus monkeys: pharmacokinetic, pharmacodynamic, and toxicologic effects. Toxicol Sci. (2006) 91:586–99. doi: 10.1093/toxsci/kfj148

PubMed Abstract | Crossref Full Text | Google Scholar

6. Furie R, Khamashta M, Merrill JT, Werth VP, Kalunian K, Brohawn P, et al. Anifrolumab, an anti-interferon-α Receptor monoclonal antibody, in moderate-to-severe systemic lupus erythematosus. Arthritis Rheumatol. (2017) 69:376–86. doi: 10.1002/art.39962

PubMed Abstract | Crossref Full Text | Google Scholar

7. Fanouriakis A, Kostopoulou M, Andersen J, Aringer M, Arnaud L, Bae S-C, et al. EULAR recommendations for the management of systemic lupus erythematosus: 2023 update. Ann Rheum Dis. (2024) 83:15–29. doi: 10.1136/ard-2023-224762

PubMed Abstract | Crossref Full Text | Google Scholar

8. van Vollenhoven RF, Petri MA, Cervera R, Roth DA, Ji BN, Kleoudis CS, et al. Belimumab in the treatment of systemic lupus erythematosus: high disease activity predictors of response. Ann Rheumatic Dis. (2012) 71:1343–9. doi: 10.1136/annrheumdis-2011-200937

PubMed Abstract | Crossref Full Text | Google Scholar

9. Parodis I, Gomez A, Emamikia S, Chatzidionysiou K. Established organ damage reduces belimumab efficacy in systemic lupus erythematosus. Ann Rheum Dis. (2019) 78:1006–7. doi: 10.1136/annrheumdis-2018-214880

PubMed Abstract | Crossref Full Text | Google Scholar

10. Wilkinson C, Henderson RB, Jones-Leone AR, Flint SM, Lennon M, Levy RA, et al. The role of baseline BLyS levels and type 1 interferon-inducible gene signature status in determining belimumab response in systemic lupus erythematosus: a post hoc meta-analysis. Arthritis Res Ther. (2020) 22:102. doi: 10.1186/s13075-020-02177-0

PubMed Abstract | Crossref Full Text | Google Scholar

11. Roth DA, Thompson A, Tang Y, Hammer AE, Molta CT, Gordon D. Elevated BLyS levels in patients with systemic lupus erythematosus: Associated factors and responses to belimumab. Lupus. (2016) 25:346–54. doi: 10.1177/0961203315604909

PubMed Abstract | Crossref Full Text | Google Scholar

12. Ramsköld D, Parodis I, Lakshmikanth T, Sippl N, Khademi M, Chen Y, et al. B cell alterations during BAFF inhibition with belimumab in SLE. EBioMedicine. (2019) 40:517–27. doi: 10.1016/j.ebiom.2018.12.035

PubMed Abstract | Crossref Full Text | Google Scholar

13. Nakano M, Ota M, Takeshima Y, Iwasaki Y, Hatano H, Nagafuchi Y, et al. Distinct transcriptome architectures underlying lupus establishment and exacerbation. Cell. (2022) 185:3375–3389.e21. doi: 10.1016/j.cell.2022.07.021

PubMed Abstract | Crossref Full Text | Google Scholar

14. Aringer M, Costenbader K, Daikh D, Brinks R, Mosca M, Ramsey-Goldman R, et al. 2019 European League Against Rheumatism/American College of Rheumatology classification criteria for systemic lupus erythematosus. Ann Rheumatic Dis. (2019) 78:1151–9. doi: 10.1136/annrheumdis-2018-214819

PubMed Abstract | Crossref Full Text | Google Scholar

15. Gladman D, Ginzler E, Goldsmith C, Fortin P, Liang M, Urowitz M, et al. The development and initial validation of the Systemic Lupus International Collaborating Clinics/American College of Rheumatology damage index for systemic lupus erythematosus. Arthritis Rheum. (1996) 39:363–9. doi: 10.1002/art.1780390303

PubMed Abstract | Crossref Full Text | Google Scholar

16. Takase Y, Iwasaki T, Doi H, Tsuji H, Hashimoto M, Ueno K, et al. Correlation between irreversible organ damage and the quality of life of patients with systemic lupus erythematosus: The Kyoto Lupus Cohort survey. Lupus. (2021) 30:1577–85. doi: 10.1177/09612033211025614

PubMed Abstract | Crossref Full Text | Google Scholar

17. Doi H, Ohmura K, Tabuchi Y, Hashimoto M, Takase Y, Inaba R, et al. Validation and verification of the Japanese version of the systemic lupus erythematosus symptom checklist for patient quality of life. Lupus. (2021) 30:1108–15. doi: 10.1177/09612033211005026

PubMed Abstract | Crossref Full Text | Google Scholar

18. Iwasaki T, Doi H, Tsuji H, Tabuchi Y, Hashimoto M, Kitagori K, et al. Phenotypic landscape of systemic lupus erythematosus: An analysis of the Kyoto Lupus Cohort. Mod Rheumatol. (2022) 32:571–6. doi: 10.1093/mr/roab020

PubMed Abstract | Crossref Full Text | Google Scholar

19. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. (2014) 30:2114–20. doi: 10.1093/bioinformatics/btu170

PubMed Abstract | Crossref Full Text | Google Scholar

20. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. (2013) 29:15–21. doi: 10.1093/bioinformatics/bts635

PubMed Abstract | Crossref Full Text | Google Scholar

21. Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinf. (2011) 12:323. doi: 10.1186/1471-2105-12-323

PubMed Abstract | Crossref Full Text | Google Scholar

22. Martin FJ, Amode MR, Aneja A, Austine-Orimoloye O, Azov AG, Barnes I, et al. Ensembl 2023. Nucleic Acids Res. (2023) 51:D933–41. doi: 10.1093/nar/gkac958

PubMed Abstract | Crossref Full Text | Google Scholar

23. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. (2014) 15:550. doi: 10.1186/s13059-014-0550-8

PubMed Abstract | Crossref Full Text | Google Scholar

24. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U.S.A. (2005) 102:15545–50. doi: 10.1073/pnas.0506580102

PubMed Abstract | Crossref Full Text | Google Scholar

25. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst. (2015) 1:417–25. doi: 10.1016/j.cels.2015.12.004

PubMed Abstract | Crossref Full Text | Google Scholar

26. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. (2019) 10:1523. doi: 10.1038/s41467-019-09234-6

PubMed Abstract | Crossref Full Text | Google Scholar

27. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol. (2019) 37:773–82. doi: 10.1038/s41587-019-0114-2

PubMed Abstract | Crossref Full Text | Google Scholar

28. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. (2015) 12:453–7. doi: 10.1038/nmeth.3337

PubMed Abstract | Crossref Full Text | Google Scholar

29. Hao Y, Hao S, Andersen-Nissen E, Mauck WM, Zheng S, Butler A, et al. Integrated analysis of multimodal single-cell data. Cell. (2021) 184:3573–3587.e29. doi: 10.1016/j.cell.2021.04.048

PubMed Abstract | Crossref Full Text | Google Scholar

30. Möckel T, Basta F, Weinmann-Menke J, Schwarting A. B cell activating factor (BAFF): Structure, functions, autoimmunity and clinical implications in Systemic Lupus Erythematosus (SLE). Autoimmun Rev. (2021) 20:102736. doi: 10.1016/j.autrev.2020.102736

PubMed Abstract | Crossref Full Text | Google Scholar

31. Carter RH, Zhao H, Liu X, Pelletier M, Chatham W, Kimberly R, et al. Expression and occupancy of BAFF-R on B cells in systemic lupus erythematosus. Arthritis Rheum. (2005) 52:3943–54. doi: 10.1002/art.21489

PubMed Abstract | Crossref Full Text | Google Scholar

32. Schmiedel BJ, Singh D, Madrigal A, Valdovino-Gonzalez AG, White BM, Zapardiel-Gonzalo J, et al. Impact of genetic polymorphisms on human immune cell gene expression. Cell. (2018) 175:1701–1715.e16. doi: 10.1016/j.cell.2018.10.022

PubMed Abstract | Crossref Full Text | Google Scholar

33. Sthoeger ZM, Bezalel S, Chapnik N, Asher I, Froy O. High α-defensin levels in patients with systemic lupus erythematosus. Immunology. (2009) 127:116–22. doi: 10.1111/j.1365-2567.2008.02997.x

PubMed Abstract | Crossref Full Text | Google Scholar

34. Arazi A, Rao DA, Berthier CC, Davidson A, Liu Y, Hoover PJ, et al. The immune cell landscape in kidneys of patients with lupus nephritis. Nat Immunol. (2019) 20:902–14. doi: 10.1038/s41590-019-0398-x

PubMed Abstract | Crossref Full Text | Google Scholar

35. Piedra-Quintero ZL, Wilson Z, Nava P, Guerau-de-Arellano M. CD38: an immunomodulatory molecule in inflammation and autoimmunity. Front Immunol. (2020) 11:597959. doi: 10.3389/fimmu.2020.597959

PubMed Abstract | Crossref Full Text | Google Scholar

36. Jacobi AM, Huang W, Wang T, Freimuth W, Sanz I, Furie R, et al. The effect of prolonged treatment with belimumab on B cells in human SLE. Arthritis Rheum. (2010) 62:201–10. doi: 10.1002/art.27189

PubMed Abstract | Crossref Full Text | Google Scholar

37. Cooper L, Xu H, Polmear J, Kealy L, Szeto C, Pang ES, et al. Type I interferons induce an epigenetically distinct memory B cell subset in chronic viral infection. Immunity. (2024) 57:1037–1055.e6. doi: 10.1016/j.immuni.2024.03.016

PubMed Abstract | Crossref Full Text | Google Scholar

38. Banerjee S, Biehl A, Gadina M, Hasni S, Schwartz DM. JAK-STAT signaling as a target for inflammatory and autoimmune diseases: current and future prospects. Drugs. (2017) 77:521–46. doi: 10.1007/s40265-017-0701-9

PubMed Abstract | Crossref Full Text | Google Scholar

39. Platanias LC. Mechanisms of type-I- and type-II-interferon-mediated signalling. Nat Rev Immunol. (2005) 5:375–86. doi: 10.1038/nri1604

PubMed Abstract | Crossref Full Text | Google Scholar

40. Frede N, Lorenzetti R, Hüppe JM, Janowska I, Troilo A, Schleyer M-T, et al. JAK inhibitors differentially modulate B cell activation, maturation and function: A comparative analysis of five JAK inhibitors in an in-vitro B cell differentiation model and in patients with rheumatoid arthritis. Front Immunol. (2023) 14:1087986. doi: 10.3389/fimmu.2023.1087986

PubMed Abstract | Crossref Full Text | Google Scholar

41. Huang W, Quach TD, Dascalu C, Liu Z, Leung T, Byrne-Steele M, et al. Belimumab promotes negative selection of activated autoreactive B cells in systemic lupus erythematosus patients. JCI Insight. (2018) 3:e122525. doi: 10.1172/jci.insight.122525

PubMed Abstract | Crossref Full Text | Google Scholar

42. Stohl W, Hiepe F, Latinis KM, Thomas M, Scheinberg MA, Clarke A, et al. Belimumab reduces autoantibodies, normalizes low complement levels, and reduces select B cell populations in patients with systemic lupus erythematosus. Arthritis Rheum. (2012) 64:2328–37. doi: 10.1002/art.34400

PubMed Abstract | Crossref Full Text | Google Scholar

43. Mariette X, Barone F, Baldini C, Bootsma H, Clark KL, Vita SD, et al. A randomized, phase II study of sequential belimumab and rituximab in primary Sjögren’s syndrome. JCI Insight. (2022) 7:e163030. doi: 10.1172/jci.insight.163030

PubMed Abstract | Crossref Full Text | Google Scholar