This study included patients diagnosed with advanced, recurrent, or metastatic NSCLC who received anti-PD-1 (nivolumab or pembrolizumab) or anti-PD-L1 therapy (atezolizumab), either as monotherapy or in combination with chemotherapy, at Kurume University Hospital (Kurume, Japan) and Kanagawa Cancer Center (Yokohama, Japan). The patient cohort partially overlapped with that of a previous study conducted by our group (21, 22). The cohort consisted of 179 patients enrolled between February 2016 and August 2019. The clinical course of the enrolled patients was followed until July 2024. The patient characteristics are summarized in Figures 1B–D , and Table 1 . Tumor PD-L1 expression was assessed via immunohistochemical staining of paraffin-embedded tumor sections using anti-PD-L1 monoclonal antibodies (clone E1L3N: Cell Signaling Technology, Danvers, MA, USA and clone 22C3: Agilent Technologies/Dako, Carpinteria, CA, USA). For most patients, PD-L1 expression was evaluated in tumor specimens collected prior to first-line therapy. This analysis was also done in patients who received ICIs as second- or subsequent-line treatment. Clinical response was evaluated based on the Response Evaluation Criteria in Solid Tumors (RECIST) version 1.1. IrAEs were defined according to the Common Terminology Criteria for Adverse Events (CTCAE) v5.0 grading scale. This study was conducted in accordance with the principles of the Declaration of Helsinki and was approved by the Institutional Review Boards of Kurume University Hospital (approval numbers: 15210 and 19240) and Kanagawa Cancer Center (approval number: 28-85). Informed consent was obtained from all enrolled patients after the nature of the study and its possible consequences were explained.

Peripheral blood samples were collected in heparin-coated tubes before and 6 weeks after therapy initiation. Plasma was separated from whole blood via centrifugation and stored at −80 °C until analysis. Plasma autoantibodies were assessed via the AlphaScreen assay using two human protein arrays generated using a wheat germ cell-free protein production system (23, 24). The HuPEX protein array, containing 19,713 human proteins, was purchased from CellFree Science Co., Ltd. (Yokohama, Japan). The Ehime-Kazusa protein array, consisting of 4,144 human proteins, was prepared in-house using a cell-free protein synthesis reagent. Gene resources for the protein array were provided by the Kazusa DNA Research Institute (Kisarazu, Japan) (25). Each cDNA clone was subcloned into the pEU-E01-FLAG-GST-K1–02 vector (26). In vitro transcription and translation were performed using the WEPRO7240 Expression Kit (CellFree Science) as previously described (27, 28). Briefly, the DNA fragments coding for each protein were amplified via PCR and used as templates for in vitro transcription. The mRNA generated by in vitro transcription was used as a template for in vitro translation (total reaction mixture, 5 µL; components, 2.5 µL mRNA; 1.67 µL WEPRO 7240 wheat germ extract; 0.14 µL creatine kinase [20 mg/mL]; and 0.11 µL RNase inhibitor) was layered below 50 µL SUB-AMIX SGC substrate solution in a 384-well plate and incubated at 20°C for 18 h. The translated protein arrays were diluted 2-fold with AlphaScreen buffer (100 mM Tris-HCl, pH 8.0; 0.01% Tween 20; 1 mg/mL BSA), aliquoted into 384-well plates, flash-frozen in liquid nitrogen, and stored at −80°C.

The AlphaScreen assay was performed as previously described, with minor modifications (27). Specificity was achieved through two independent high-affinity binding steps (1): glutathione-coated donor beads that selectively bind GST-tagged recombinant antigens, and (2) Protein G–conjugated acceptor beads that specifically recognize human IgG, thereby minimizing nonspecific interactions. For primary screening, autoantibody reactivity was assessed after incubating 10 paired plasma samples from patients with NSCLC (five pre- and post-ICI treatment, respectively) with a library of 23,857 human proteins. Plasma samples were diluted 1:667 in AlphaScreen reaction buffer (100 mM Tris-HCl, pH 8.0, 0.01% [v/v] Tween 20, 0.1% [w/v] bovine serum albumin), and 20 µL diluted plasma was dispensed into each well of OptiPlate-384 plates (Revvity, Yokohama, Japan) using a Liquidator96 pipetting system. Next, 1 µL of each FLAG-GST tagged protein was transferred from 384-well stock plate to the reaction plate using JANUS automated dispensing workstation equipped with a NanoHead 384-channel microsyringe head (Revvity). Subsequently, 9 µL detection mixture (containing 0.06 µL AlphaScreen GSH Donor Beads [Revvity] and 0.06 µL Protein G-conjugated AlphaScreen acceptor beads in reaction buffer) was added to the reaction plates with a FlexDrop dispenser (27). After incubation at 25°C for 1 h in the dark, AlphaScreen signals were detected using an EnVision plate reader (Revvity).

For secondary screening, 358 paired plasma samples from 179 patients with NSCLC (pre- and post-treatment) were tested for the 59 selected proteins, each in quadruplicate. All reactions were performed on AlphaPlate-1536 plates (Revvity). Three microliters of protein dilution containing 0.05 µL FLAG-GST tagged protein was dispensed into an AlphaPlate 1536 plate (Revvity) via a Multidrop Combi nL (Thermo Fisher Scientific, Tokyo, Japan). All 358 plasma samples were diluted 1:40 and dispensed into a 384 well plate. Subsequently, 0.2 µL of the 1:40 diluted plasma from each well of these 384-well plates was dispensed into the 1536-well reaction plates via a JANUS workstation with a NanoHead. Finally, 1.8 µL detection mixture (containing 0.01 µL AlphaScreen GSH Donor Beads and 0.01 µL Protein G-conjugated AlphaScreen acceptor beads) was added to the reaction plates via a Multidrop Combi nL. After incubation at 25°C for 1 h in the dark, AlphaScreen signals were detected with an EnVision plate reader. The secondary screening minimized inter-plate variability by ensuring that paired plasma samples (pre- and post-ICI) from the same patient were analyzed on the same plate. To further minimize experimental variability, all assays were conducted on a single 1536-well plate uniformly coated with the same antigen. Plasma samples from all patients were dispensed in quadruplicate. This design enabled direct assessment of treatment-induced changes both within individual patients and across the cohort, while controlling for inter-plate variability.

All analyses were performed using R (version 4.4.1; https://www.r-project.org). Univariate binary logistic and Cox proportional hazards regression analyses were performed using the glm function (R Basic package). P-values were derived using Wald z-tests and subsequently adjusted for multiple comparisons using the Benjamini–Hochberg false discovery rate (FDR) procedure. To mitigate overfitting in multivariate predictive modeling, we applied elastic net (EN) regularized regression and random forest (RF) regression for binary outcome variables (irAE, ir-pneumonitis, and treatment response), and Cox proportional hazards regression with elastic net regularization (EN-Cox) and random survival forest (RSF) for survival outcomes (PFS and OS) within a machine learning framework (21). The cohort was randomly stratified into training (80%) and validation (20%) groups. Model development and hyperparameter optimization were conducted using nested cross-validation (CV) exclusively on the training set. The final models were subsequently retrained on the entire training set and independently evaluated on the validation set, maintaining a strict separation between the development and evaluation datasets to prevent data leakage. For EN regression, hyperparameter tuning (α = 0.5 for balanced L1/L2 regularization) was performed using 10-fold CV via the cv.glmnet implementation (glmnet package), optimizing the regularization parameter λ to minimize classification error. The λ value yielding minimum CV error was selected for final model training, followed by external validation on the holdout set. For the RF regression, feature selection was performed via sequential variable inclusion based on importance rankings, with variable importance assessed using the mean Gini index (randomForest package). The optimal feature subset (important factors) was identified by maximizing predictive accuracy during stepwise variable addition. Model performance was assessed using leave-one-out cross-validation (LOO-CV) of the training set with discriminatory power quantified by the area under the receiver operating characteristic curve (AUC). All RF implementations utilized an ensemble of 1,000 decorrelated decision trees with the minimum node purity set to 1, following established methodologies in high-dimensional regression (29). Optimization of the model cutoff value was performed by maximizing Youden’s J statistic (sensitivity + specificity - 1), implemented through the coords function in R. Permutational multivariate analysis of variance (PERMANOVA) was performed using the Adonis function (vegan package). Nonparametric comparisons were performed using the Wilcoxon signed-rank test via the wilcox.test function. Survival analyses included Kaplan–Meier curve generation with log-rank testing, alongside Cox proportional hazards modeling, executed using the survival and survminer packages in R.

As shown in Figure 1 and Table 1 , of the 179 patients with NSCLC included in this study, 128 (71.5%) received anti-PD-1/PD-L1 monotherapy and 51 (28.5%) received combination therapy with chemotherapy. The best overall response was partial response (PR) in 56 patients (31.3%), stable disease (SD) in 37 patients (20.7%), and progressive disease (PD) in 86 patients (48.0%). IrAEs occurred in 69 patients (38.6%) and ir-pneumonitis was observed in 21 patients (11.7%). The median PFS in the cohort was 131 days and the median OS was 468 days. Figure 1C illustrates PFS across the clinical subgroups. No significant differences in PFS were observed between patients who received combination therapy and those who received monotherapy (hazard ratio [HR] = 0.95, p = 0.786). Patients in the responder group (PR + SD) demonstrated significantly prolonged PFS compared to those in the PD group (HR = 0.12, p < 0.001). Similarly, patients who developed irAEs had a longer PFS than those who did not (HR = 0.52, p < 0.001). However, patients with ir-pneumonitis did not show a significant difference in PFS compared to those without ir-pneumonitis (HR = 0.74, p = 0.225). Figure 1D shows OS data across the same subgroups. Consistent with the PFS results, the OS did not differ significantly between the ICI monotherapy and combination therapy (HR = 0.92, p = 0.683) groups. The responder group exhibited significantly improved OS compared with the PD group (HR = 0.27, p < 0.001). Additionally, patients with irAEs had a longer OS than those without irAEs (HR = 0.69, p = 0.0377). In contrast, patients who developed ir-pneumonitis showed no OS benefit from treatment (HR = 1.19, p = 0.498).

For primary screening, we used AlphaScreen to quantify plasma IgG antibodies against approximately 24,000 biotinylated proteins in pre- and post-treatment plasma from five patients with NSCLC who developed ir-pneumonitis during ICI treatment. Based on pre/post signal ratios, 37 antibodies showed significant changes (p < 0.05, Student’s t-test; fold changes < 0.7 or > 1.4). To offset the limited sample size, 22 additional antibodies were included due to notable fold changes and prior associations with immune-related pathways, cancer biology, or ICI responses reported in the literature. In total, 59 antibodies were selected for further analysis ( Supplementary Table S1 ).

In secondary screening, these 59 candidates were quantified in pre- and post-treatment plasma from a formal cohort of 179 patients. A volcano plot ( Figure 2A ) illustrates changes in CAAB levels before and after anti-PD-1/PD-L1 therapy: 53 decreased (average log₂-transformed fold-change < 0), 30 significantly (p < 0.05), including BASP1, ITGAE, CPB1, and CD200 autoantibodies. Six increased (average log₂-transformed fold-change > 0), with two autoantibodies against CASP10 and SCGN reached statistical significance.

We calculated the log₂-transformed fold-change for each CAAB, designating this dataset as the dynamic repertoire. Clinical relevance of the CAAB dynamic repertoire by computing multivariate effect sizes for treatment regimen, treatment line, PD-L1 expression, occurrence of irAE and ir-pneumonitis, therapeutic response, PFS, and OS using PERMANOVA (30). As shown in Figure 2B , the CAAB dynamic repertoire exhibited a significant association only with chemotherapy exposure among all clinical factors assessed (highest R² = 2.5%, p < 0.05), indicating marked sensitivity to anticancer agent administration. Regarding clinical baseline characteristics, CAAB dynamics demonstrated a stronger effect size with PD-L1 expression than with the line of therapy, although without significance. In the present cohort, the therapeutic decisions (monotherapy versus combination therapy) were directly influenced by PD-L1 expression ( Figure 2C ) (31). Among the five treatment outcome measures—irAEs, ir-pneumonitis, RECIST-assessed response, PFS, and OS—the effect size of the association with CAAB dynamics decreased in the following order: irAEs, ir-pneumonitis, response, PFS, and OS, although none reached statistical significance ( Figure 2B ). To further elucidate which CAABs exhibited the strongest association with concomitant anticancer agent administration during ICI therapy, we performed univariate logistic regression analyses stratified by anticancer agent use. CAABs with a significance threshold of p < 0.05 in the regression analyses are summarized in Figure 2D . Nine CAABs retained statistical significance after FDR correction: antibodies against ZNF551, S100A7A, KRT8, CASP10, IRF4, CCDC3, DBX1, SCGN, and KLHL18. Several other CAABs exhibited notable associations, including PKDCC, SOX1, UBTF, BIRC3, WDR24, IL7, BASP1, CPB1, IL4, FBXL6, LAMA4, CXCL1, CSF2, GATAD2B, PFKFB4, and F13B. Plasma levels of these autoantibodies were significantly modulated by anticancer agent coadministration. Violin plots ( Figure 2D ) illustrate distinct expression patterns. ZNF551-targeting CAABs decreased after monotherapy but increased with combination therapy. Similarly, S100A7A and KRT8 CAABs decreased with monotherapy, but remained stable before and after combination therapy. CASP10 autoantibody levels increased in both groups, albeit less markedly in monotherapy.

Through the above analyses, we identified that the treatment regimen was the primary significant factor influencing the dynamic repertoire of CAAB during ICI therapy. Given the non-negligible influence of anticancer agents, we further examined the association between CAAB fluctuation patterns and treatment outcomes in a subgroup of patients receiving ICI monotherapy. To identify the autoantibodies specifically associated with each treatment outcome (irAEs, ir-pneumonitis, treatment response, PFS, and OS), we performed individual univariate logistic/Cox regression analyses for each factor. Figure 3A shows that in patients receiving ICI monotherapy, CAABs targeting DBX1, BIRC3, BASP1, STAT4, and SPATC1L were significantly associated with the occurrence of irAEs. The irAE-positive subgroup exhibited more pronounced post-treatment reductions in autoantibody titers than the controls. Figure 3B presents the univariate regression analyses of pneumonitis-associated CAABs, revealing significant correlations for autoantibodies targeting CXCL2, ROPN1, SPATC1L, FURIN, and DBX1. Consistent with the irAE pattern, patients with ir-pneumonitis showed greater treatment-related decreases in CAABs, with CXCL2 showing the greatest reduction. Regarding the therapeutic efficacy, response-associated CAABs included autoantibodies against SNCA, CCDC3, ECSCR, STAT4, and BIRC3 ( Figure 3C ). The treatment responders displayed a significantly greater absolute reduction in CAAB levels after ICI administration. PFS- and OS-associated CAABs included autoantibodies against F13B, ECSCR, and CXCL1 ( Figures 3D, E ). The Venn diagram in Figure 3F illustrates both shared and unique CAABs across different treatment outcomes. Autoantibodies against BIRC3 and STAT4 were common to irAEs and treatment responses, whereas SPATC1L and DBX1 were associated with both irAEs and ir-pneumonitis. Treatment response was characterized by distinct autoantibodies targeting SNCA, CCDC3, and ECSCR, with ECSCR additionally shared with PFS. Autoantibodies against CXCL1 were specifically associated with OS, while ir-pneumonitis showed unique associations with CXCL2, ROPN1, and FURIN, highlighting the differential CAAB dynamics underlying these outcomes.

Next, we assessed the potential of CAAB as a predictive biomarker for treatment outcomes associated with anti-PD-1/PD-L1–based ICI therapy. The predictive performance of CAAB dynamics for five treatment outcomes across the overall cohort and six subgroups defined by three stratification variables is summarized in Table 2 . As shown, models using the entire cohort without stratification demonstrated limited predictive accuracy. In contrast, stratified analyses revealed notable performance: for irAE prediction, the EN model in the high PD-L1 subgroup attained AUCs of 0.81 and 0.68 for the training and test sets, respectively ( Figure 4 ); for ir-pneumonitis prediction, the RF model in the monotherapy subgroup achieved AUCs of 0.83 and 0.88 ( Figure 5 ); for response prediction, the EN model in the combination therapy subgroup exhibited AUCs of 0.90 and 0.79; and in the first-line treatment subgroup, the RF model achieved AUCs of 0.79 and 0.89 ( Figure 6 ). By contrast, the predictive capacity of CAAB dynamics for PFS and OS remained limited.

In the high PD-L1 subgroup, the EN model selected five CAABs (RNASE2, FGFR3, PKDCC, DBX1, and CD40). The optimized cutoff value of 0.35 yielded a specificity of 0.86 and a sensitivity of 0.72 for predicting irAEs ( Figure 4 ). the RF model with the top three features (PARK2, KLHL18, CXCL2) achieved the best performance, with a cutoff of 0.09 (specificity 0.86, sensitivity 0.75) for predicting ir-pneumonitis ( Figure 5 ). For treatment response in the combination therapy subgroup, the EN model highlighted nine CAABs (HTATSF1, NKIRAS1, IL17A, ITGBL1, LAMA4, IRF4, BACH2, HES1, PKDCC), with a cutoff of 0.49 (specificity 0.79, sensitivity 0.91) ( Figures 6A–C ). In the first-line treatment subgroup, the RF model selected seven CAABs (CCDC3, STAT4, ECSCR, PFKFB4, F13B, SNCA, FURIN), with a cutoff of 0.55 (specificity 0.92, sensitivity 0.65) ( Figures 6D–F ). Both models demonstrated favorable predictive performance for treatment response within their respective subgroups.

These results indicate that patient background and treatment options distinctly shape CAAB dynamics, and that capturing subset-specific dynamic signatures may enable effective prediction of treatment outcomes across different clinical contexts or therapeutic regimens.