The retrospective pharmacovigilance study collected data from the FAERS Quarterly Data Extraction website between the first quarter of 2004 and the fourth quarter of 2023. The FAERS datasets contained demographic information, drug information, and reaction information. Each report contained at least one suspected drug and one or more ADE terms. The downloaded data were imported into SQL server software, and we established our own convenient and directly usable database. Data cleaning procedures were as follows: when the CASEID (the number used to identify FAERS cases) and FDA_DT (the date when the FDA received the case) were the same, the data unit with the latest FDA_DT was kept in the DEMO table for de-duplication, and then the “deleted” cases were removed. The entire process of de-duplication was based on recommendations by the FDA.

The ADE terms in the FAERS database were based on the International Medical Terminology Dictionary for Regulatory Activities (MedDRA). The DIA-related preferred terms (PTs) in the MedDRA structure were investigated according to the standardized MedDRA query (SMQ) narrow search “agranulocytosis.” The resultant 15 PTs and PT_code are listed in Supplementary Table S1.

The generic names, brand names, synonymous names, or abbreviations contribute to the various drug names in the FAERS database. We used MedEx software (MedEx UIMA 1.3.8, Vanderbilt University, US) to convert different forms of drug names into the “generic name” for standardization. The potential target drugs were counted on the basis of the “generic name,” and the primary suspected (PS) drugs were included for further analysis.

As agranulocytosis cases were identified, drugs with case number >10 were included for further analysis. Frequency difference between specific drug event combination (DEC) and the background was investigated following the disproportionality principle. It was considered a significant signal when the DEC frequency was higher than the background frequency and within the set threshold, which revealed a statistical correlation between the target drug and the target ADE.

The reporting odds ratio (ROR), a frequency-based method, was one of the most commonly used disproportionality algorithms. A greater ROR value indicates a stronger signal, signifying a robust statistical relationship between the target drug and ADE. An ROR signal was identified when at least three cases were reported with the lower limit of 95% confidence interval (CI) of ROR >1. In addition to the ROR approach, this study incorporated the information component (IC) method, a Bayesian method, to detect signals of disproportionality. An IC signal was detected when the lower limit of 95% CI > 0. Both the methods were regularly utilized for the identification of disproportionality signals, with the IC method demonstrating superior sensitivity and specificity, particularly in datasets with smaller sample sizes. In this study, the DIA signal was identified through the generation of both the ROR and IC signals.

FAERS reports containing patient information (gender, age, and weight) were extracted, and only those with complete data were analyzed. Age >120 years and weight >400 kg were defined as outliers and excluded from the analysis.

A univariable logistic analysis for suspicious drugs was conducted using the following parameters: 95% CI lower limit of ROR > l, case number >100, and p-value <0.01. As a result, drugs with p <0.01 were further used for least absolute shrinkage and selection operator (LASSO) regression (Zhang et al., 2023). Multivariate logistic regression was able to determine the existence of DIA risk factors using the drugs screened by cross-validation LASSO (CV-LASSO) combined with basic patient information as independent variables. After constructing a multivariate logistic regression model, we used the receiver operating characteristic (ROC) curve to calculate the area under the curve (AUC) value for assessing the accuracy of the model’s predictions. When 0.5< AUC <1, the larger AUC value was regarded as a better classifier of the model.

The medications would be systematically classified for further analysis according to the multivariate logistic regression results.

The onset time was defined as the interval between EVENT_DT (date of adverse event occurrence) and START_DT (start date for drug use). Moreover, the input errors (EVENT_DT earlier than START_DT) reports or inaccurate date entries were excluded. Therefore, the drug-induced time was the difference between the two aforementioned time points and depicted as the cumulative curve and violin plot.

R software (version 4.2.1) was used for logistic regression and time-to-onset analysis.

We systematically searched PubMed, Embase, The Cochrane Library, and Web of Science from the inception to April 2025. For pharmacovigilance, case reports serve as a critical source for identifying rare and serious ADE. Therefore, SRs of randomized clinical trials (RCTs), observational studies, and case reports were included. The whole search strategy was provided in Supplementary Table S2. In brief, we utilized the controlled term “agranulocytosis” combined with free terms such as “neutropenia.” No language limits were applied in the search strategy.

Regarding the Population, Intervention, Comparison, Outcome, and Study design (PICOS) framework, there were no limitations in population, intervention, and comparison. Consistent with the FAERS study, we included SRs that clearly investigated the association between any intervention and agranulocytosis. The excluded criteria were as follows: not full-text report, grade-1/2 neutropenia, grade not reported, congenital neutropenia, neutropenia developed before therapy, and reported outcome not identified as severe neutropenia or agranulocytosis.

Two authors (YL and KL) screened each record independently with discrepancies, and conflicts were settled by the third author (BW). Relevant data such as the author, year, time range, number of studies included, number of patients with DIA, interventions, and outcomes were extracted and listed in a Microsoft Excel document by YL and meticulously reviewed by KL.

The methodological quality of the included SRs was assessed in accordance with the Assessment of Multiple Systematic Reviews 2 (AMSTAR-2) (Shea et al., 2017). Following the guidance document of AMSTAR-2, the 16 items in the questionnaire were completed by YL and KL independently, and then the overall quality of each study was generated based on critical weaknesses. Meanwhile, any disagreement regarding the AMSTAR-2 scale was resolved by BW to reach consensus.

To obtain safety information comprehensively, SRs of all research types were included, so quantitative synthesis was not applicable in this study. As part of qualitative evaluation, the synthesis was structured thematically around drug classes, aligning with the classification framework used in the disproportionality analysis to facilitate cross-validation. We systematically extracted and synthesized key findings from the included SRs, including the cumulative incidence, effect size of relative risk (RR) with their CI, heterogeneity, and discrepancies.

A total of 160,487 agranulocytosis cases with PS drugs were identified, accounting for 0.92% of all ADE cases in FAERS. The PS drugs of each report were consolidated, resulting in 1,853 PS drugs. Then, 598 PS drugs with more than 10 cases were included and analyzed (Figure 1A). The annual numbers of reports from 2004 to 2023 are shown in Figure 2.

The average age of the reported cases was 54.79 ± 21.69 years, and the other characteristics of the patients, including basic information such as gender, outcome, the reporter’s profession, and reporter’s region, are shown in Supplementary Table S3. The number of female patients (71,422, 44.50%) was higher than that of male patients (68,183, 42.49%). A total of 137,500 cases (85.68%) were reported by health professionals. The top region of reporting was Europe (58,051, 36.17%). There were 19,679 (12.26%) death cases and 10,384 (6.47%) life-threatening cases.

A total of 251 drugs with more than ten cases were identified with DIA signals. We classified the positive DIA signal drugs according to the WHO Anatomical Therapeutic Chemical (ATC) system (Supplementary Table S4). The top-five drug categories warranting attention were monoclonal antibodies and antibody–drug conjugates (ATC code: L01F, 23,656 cases, 35 drugs), immunosuppressants (ATC code: L04A, 21,099 cases, 21 drugs), other antineoplastic agents (ATC code: L01X, 17,345 cases, 33 drugs), antimetabolites (ATC code: L01B, 14,869 cases, 13 drugs), and plant alkaloids and other natural products (ATC code: L01C, 9,230 cases, 11 drugs), as shown in Supplementary Table S5.

According to the WHO ATC classification, we reported the top-50 drugs by the case number (Table 1) and the top-50 drugs by the signal strength (Table 2) in the main text. The top-five drugs by case number were methotrexate (6,462 cases), lenalidomide (5,722 cases), rituximab (5,691 cases), doxorubicin (4,391 cases), and carboplatin (4371 cases). The top-five drugs by signal strength were toripalimab (ROR = 73.93; 95% CI 45.26–120.77), mesna (ROR = 46.16; 95% CI 38.52–55.32), cytarabine (ROR = 43.32; 95% CI 41.55–45.17), dexrazoxane (ROR = 42.55; 95% CI 33.50–54.05), and clofarabine (ROR = 38.63; 95% CI 35.23–42.36).

Notably, since we defined a positive signal when both the ROR and IC signals were generated, we excluded 52 drugs. These drugs had positive ROR signals but negative IC signals. For example, tislelizumab (ROR = 4.45, 95% CI 2.43–8.13; IC = 2.11, 95% CI -0.17–3.67), diazoxide (ROR = 4.26, 95% CI 2.39–7.58; IC = 2.05, 95% CI -0.11–3.57), norfloxacin (ROR = 3.58, 95% CI 1.96–6.53; IC = 1.81, 95% CI -0.38–3.44), and other drugs could not confirm the association with agranulocytosis.

To systematically evaluate the drug–event associations, we initially performed univariable logistic regression on 27 pharmacologically active compounds identified through univariate screening. Subsequently, LASSO regression with 10-fold CV-LASSO was applied to these candidate drugs and potential confounding factors (age, weight, and gender). The optimal penalty parameter λ selected 26 significantly associated drugs, along with age and weight as meaningful covariates (Figure 3). As shown in Figure 4, a multivariable logistic regression model incorporating these selected variables demonstrated statistically significant associations for all included parameters (p < 0.01). Deferiprone, thiamazole, cyclophosphamide, and clozapine are strongly associated with DIA. The accuracy of the model was validated by ROC curve analysis, yielding an AUC of 0.761, indicating satisfactory predictive accuracy (Figure 5).

Based on the multivariate logistic regression results, the medications were systematically classified into seven subgroups for comparative analysis of the temporal patterns in adverse reaction occurrence.

Then, we performed a statistical analysis of the association between the duration of drug use and DIA. The results in Supplementary Table S6 and Figure 6 reveal significant differences in the duration of use. Clozapine, an antipsychotic drug, had the longest onset time (mean: 1,121.3, Q1: 51.0, Q3: 1,331.0). In contrast, antibiotic drugs such as amoxicillin and azithromycin had relatively shorter usage durations, averaging 11.0 days and 8.1 days, respectively, with lower associations with agranulocytosis. Chemotherapeutic agents, such as carboplatin, had an average duration of 18.0 days, with a certain level of association with agranulocytosis.

A total of 742 records were identified, which included 261 from Embase, 183 from PubMed, 21 from Cochrane Library, and 277 from Web of Science. After 377 duplicates were removed, 365 records were left to be screened. Of these, we excluded 346 records at title/abstract screening, eight at full-text, and four at grade 1/2 or not reported, resulting in seven studies (Andersohn F. et al., 2007; Lyman et al., 2016; Kassem et al., 2018; Omar et al., 2020; Holz J.M. et al., 2022; Jiang X. et al., 2022; Moinard-Butot F. et al., 2023) being included for further review (Figure 1B). The records excluded along with the specific reasons are listed in Supplementary Table S7.

Table 3 summarizes the key characteristics of included SRs examining DIA events, which spanned a broad time range from 1966 to 2023. The included SRs investigated a range of interventions, encompassing tyrosine kinase inhibitors, CDK4/6 inhibitors, immune checkpoint inhibitors (ICIs), chemotherapeutic regimens, antibiotics, and other agents such as antithyroid drugs and psychotropics. The majority of the outcomes described in these studies was agranulocytosis, with or without severe complications. Among these SRs, three included only RCTs, another three included solely case reports, and just one incorporated both RCTs and observational studies. Sample sizes varied greatly from 83 to 8,593 patients, which is largely attributed to the difference in the types of included research.

The quality of the included SRs assessed by the AMSTAR-2 scale was moderate in one, low in four, and critically low in two. Four SRs failed to satisfy the critical domains that pertain to the consideration of risk-of-bias (RoB) in individual studies during the interpretation and discussion of review results (Supplementary Table S8).

Five studies focused on cancer treatment, three of which investigated protein kinase inhibitors (PCIs). As the first-generation anaplastic lymphoma kinase (ALK) inhibitor, crizotinib-induced grade 3/4 neutropenia was reported in 0%–30% of patients with lung cancer, while the rate of the second-/third-generation was reported as 1.1% and <1% (Moinard-Butot F. et al., 2023), respectively. The summary incidences of grade 3/4 neutropenia induced by sunitinib in patients with gastrointestinal stromal tumors was 9.3% (95% CI, 5.6–13.7%) (Jiang X. et al., 2022). In addition, the RR of sunitinib was reported as 10.39 (95% CI, 1.53–70.72; p = 0.017) compared to the placebo, with heterogeneity not being observed (p < 0.01). Another SR reported that the incidence of CDK 4/6 inhibitor-induced grade 3/4 neutropenia in breast cancer patients ranged from 21.1% to 66% with the calculated RR at 44.00 (95% CI 24.72–78.33; p < 0.0001) (Kassem et al., 2018). Subgroup analyses showed no significant differences in heterogeneity attributed to the dosage of sunitinib and the type of CDK 4/6 inhibitors (palbociclib vs ribociclib vs abemaciclib).

Agranulocytosis and its complications were also identified as an ADE caused by ICIs through SRs of case reports, but they were less commonly reported than all grades of neutropenia (Omar et al., 2020).

Regarding chemotherapy-induced grade 3/4 neutropenia, the median event rate varied significantly with regimen intensity, ranging from single agents at 5.0% to ACVBP-like at 18.0% (Lyman et al., 2016). Notably, the addition of rituximab did not significantly alter the risk characteristics of chemotherapy-related myelosuppression.

In a systematic review of 62 case reports involving 83 patients, antibiotic-induced neutropenia predominantly presented as grade 3/4 severity, accounting for 85% of cases (Holz J.M. et al., 2022). Among these cases, vancomycin and ceftaroline were reported as the most abundant drugs. The lowest ANC caused by non-penicillin β-lactams was significantly lower than that of other classes (p = 0.008, t-statistic = 9.604), with a further significant difference noted between β-lactams and glycopeptides (p = 0.006, t-statistic = 3.099).

In another systematic review of 672 case reports, antibiotic-induced agranulocytosis demonstrated the highest proportion among all the included DIA cases, with 113 cases (16.8%) determined as definite/probable casual relationships (Andersohn F. et al., 2007). This study also showed that antithyroid drugs, psychotropic drugs, and antiarrhythmics were commonly linked to agranulocytosis, with 87, 85, and 39 cases reported for each of them, respectively.