We downloaded GSE142255 and GSE156382 from the Gene Expression Omnibus (http://www.ncbi.nlm.nih.gov/geo/). The GSE142255 dataset, which was extracted on the GPL17586 platform (Affymetrix Human Transcriptome Array 2.0), comprises transcriptome analysis of whole blood samples from 17 ACLF patients, 7 acute decompensated cirrhosis patients, 7 stable cirrhosis, 5 fibrosis patients and 7 healthy individuals (6). Meanwhile the dataset GSE156382, which was derived on the GPL20844 platform (Agilent-072363 SurePrint G3 Human GE v3 8x60K Microarray 039494), includes circulating neutrophil samples from 6 chronic liver disease, 12 ACLF samples and 5 healthy samples (13). Further analysis were all performed on the data from ACLF patients and healthy controls. GSE142255 was used as the training set, while GSE156382 was utilized as verifying datasets.

The DEGs between patients and normal controls in the GSE142255 dataset were screened by R software “limma” package. DEG identification threshold values were set at | logFC) |>1 and P < 0.05. Volcano maps and heatmaps were plotted using the R software packages ggplots and pheatmap. Through published literature, 137 neutrophil-related genes (Ne-related genes, NRGs) were attained to screen for neutrophil-related DEGs (14). These genes were selected based on their roles in neutrophil function, infection and immune response.

To specifically investigate the role of neutrophil in ACLF, we focused on genes that were both DEGs in ACLF patients and related to neutrophil. WGCNA was performed on the samples and data in the training dataset using the R package WGCNA to identify key modules, select an appropriate soft threshold (β), and adjust the scale-free topology fitting index to exceed 0.9. Based on Pearson correlation analysis, the module most relevant to disease traits was selected as the key module, and genes from this key module were extracted as ModelGenes.

To obtain differentially expressed NRGs related to ACLF (ACLF-NRGs), we took the intersection of DEGs, NRG, and module genes (ModelGenes) related to ACLF. We aimed to pinpoint the most relevant genes that could serve as specific markers or therapeutic targets related to neutrophil function in the context of ACLF. This targeted approach allows us to highlight the molecular mechanisms and pathways directly involved in neutrophil function, which might be overlooked if all DEGs were analyzed collectively. Network module analysis is vulnerable to the impact of outlier samples; therefore, we employed a hierarchical clustering tree to identify and eliminate these outliers from the dataset.

ACLF-NRGs were mapped onto human chromosomes by using the RCircos package of R language. We performed Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) functional enrichment analysis on ACLF-NRGs using the ClusterProfiler package, with p.adj<0.05 and count ≥ 2 as the threshold for screening. The Cytoscape STRING database was used to analyze the protein-protein interaction (PPI) networks of ACLF-NRGs.

Random Forest (RF) analysis was performed on ACLF-NRGs by using the random Forest package. The optimal number of classification trees and gene importance ranking were obtained. The selection criteria for this analysis was to include genes with mean decrease accuracy that is not zero. To eliminate the recursive features of the obtained differential genes, the SVM algorithm in R package e1071was employed. The RF method was used to sort and calculate the above feature genes, determining the significance and ranking of each gene, along with the error rate and accuracy for each iterative combination. The combination yielding the lowest error rate was identified as the optimal configuration, and the corresponding gene was designated as the feature gene. The intersection of genes obtained by RF and SVM-RFE algorithms was used to obtain feature genes.

Lasso analysis was conducted on feature genes employing the glmnet package to obtain error maps, gene co-efficient maps for cross validation and diagnostic genes. The optimal model parameter λ=0.008889174 was selected using 10 fold cross validation (the λ corresponding to the smallest error mean). These diagnostic genes were used to construct diagnostic models in the ACLF and Control groups to evaluate the diagnostic ability of those key genes for ACLF. The Receiver Operating Characteristic (ROC) analysis was conducted on this model to assess its judgment accuracy using training set GSE142255. The lasso diagnostic model constructed by diagnostic genes was validated in datasets GSE156382.

To better understand the molecular mechanisms of diagnostic genes acting on ACLF, signaling pathways affected significantly by diagnostic genes were identified using Ingenuity Pathway Analysis (IPA). We performed enrichment analysis on all differentially expressed genes as input, displayed the pathways in which the diagnostic genes were involved, and observed whether the enriched pathways contained disease-related pathways. “c2. cp. kegg. v11.0. symbols” was downloaded as the internal reference gene set from MSigDB database. Then Gene Set Enrichment Analysis (GSEA) was conducted on the diagnostic genes.

Through published literature, ferroptosis-related genes (FRGs) and cuprotosis-related genes (Cuprotosis) were obtained. ssGSEA was used to analyze the ferroptosis ssGSEA score (FRGs Score) and cuprotosis ssGSEA score (Cuprotosis Score) of the ferroptosis-related gene set and cuprotosis-related gene set on GSE142255 samples.

Enzyme-Linked Immunosorbent Assay were conducted on 116 patients recruited as outpatients or inpatients at the First Affiliated Hospital of Zhejiang University between August 2015 and September 2017. These patients were divided into four groups: 82 chronic hepatitis B (CHB) patients with ACLF (HB-ACLF), 8 CHB patients with acute decompensated cirrhosis (HB-AD), 11 age-/sex-matched CHB and 15 age-/sex-matched healthy subjects. RT-qPCR and western blot analysis was carried out on white blood cells from a different cohort of 54 subjects, who were recruited from inpatients of the First Affiliated Hospital of Wenzhou Medical University between July 2023 and December 2024. The Asian-Pacific Association for the Study of the Liver (APASL) defines ACLF as follows: “acute hepatic insult manifesting as jaundice (≥5 mg/dL) and coagulopathy, complicated within 4 weeks by ascites and/or encephalopathy in a patient with previously diagnosed or undiagnosed chronic liver disease.” (15) The AD patients presented with cirrhosis acutely complicated with hepatic encephalopathy, ascites and upper gastrointestinal bleeding. The diagnosis of cirrhosis was based on liver biopsy, radiological evidence, or clinical manifestation of liver decompensation. CHB is diagnosed when serum HBsAg is positive for at least six months, coupled with cirrhosis or liver fibrosis or long-term liver inflammation confirmed by histology or imaging or laboratorial or clinical evidence; Control subjects were 18 age-/sex-matched healthy individuals. We excluded patients who were infected with HIV, younger than 18 years old, pregnant, using immunotherapy, had a history of autoimmune diseases, cancer or complicated with serious diseases of other major organs. Each subject provided written consent. If participants were unable to give informed consent, their next of kin provided written consent. The Helsinki Declaration was followed in the conduct of this study. Approval for this study was obtained from the Ethic Committee of the First Affiliated Hospital of Zhejiang University and the First Affiliated Hospital of Wenzhou Medical University. A summary of the baseline characteristics of the patients can be found in the Supplementary Tables S1 and S2 .

Human blood samples were lysed with Red Blood Cell Lysis Buffer to obtain white blood cells. RNA extraction was conducted in accordance with the protocol provided by the TRIzol reagent manufacturer (Takara, Dalian, China). Then, the RNA was reverse transcribed into cDNA using the Primescript reverse transcription reagent kit for qPCR (Takara, Dalian, China). Finally, reverse transcription quantitative polymerase chain reaction (RT-qPCR) experiments were performed utilizing the SYBR Green Premix Pro Taq HS qPCR Kit (Takara, Dalian, China). The primers for MMP9, S100A12, CCL5 and beta-actin were offered by Meigu Biotechnology Co., Ltd. (Wenzhou, China).

We conducted Western blot (WB) analysis as previously described (16). Briefly, 20 μg of protein from the human blood white blood cells was subjected to SDS-PAGE using a 12% gel. After trans-membrane, the membranes were blocked with rapid blocking buffer for ten minutes at room temperature. Subsequently, the membranes were incubated overnight at 4°C with primary antibodies targeting MMP9, S100A12 and β-actin followed by incubation with secondary antibodies. At last, enhanced chemiluminescence were used to detect the protein levels. All the antibodies were purchased from Proteintech, USA.

After enrolling in the study, whole blood samples were collected from participants. As soon as plasma was collected after centrifugation, it was stored at -80°C immediately. In accordance with the manufacturer’s instructions, enzyme-linked immunosorbent assays from R&D and biolegend were used to determine MMP9 and CCL5 levels in plasma samples (10μL) respectively.

Small molecule drugs for diagnostic genes were predicted using the Drug Bank database. The NetworkAnalyst database was searched for transcription factors (TFs) that regulate diagnostic genes (https://www.networkanalyst.ca/). Then the upstream kinases of the identified TFs were predicted by X2Kgui, and miRNAs associated with the diagnostic genes were predicted by using miRDB database. At last, the Kinase-TF-mRNA-miRNA network was constructed by Cytoscape software.

GraphPad Prism 6.0 or IBM SPSS statistics 20 were used to generate and draw the statistical analysis and charts. Image J software was employed to determine the gray value of the WB. There are four ways to express data: means ± standard deviation (SD), means ± standard error of the mean (SEM), number (percentage) and median (range). For the comparison of two independent groups, Mann-Whitney U tests were used. The correlation between FRGs Score or Cuprotosis Score and the expression of the diagnostic genes were evaluated by Pearson correlation analysis. Kaplan Meier survival curves were also developed to evaluated whether plasma MMP9 levels were highly correlated with mortality. Receiver operating characteristic curves (ROC curve) were employed to calculate diagnostic capability. In all statistical tests, p less than 0.05 was considered statistically significant.

Study protocol was shown in Figure 1 . The mRNA expression profiles in the GSE142255 and GSE156382 datasets were obtained from the GEO database. Samples in the GSE142255 dataset are whole blood samples, including 7 healthy samples and 17 ACLF samples (6); The samples in the GSE156382 dataset are neutrophil samples from blood, including 5 healthy samples and 12 ACLF samples (13); Through the published literature, 137 NRGs were identified (14). Both iron death and copper death-related genes were obtained from a comprehensive literature review, among which 268 ferroptosis-related genes were obtained (17). Additionally, ten cuproptosis-related genes were identified, namely LIAS, FDX1, CDKN2A, DLD, DLAT, LIPT1, GLS, MTF1,PDHB and PDHA1 (18).

We analyzed the differentially expressed genes between ACLF and the healthy control (HC) samples in the GSE142255 dataset using the R software limma package. To identify novel diagnostic biomarkers associated with ACLF, a total of 288 DEGs were screened. DEGs that had log fold changes greater than 1 and p-values below 0.05 after adjusting for false discovery rate (FDR) were selected. Among these DEGs, 150 genes were found to be up-regulated ( Supplementary Table S3 ), while 138 genes were identified to be down-regulated ( Supplementary Table S4 ). The distribution of these 288 DEGs is illustrated in the accompanying volcano plots and heat maps and 10 top upregulated and downregulated genes were selected for display ( Figures 2A, B ).

To identify key modules, the R package WGCNA was utilized to conduct weighted correlation network analysis on the samples and data within the training set. The network module analysis is susceptible to the influence of outlier samples. We first constructed a hierarchical clustering tree for the sample (N=24), as shown in ( Supplementary Figure 1A ). The results showed the presence of outlier samples, and after removing two outlier samples ( Supplementary Figure 1B ), there were no obvious outlier samples. Therefore, the filtered samples were used for subsequent analysis. A scale-independent topological network was established, characterized by a soft-thresholding power of nine and a scale-free R² value of 0.9, along with a mean connectivity network ( Figures 3A, B ). We segmented the dendrogram at pertinent transition points using dynamic hybrid cutting to construct a hierarchical clustering tree. The leaves of the trees represent individual genes, while the branches of the dendrogram represent clusters of genes with similar expression profiles. Branches containing genes with analogous expression patterns were designated as gene modules. In total, six gene modules were produced through the amalgamation of analogous modules ( Figure 3C ). Among these six modules, the gray module is meaningless. As shown in Figure 3D , the blue module had the highest and significant correlation with disease phenotype (Cor=-0.8261517 p value=2.17e-06). Therefore, the 6332 genes included in this module would be marked as ModuleGenes and used for subsequent analysis. The correlation between gene modules and blue modules was shown in Figure 3E .

To obtain differentially expressed NRGs related to ACLF (ACLF-NRGs), an intersection of DEGs, neutrophil-related genes (NRGs), and module genes (ModelGenes) associated to ACLF was identified. As shown in Figure 3F , a total of 14 intersection genes (ACLF-NRGs) were obtained. Details of the 14 ACLF-NRGs were shown in the Supplementary Materials ( Supplementary Table S5 ).

14 intersecting genes (ACLF-NRGs) were mapped onto human chromosomes by using the RCircos package in the R programming language ( Figure 4A ). As illustrated in Figure 4B , Gene Ontology (GO) analysis identified and examined 13 cellular components (CC), 118 biological processes (BP), and 6 molecular functions (MF) associated with these intersecting genes. The top five enrichment terminologies pertaining to BP, CC, and MF are presented in Supplementary Table S6 . The three most significantly enriched biological process (BP) terms were related to the immune response: neutrophil degranulation, regulation of neutrophil activation, and neutrophil activation involved in the immune response. Correspondingly, the top three highly enriched cellular component (CC) terms identified through Gene Ontology (GO) analysis were specific granule, tertiary granule membrane, and tertiary granule lumen. Concurrently, the most enriched molecular function (MF) terms were protein serine/threonine/tyrosine kinase activity, calcium-dependent protein binding, and protein kinase activator activity. Comprehensive details of the Gene Ontology (GO) analysis reports can be found in the Supplementary Materials ( Supplementary Table S7 ).

A Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis of the intersecting genes identified 15 immune-associated pathways, encompassing the TNF signaling pathway, Staphylococcus aureus infection, prolactin signaling pathway, epithelial cell signaling in Helicobacter pylori infection, and leukocyte transendothelial migration ( Figure 4C ). Comprehensive details of the KEGG analysis are provided in the Supplementary Materials ( Supplementary Table S8 ).

The analysis of Protein-protein interaction (PPI) networks for ACLF-NRGs was conducted using the Cytoscape STRING database, resulting in the construction of a PPI network. LogFC of each gene was obtained from differential analysis to enrich the network graph information ( Figure 4D ).

A Random Forest (RF) analysis was conducted on 14 intersecting genes utilizing the Random Forest package. The optimal number of classification trees ( Supplementary Figure 2A ) and the gene importance ranking ( Supplementary Figure 2B ) were obtained. The selection criteria for this analysis were to include genes with mean decrease accuracy that is not zero. The scores of all 14 genes were not zero, so they were all included in the subsequent analysis.

To eliminate the recursive features of the obtained differential genes, the SVM algorithm in R package e1071was employed. RF method was used to sort and calculate the above feature genes, obtaining their significance and ranking, along with the accuracy of each combination of iterations and error rate. The combination with the lowest error rate was selected as the optimal set, and the corresponding genes were designated as feature genes. As shown in the figure below, the highest accurate subset (genes included in the blue dots) of GSE142255 contained four feature genes: FCAR, MMP9, S100A12, CCL5 ( Supplementary Figure 2C ).

The intersection of genes obtained by RF and SVM-RFE algorithms was used to determine feature genes, as shown in the Figure 5A . A total of four biomarkers were obtained, namely FCAR, MMP9, S100A12, and CCL5 ( Figure 5A ).

Subsequent Lasso analysis on these four genes identified three diagnostic genes: MMP9, S100A12, and CCL5. Using the training set data GSE142255, these three diagnostic genes were employed to construct diagnostic models for the ACLF and Control groups, formulated as follows: CCL5 * 1.86 + MMP9 * (-3.21) + S100A12 * (-1.06). This model was then employed to assess the diagnostic capability of the identified genes for ACLF. The Receiver Operating Characteristic (ROC) analysis was conducted on this model, yielding an Area Under the Curve (AUC) value of 1.0, which signifies perfect predictive accuracy and suggests its potential utility in disease diagnosis ( Figure 5B ).

A nomogram was subsequently developed using the diagnostic genes, facilitated by the RMS package ( Figure 5C ). The differential expression levels of MMP9, S100A12, and CCL5 between ACLF and HC are presented in Figure 5D , with significant p-values less than 0.0001. Among the three diagnostic genes, MMP9 and S100A12 were up-regulated, while CCL5 was down-regulated in ACLF. Furthermore, in GSE168049 dataset, compared to the survival group, S100A12 and MMP9 were up-regulated while CCL5 was down-regulated in the death group, indicating the three diagnostic genes were potential prognostic factors ( Figure 5E ).

The lasso diagnostic model constructed by three diagnostic genes was validated using datasets GSE156382, referred to as verifying dataset. In the verifying dataset, the results of the ROC analysis are shown in Figure 6A . The model achieved an AUC value greater than 0.9, demonstrating high diagnostic accuracy. Furthermore, the expression levels of the three diagnostic genes were consistent with those in the training set ( Figure 6B ).

In order to further understand the role of diagnostic genes play in the molecular mechanisms of ACLF, IPA was employed to identify the signaling pathways significantly influenced by these genes. The top 15 pathways involved in diagnostic genes were displayed ( Figure 7A ), including four disease-related pathways: pathogen induced cytokine storm, glucocorticoid receptor signaling, natural killer cell signaling and P38 MAPK pathway (19–21). The biological function heatmap enriched by diagnostic genes is presented in Figure 7B . Moreover, the upstream regulatory factors and downstream target genes of the three diagnostic genes are presented in Figure 7C .

“c2. cp. kegg. v11.0. symbols” were downloaded as the internal reference gene set from MSigDB database. Gene Set Enrichment Analysis (GSEA) were performed on three specific genes: MMP9, S100A12, and CCL5. The KEGG pathways with the top three positive and negative Normalized Enrichment Scores (NES) were displayed ( Supplementary Figure 3 ).

Livers from ACLF patients exhibited key features of ferroptosis, implicating its role in the pathogenesis of ACLF (22). To investigate the correlation between Ferroptosis or cuprotosis and the three diagnostic genes, Pearson correlation analysis was performed. Through published literature, ferroptosis-related genes (FRGs) and cuprotosis-related genes (Cuprotosis) were obtained. Single-sample Gene Set Enrichment Analysis (ssGSEA) was employed to evaluate the ferroptosis ssGSEA score (FRGs Score) and the cuprotosis ssGSEA score (Cuprotosis Score) for the ferroptosis-related and cuprotosis-related gene sets in the GSE142255 dataset. Among the three diagnostic genes analyzed, CCL5 exhibited a significant negative correlation with the FRGs Score, whereas MMP9 and SA1002 demonstrated significant positive correlations with the same score ( Figure 8A ). Conversely, the correlation patterns between these three diagnostic genes and the Cuprotosis Score were found to be opposite to those observed with the FRGs Score ( Figure 8B ). Additionally, three diagnostic genes demonstrated a stronger association with the FRGs Score compared to the Cuprotosis Score.

qRT-PCR, longside WB or ELISA was performed on peripheral white blood cells or plasma samples from HC, CHB and ACLF to validate the expression levels of the three diagnostic genes. A novel diagnostic criterion for HB-ACLF was established, known as the Chinese Group on the Study of Severe Hepatitis B-ACLF (COSSH-ACLF) (23). According to this criterion, patients with CHB exhibiting a total bilirubin (TB) level of ≥12 mg/dL and an international normalized ratio (INR) of ≥1.5 should be diagnosed with ACLF, irrespective of the presence of cirrhosis. This criterion includes patients with higher TB levels than those defined by APASL. Since all the patients included were HB-ACLF except for one case, the qRT-PCR ACLF cohort was further stratified into two groups: patients with TB ≥12 mg/dL were categorized under the COSSH-ACLF group, while those with TB <12 mg/dL were classified under the APASL-ACLF group. The mRNA expression levels of the three diagnostic genes were significantly elevated in the COSSH-ACLF group compared to the HC, CHB and APASL-ACLF groups ( Figure 9A ). Additionally, there were an increase in S100A12 and CCL5 expression in COSSH-ACLF than AD. Among the three key genes, only the expression of S100A12 were upregulated in AD, compared with HC and CHB. No statistically significant differences were observed among the HC, CHB and APASL-ACLF groups ( Figure 9A ). However, at the protein level, only MMP9 was confirmed to be more highly expressed in the ACLF group compared to the CHB group ( Figures 9B, C ). No significant differences were found between the CHB and ACLF groups in the protein expression of S100A12 and CCL5 ( Figures 9B–D ).

Among the three key genes, MMP9 expression exhibited the most significant difference between the death and survival groups ( Figure 5E ). Consequently, additional ELISA analyses were conducted to assess MMP9 levels. Specifically, MMP9 ELISA was conducted on plasma samples from 15 HC, 11 patients with CHB, 8 patients with HB-AD, and 40 patients with HB-ACLF to validate the differential expression levels of MMP9. The results indicated a significant elevation of MMP9 levels in the ACLF group compared to the other three groups ( Figure 10A ). MMP9, also known as type IV collagenase, is implicated in fibrosis. There was no significant difference in plasma MMP9 levels between patients with cirrhosis and those without ( Figure 10B ). It was interesting to note that plasma MMP9 levels in AD patients were significantly lower than those in HC and CHB patients. CHB and HC did not differ significantly in MMP9 levels ( Figure 10A ). Further ELISA analysis on additional cases of ACLF revealed that plasma MMP9 levels were significantly higher in patients who succumbed to the condition compared to those who survived ( Figure 10C ). ACLF patients with elevated MMP9 levels (>175.8 ng/ml) exhibited higher short-term mortality rates within both 30 and 90 days ( Figure 10D , p<0.001).

Drug Bank database was used to predict small molecule drugs for diagnostic genes. A total of 21 drugs were obtained targeting the three diagnostic genes. Specifically, CCL5 and S100A12 each yielded two predictive results, while MMP9 yielded 17 predictive results ( Figure 11A ). Among the 17 predictive drugs identified for targeting MMP9, the following were approved: Minocycline, Captopril, Glutathione, Zinc, Zinc acetate, Zinc chloride, and Zinc sulfate. Additionally, Amlexanox and Olopatadine were the two approved predictive drugs targeting S100A12. No approved drugs were identified for CCL5; however, two experimental drugs, Heparin Disaccharide I-S and Heparin Disaccharide III-S, were noted. Detailed information regarding these drugs was provided in the Supplementary Materials ( Supplementary Table S9 ). The approved drugs would be discussed below.

Transcription factors (TF) that regulate diagnostic genes were retrieved from NetworkAnalyst database (https://www.networkanalyst.ca/). Then the upstream kinases of TF were predicted by X2Kgui and miRNA of diagnostic genes were predicted by using miRDB database. At last, The Kinase-TF-mRNA-miRNA network was constructed by Cytoscape software ( Figure 11B ).

The prediction of these drugs aims to screen for drugs that may affect the pathological process of ACLF by regulating the expression or function of these genes. Through drug screening and regulatory network analysis, we provided preliminary insights into potential treatment strategies and targets. Moreover, The Kinase-TF-mRNA-miRNA network also provided the potential regulatory mechanisms of these genes, establishing a basis for mechanistic research on these key genes in ACLF.