CD86 Is Associated with Immune Infiltration and Immunotherapy Signatures in AML and Promotes Its Progression

Background Cluster of differentiation 86 (CD86), also known as B7-2, is a molecule expressed on antigen-presenting cells that provides the costimulatory signals required for T cell activation and survival. CD86 binds to two ligands on the surface of T cells: the antigen CD28 and cytotoxic T lymphocyte-associated protein 4 (CTLA-4). By binding to CD28, CD86—together with CD80—promotes the participation of T cells in the antigen presentation process. However, the interrelationships among CD86, immunotherapy, and immune infiltration in acute myeloid leukemia (AML) are unclear. Methods The immunological effects of CD86 in various cancers (including on chemokines, immunostimulators, MHC, and receptors) were evaluated through a pan-cancer analysis using TCGA and GEO databases. The relationship between CD86 expression and mononucleotide variation, gene copy number variation, methylation, immune checkpoint blockers (ICBs), and T-cell inflammation score in AML was subsequently examined. ESTIMATE and limma packages were used to identify genes at the intersection of CD86 with StromalScore and ImmuneScore. Subsequently, GO/KEGG and PPI network analyses were performed. The immune risk score (IRS) model was constructed, and the validation set was used for verification. The predictive value was compared with the TIDE score. Results CD86 was overexpressed in many cancers, and its overexpression was associated with a poor prognosis. CD86 expression was positively correlated with the expression of CTLA4, PDCD1LG2, IDO1, HAVCR2, and other genes and negatively correlated with CD86 methylation. The expression of CD86 in AML cell lines was detected by QRT-PCR and Western blot, and the results showed that CD86 was overexpressed in AML cell lines. Immune infiltration assays showed that CD86 expression was positively correlated with CD8 T cell, Dendritic cell, macrophage, NK cell, and Th1_cell and also with immune examination site, immune regulation, immunotherapy response, and TIICs. ssGSEA showed that CD86 was enriched in immune-related pathways, and CD86 expression was correlated with mutations in the genes RB1, ERBB2, and FANCC, which are associated with responses to radiotherapy and chemotherapy. The IRS score performed better than the TIDE website score. Conclusion CD86 appears to participate in immune invasion in AML and is an important player in the tumor microenvironment in this malignancy. At the same time, the IRS score developed by us has a good effect and may provide some support for the diagnosis of AML. Thus, CD86 may serve as a potential target for AML immunotherapy.


Introduction
Acute myeloid leukemia (AML) is a common hematological disease characterized by the clonal proliferation, abnormal diferentiation, and cell death evasion of bone-marrow-derived hematopoietic stem and progenitor cells [1]. Tese cells proliferate in the peripheral blood and infltrate the bone marrow. Te tumor microenvironment in AML is characterized by immunosuppression, which promotes immune tolerance and the immune escape of malignant cells [2]. Te main components of the AML bone marrow microenvironment (BMM) include T cells, B cells, and NK cells [3]. Te immune imbalance of T helper cells (T cells) is a major contributor to the sudden progression of AML [4].
T-cell-mediated cellular immunity is primarily achieved by the specifc binding of antigen peptides to the major histocompatibility complex (MHC) (frst signal) and the binding of costimulatory molecules located on the surface of antigen-presenting cells (APCs) to their receptors (second signal) [5]. Te absence of costimulatory molecules leads to immune unresponsiveness, which promotes tumor escape in AML. Owing to advancements in research, immunotherapies that retarget efector cells (T cells, NK cells) have been developed, and these have become the key for AML treatment [6]. Typically, tumors suppress the immune system, resulting in the impairment of T-cell function. Te goal of immunotherapy is to eliminate this impairment. Studies have shown that vaccines for AML/dendritic cell fusion can amplify T-cell populations and prevent AML recurrence [7]. Terefore, efective immunotherapy approaches that target specifc proteins are the key to AML treatment.
B7-2, better known as CD86, is a member of the B7 family [8]. CD28 and cytotoxic T lymphocyte antigen 4 (CTL1-4) are regulated. CD86 can bind to CD28, leading to signal production and the recognition of antigenic peptides by T-cell receptors (TCRs), which leads to T-cell proliferation and IL-2 production [9]. CD86 has been reported to be overexpressed in samples from AML patients [10]. CD86 is a marker for monocytes and dendritic cells and is involved in the progression of AML [11]. Improvements in sequencing technology have promoted extensive research on molecular networks using gene sequencing data from public databases [12]. However, the correlation between CD86 and immunomodulators (chemokines, receptors, and MHC proteins), immunotherapy results, and immune checkpoint proteins in AML has not been reported. Terefore, it is very important to explore the associations among CD86-related molecules, immune infltration, and immunotherapy.

Data.
Te Cancer Genome Atlas (TCGA) data: Pancarcinoma (33 species) RNA sequencing (RNA-SEQ) data (FPKM values) were downloaded from the UCSC Xena data portal (https://xenabrowser.net/). Tey were converted to TPM format, and somatic mutation data and survival information were downloaded. Log2 transformation was performed on the RNA-SEQ data, and somatic mutation data were analyzed using MuTect. Copy number variation (CNV) data processed using GISTIC were downloaded from the UCSC Xena data portal. Further, the methylation data were downloaded from the LinkedOmics data portal.
Ten, information was obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/) LAML GEO queue, which contains detailed survival data. Te information including data from the GSE10358, GSE37642 (including data from the GPL570 and GPL96 platforms), GSE146173, GSE106291, and GSE12417 databases (including data from the GPL97 and GPL96 platforms). Te sample data for leukemia was retained. Further, three GEO databases containing information on responses to immunotherapy were downloaded: GSE78220 (melanoma), GSE135222 (nonsmall-cell lung cancer), and GSE91061 (melanoma). Moreover, complete expression data and detailed clinical information of patients from the IMvigor210 study (AML immunotherapy-related data) were obtained from https:// research-pub.Gene.com/imvigor210corebiologies/under the 3.0 license.

Analysis of Immunological Characteristics of AML.
First, using the web portal TISIDB (https://cis.hku.hk/ TISIDB) [13], genes related to the immune response, including those encoding immune stimulants, MHC proteins, immune receptors, and chemokines, were identifed. 'ggplot2' in R software was used for visualization, and the R package 'Corrlot' was used to calculate the Spearman correlation coefcients between the expression of CD86 and that of the abovementioned genes. In order to calculate the correlation between CD86 expression and that of various oncogenes expressed in the tumor microenvironment, single-sample gene enrichment analysis (ssGSEA) was performed and the correlation between CD86 and immune cell scores was calculated. Te association of CD86 with immune risk scores (IRSs) and the infammatory coefcient of T cells was calculated using the generalized T-cell infammation score formula [14].
Ten R package 'limma' was used to analyze the differences in the expression of chemokines, immunostimulators, MHC proteins, and immune receptors based on high vs. low CD86 expression. CIBERSORT, MCPcounter, TIMER, Quantiseq, and Xcell were used to examine the immune-infltrating cells in AML. Te correlation between CD86 and common immune checkpoint blockers (ICBs) was calculated. Further, StromalScore and ImmuneScore were calculated for AML samples using the R package 'ESTI-MATE.' 'Limma was used to identify the DEGs in the high vs. low CD86 expression, StromalScore, and ImmuneScore groups. Ten, 'ggplot2' was used to draw volcano maps and heat maps of the DEGs. A total of 308 up-regulated genes and 16 down-regulated genes were identifed through this analysis.

Immune Risk Score (IRS) Calculation.
IRSs were calculated based on the time of patient enrollment. Te 324 DEGs were randomly sampled from TCGA to establish the training and validation sets at a 1 : 1 ratio. Te R package 'SurvMiner' was used to conduct univariate Cox regression analysis for the DEGs, and the optimal characteristic genes were identifed according to the Least Absolute Shrinkage and Selection Operator (LASSO) method. Multivariate Cox regression analysis was performed, and based on the median IRS, the sample was divided into groups. Te Kaplan-Meier method was used to compare survival outcomes between these groups. Univariate Cox analysis was used to screen PRGs with a prognostic value. Te P value threshold for signifcance was set at 0.05, and 17 survival-related genes were selected for further analysis. LASSO-penalized Cox regression analysis and GLMNET R software package were used to establish a prognostic model to reduce overftting. Finally, six genes and their coefcients were retained to determine the penalty parameter (λ) with the minimum criterion. Te risk score was obtained using the formula IRS � n i�1 β * xi, where β � the regression coefcient. AML patients were divided into two groups: high-risk group and low-risk group. Te 'SurvMiner' R software package was used to compare survival status between the two risk groups, and 'Survival' and 'timeROC R software packages were used for receiver operating characteristic (ROC) curve analysis. In addition, univariate and multivariate Cox regressions were used to determine the independent prognostic value of the three genes. To verify the validity of the model, analyses were performed using data from the GEO internal test queue or ICGC external validation queue. Median risk scores were obtained using the GEO training cohort, whereas patients in the GEO test cohort were divided into low-and high-risk groups.

GO, KEGG and PPI Analysis.
Based on CD86 expression, StromalScore, and ImmuneScore, the patients were divided into two groups. Using the limma package and subsequent fltering based on a |Log2FC| ≥ 1 and FDR < 0.05, DEGs were identifed in the high vs. low CD86 expression, StromalScore, and ImmuneScore groups. Te 'cluster analyzer' R package was used for Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis. Normalized P values < 0.05 and an FDR q < 0.05 were considered statistically signifcant.
PPI analysis was performed for the 324 DERs using STRING, Cytoscape was used for visualization, and the MCODE plug-in was used to identify critical clusters.

Western Blot Assay.
After cell digestion and centrifugation, RIPA lysis bufer (Beyotime Biotech) was added, and the cells were lysed on ice for 30 min. Ten, cells were centrifuged at 12000 rpm for 30 min. Te supernatant was removed, and protein levels were quantifed using the BCA kit (Beyotime Biotech). Te proteins were separated using SDS-PAGE and electrotransferred to PVDF membranes. Te CD86 primary antibody (Proteintech) was incubated overnight at 4°C after 2 hours of rapid blocking solution (BSA; Beyotime Biotech). On the following day, the corresponding secondary antibody was added. Protein bands were detected using the ECL exposure solution.
2.9. Statistical Analysis. Data were plotted using R package (V 4.0.0). Te T test and Utest were used to compare variables between two groups. Categorical variables were evaluated using the Chi-square test. Pearson and Spearman coefcients were used for correlation analysis. Te Kaplan-Meier method was used to plot survival outcomes, and the logarithmic rank sum test was used to analyze statistical diferences. P < 0.05 was considered statistically signifcant.

CD86 Is Overexpressed in Many Cancers and Is
Associated with the Prognosis and Immune Response of AML. Using TCGA data on the expression profles of 33 cancers, CD86 expression was examined. Te fndings showed that CD86 was highly expressed in most of the cancers, such as breast cancer, cholangiocarcinoma, colorectal cancer, esophageal cancer, glioma, renal clear cell carcinoma, renal papillary cell carcinoma, lung adenocarcinoma, lung squamous cell carcinoma, pancreatic cancer, rectal adenocarcinoma, gastric cancer, thyroid cancer, and endometrial cancer (Supplementary Figure 1). CD86 was also overexpressed in AML (Supplementary Figure 1). Ten, based on the median expression value of CD86, patients were divided into high-and low-expression groups.
Kaplan-Meier analysis was performed to examine high vs. low CD86 expression in various cancers using TCGA data, and log-rank tests were used for survival analyses. Te results showed that low CD86 expression was associated with bladder urothelial carcinoma, cervical squamous cell carcinoma, and endocervical adenocarcinoma. In AML, low expression of CD86 had a statistically signifcant better prognosis (Supplementary Figure 2). Te results of univariate Cox regression analysis were then used to create a forest map, which showed that CD86 expression was statistically signifcant in various cancers (Supplementary Figure 3A). Subsequently, using the TISIDB website established by Ru et al., four gene sets-chemokines, immunostimulators, MHC proteins, and receptors-were downloaded (Supplementary Table 1). Spearman correlation coefcients were used to analyze the association between CD86 and these four gene sets in diferent cancer types (Figure 1(a)). Subsequently, the correlation between key pancarcinoma molecules (including CTLA4, PDCD1LG2, IDO1, and HAVCR2) and CD86 was calculated. Tese genes were found to be positively correlated with CD86 in AML (Figures 1(b)-1(d)). ssGSEA method was then used to evaluate the scores of 28 immune cell types in diferent cancer types, and then calculated the correlation between CD86 and them. Te results showed that CD86 expression was positively correlated with 28 types of immune cells (Figure 1(f)).

Single Nucleotide Variation (SNV), Gene Copy Number Variation (CNV), and Methylation Analysis of CD86 in AML.
Site mutations are a key pathogenic factor causing abnormal proliferation in AML. To investigate whether CD86 is mutated in AML, SNV, and gene CNV data for AML were analyzed. Te results showed that CD86 was not mutated in AML. Te AML samples were divided into two groups according to a CD86-expression-based cutof. Te group with high CD86 expression had a higher risk, indicating that the high CD86 level was a risk factor for leukemia (Figure 2(a)). Ten, the 10 genes with the highest mutation frequencies in the high-vs. low-expression groups were plotted. Accordingly, we found that DNMT3A, FLT3, NPM1, IDH2, and other genes had a relatively high mutation frequency in the low expression group (Figure 2(b)). Differences in tumor mutation load (TMB) were examined in the CD86 high-vs. low-expression groups, but the results revealed no signifcant diferences (Figure 2(c)). Te amplifcation and deletion of CD86 was examined. However, most samples showed no copy number changes in the CD86 gene ( Figure 2(d)). Te expression of the CD86 gene was compared across diferent groups. Meanwhile, the correlation between the expression of CD86 and the degree of methylation was calculated and plotted. CD86 expression showed a signifcant negative correlation with CD86 methylation (Figure 2(e)). All previous experiments were conducted using public databases. To validate whether CD86 is associated with AML, we examined CD86 expression in vitro. QRT-PCR and Western blot were used to detect CD86 expression in SKM-1 (myelodysplastic syndrome), OCI-AML2 (human myeloid leukemia cell), SH-1 (human myeloid leukemia cell), HL-60 (human myeloid leukemia cell), MEG01 (human megakaryoblastic leukemia cell), and K562 cells (human myeloid leukemia cell). Te results showed that CD86 was overexpressed in OCI-AML2, THP-1, SH-1, and K652 cells (Figures 2(f ) and 2(g)). Tese results demonstrated that while CD86 was not mutated in AML and was not related to the TMB, the degree of CD86 methylation decreased with an increase in CD86 expression.

Immune Status of CD86 High-vs. Low-Expression Groups in AML.
To further understand the association between CD86 expression and immunoassay sites in AML, the diferences in chemokine, immunostimulator, MHC protein, and immune receptor expression were compared between the high vs. low CD86 expression groups (Supplementary Table 2). A heat map was drawn to represent the DEGs (Figure 3(a)). Te distribution of 28 types of immune cells in the high vs. low CD86 expression groups was analyzed. Te results showed that for 24 types of immune cells, the group with the high expression of CD86 had a higher immune score (Figure 3(b)). To further understand the correlation between CD86 expression and tumorinfltrating immune cells (TIICs) in AML, CIBERSORT, MCPcounter, TIMER, Quantiseq, and Xcell were used. Immune infltration analysis was performed, and correlation between CD86 expression and immune scores was calculated. Further, given that CD8+ T cell recruitment and dendritic cell, macrophage, NK cell, and T1 activation are required during the migration of immune cells to tumors, the marker genes of these cell types were analyzed in the CD86 high-vs. low-expression groups (Supplementary Table 3). Te heat map is shown in Figure 3(c). In addition, the correlation between CD86 and immune checkpoints was calculated. Te results indicated that CD86 was positively correlated with these aforementioned immunoassay sites (Figure 3(d)).

IRS Model Construction and Verifcation.
After a series of analyses, 324 DERs were identifed. Subsequently, 65 prognostic genes were obtained through random sampling based on TCGA samples (training: test � 1 : 1) and univariate Cox regression (P < 0.05, Supplementary Table 5). Ten, the LASSO method was used to select the best genes, and six genes were obtained according to the minimum lambda cutof of 0.1452 (Figure 7(a)). Multivariate Cox regression analysis was performed using these six genes, and the risk coefcients of related genes were obtained and represented by a forest map (Figure 7(b)). Ten, the risk score of each sample in the TCGA training and verifcation datasets was calculated. Te samples were divided into two groups (high vs. low expression) based on the best cutof, and Kaplan-Meier curves were drawn. Further, ROC curve analysis was also performed. Te results showed that the low-expression group had a better survival prognosis (Figures 7(c) and 7(d)). Subsequently, IRS model validation was performed using all TCGA datasets, GSE10358 datasets, and GSE37642 (GPL570) datasets. Te results showed that patients in the low-expression group had a good prognosis (Figures 7(e)-7(g)). In order to further verify the accuracy of the IRS, evaluations were performed using the following GEO datasets: GSE146173, GSE106291, GSE37642 (a subset of the GPL96 platform), GSE12417 (a subset of the GPL97 platform), and GSE12417 (a subset of the GPL96 platform).
Here too, the results revealed a better prognosis in the lowrisk group (Supplementary Figures 7A-7E).

Association between IRS and Immunity.
Meanwhile, based on the TCGA dataset, we compared the diferential expression of high and low IRC expression groups and concentration of chemokine, immunostimulator, MHC, and receptor genes. Tese were represented by heat maps (Figure 8(a)). Te diference in CD86 expression between the high-and low-expression groups was detected (Figure 8(b)). Analyses of infammation scores for pan-cancer T cells revealed signifcantly higher scores in the high-expression group (Figure 8(c)). Subsequently, we plotted the correlation between IRC and 28 types of immune cells using ssGSEA method. Te results suggested that the high-expression group was enriched for a variety of immune cells (Figure 8(d)). Diferences in immunoassay sites and IRC groupings were also examined (Figure 8(e)). Te results suggested that a high IRS is correlated with immune cells in AML.

Performance Comparison between IRS and TIDE.
To verify the efect of the IRS model constructed by us, we collected data from the IMvigor210, GSE91061, GSE78220, and GSE135222 datasets after immunotherapy. We used our method to calculate the IRS, and the TIDE website was used to evaluate the TIDE score (https://tide.dfci.harvard.edu/) for immune treatment efects. Te predictive value of the IRS and TIDE for the response to treatment was then compared. Survival prediction curves and Kaplan-Meier curves (median cutof) were used for analysis. Our IRS score was found to be better than the TIDE score (Figures 9(a)-9(k)).
A molecule can be a central target for cancer immunotherapy depending on its specifc expression in the tumor microenvironment. CD86 (B7-2), a member of the B7 family of proteins, is one of the surface proteins of APCs [22]. Te B7 family has been implicated in the progression of AML. Te levels of CD80 (B7-1) are elevated in AML [23]. Moreover, programmed cell death ligand (PD-L1, B7H-1) is abnormally expressed in AML patients and is directly associated with a poor prognosis [24]. T cells can be activated to exert immune efects only when CD86 is expressed on APC membranes and binds to CD28 on the surface of T cells [11]. Using data from public databases, we found that CD86 is overexpressed in many cancers, and especially in AML. We also demonstrated this in AML cell lines. In AML, a high expression of CD86 was found to be associated with a poor prognosis. Further, clinical data from GEO and TCGA datasets show that high CD86 expression is directly associated with a poor prognosis in AML.
In AML, mutation sites not only afect disease classifcation but also afect risk stratifcation and chemotherapeutic resistance. For example, FLT3 mutations are detected in about one-third of AML patients, and these mutations are directly related to the poor prognosis of AML [33]. However, interestingly, the mutation rates of DNMT3A, FLT3, NPM1, and IDH2 were higher in the low CD86 expression group in our study [34,35]. Te mutation rate of RUNX1 was higher in the CD86 group, which could be because of the number of samples. Our study also showed that CD86 expression was negatively correlated with DNA methylation. Tis was noteworthy because methylation has been found to predict chemotherapy outcomes in AML [36]. Meanwhile, we predicted that mutations in RB1, ERBB2, and FANCC increased as CD86 expression increased, suggesting that CD86 may be related to radiotherapy and chemotherapy resistance in AML. However, further verifcation is still needed.
Te IRS is a genetic prognostic model calculated using a formula to assess the risk of a disease. An IRS can predict the survival and prognosis of AML patients undergoing chemotherapy. Immune risk scores can be used to predict the beneft of adjuvant chemotherapy in diferent risk groups of patients. Wang Yun et al. [37] constructed diferent algorithms to evaluate the prognostic models of AML immune components after receiving diferent degrees of radiotherapy, and the results were highly accurate. However, in cases of AML, the IRS is currently inaccurate and inconsistent. Hence, we developed an IRS model to predict the overall prognosis of AML. Verifcation with external datasets showed that our model is superior to the TIDE score. Tis complements the enrichment of AML risk scores.
Nevertheless, there are some limitations to our study. First, all our samples were obtained from public databases, and a large number of patient samples are still needed for follow-up verifcation. Second, no in vivo experiments or mechanistic studies were performed. Tis area needs to be explored further.

Conclusion
Tis study found that CD86 is involved in the progression of AML and is closely related to the BMM in AML. Te expression of CD86 could be used to predict immunotherapy efcacy. Terefore, the development of CD86-targeting drugs could lead to advancements in AML treatment.

Data Availability
Te data used to support the fndings of this study are included within the article.

Conflicts of Interest
Te authors declare that they have no conficts of interest.
Diferential genes in StromalScore, ImmuneScore and CD86 high-and low-expression groups. Supplement Table 5. 65 genes screened out by Cox regression analysis. (Supplementary Materials)