Between January 2022 and December 2023, a total of 1,048 ART-treated individuals with virologic failure were identified in Aksu Prefecture based on routine clinical monitoring. Virologic failure was defined as a confirmed HIV-1 diagnosis, at least 6 months of ART, and an HIV RNA level >1,000 copies/mL. Among these individuals, HIV-1 pol sequences were successfully obtained from 704 samples (67.2%), while sequencing failure in the remaining 344 samples was mainly due to low viral RNA levels, RNA degradation, or PCR amplification failure. If an individual experienced more than one episode of virologic failure during the study period, only the earliest eligible sample was included. All subsequent samples from the same individual were excluded to avoid duplicate representation. After removing duplicate samples, 675 ART-treated individuals with valid pol sequences constituted the final analytic population.

Sociodemographic and clinical information was collected for all participants. Sociodemographic variables included sex, ethnicity, age, marital status, educational attainment, and mode of HIV transmission. Clinical data included duration of ART, current ART regimen, and drug resistance results, which were extracted from participants’ medical records.

Five milliliters of EDTA-anticoagulated whole blood was collected from each participant. One tube of whole blood was processed within 6 h for CD4⁺ T-lymphocyte testing at the local Centers for Disease Control and Prevention (CDCs), with analysis completed within 24 h. The remaining whole-blood samples were centrifuged at 3,000 rpm for 15 min, after which plasma was separated and stored at −80 °C. All plasma samples were subsequently transported under frozen conditions to the Aksu Regional Center for Disease Control and Prevention for HIV drug resistance testing.

Viral RNA was extracted from plasma samples using the QIAamp Viral RNA Mini Kit (Qiagen, Germany). A nested PCR approach was used to amplify a 1.3-kb fragment of the HIV-1 pol region (corresponding to positions 2,166–3,462 of the HXB2 reference strain), according to a previously described protocol (Wang et al., 2018). PCR products were visualized by electrophoresis on a 1% agarose gel, purified, and sequenced by Beijing Dehongchangyuan Biotechnology Co., Ltd. Sequence fragments were assembled using Sequencher version 5.4.6.

All sequences were aligned with HIV-1 reference sequences retrieved from the Los Alamos National Laboratory HIV Sequence Database using MEGA version 11.0.13 for multiple sequence alignment and preliminary phylogenetic analysis. HIV-1 subtypes were determined using both the HIV BLAST tool and the China HIV Gene Sequence Data Management and Analysis System¹ (Feng et al., 2023). Drug resistance mutations were identified using the Stanford HIV Drug Resistance Database (HIVdb) algorithm (version 9.8, accessed in March 2025). All sequences were screened for APOBEC-mediated hypermutation using the Hypermut 2.0 tool and were further evaluated for sequence quality. No sequences were excluded based on these quality control procedures.

Resistance levels were assessed for a total of 20 antiretroviral drugs, including seven nucleoside reverse transcriptase inhibitors (NRTIs), five non-nucleoside reverse transcriptase inhibitors (NNRTIs), and eight protease inhibitors (PIs). The evaluated drugs comprised doravirine (DOR), efavirenz (EFV), nevirapine (NVP), etravirine (ETR), rilpivirine (RPV), fosamprenavir (FPV), abacavir (ABC), zidovudine (AZT), stavudine (D4T), didanosine (DDI), emtricitabine (FTC), lamivudine (3TC), tenofovir (TDF), atazanavir/ritonavir (ATV/r), darunavir/ritonavir (DRV/r), fosamprenavir/ritonavir (FPV/r), indinavir/ritonavir (IDV/r), nelfinavir (NFV), lopinavir/ritonavir (LPV/r), saquinavir/ritonavir (SQV/r), and tipranavir/ritonavir (TPV/r). Drug resistance was categorized into five levels based on the Stanford HIVdb scoring system: susceptible (0–9), potential low-level resistance (10–14), low-level resistance (15–29), intermediate-level resistance (30–59), and high-level resistance (≥60). Acquired drug resistance (ADR) was defined as the presence of at least one antiretroviral drug with a resistance score≥15, excluding cases classified solely as potential low-level resistance. ART-treated individuals with DRMs were defined as those harboring any mutation with a Stanford HIVdb resistance score≥10.

An approximate maximum-likelihood (ML) phylogenetic tree was reconstructed using FastTree version 3.0, with reference sequences representing major HIV-1 subtypes included as outgroups. The generalized time-reversible model with gamma-distributed rate variation and a proportion of invariant sites (GTR + G + I) was applied for nucleotide substitution. Node support was evaluated using the Shimodaira-Hasegawa (SH) test implemented in the software. The resulting phylogenetic tree was further analyzed using Cluster Picker version 1.2.5 to identify genetic clusters, defined as clades with bootstrap support ≥95% and pairwise genetic distance ≤3%. All 675 sequences were included in the phylogenetic reconstruction. The final tree was visualized using iTOL version 6.

Pairwise genetic distances between sequences were calculated using HyPhy version 2.2.4 under the TN93 model, and the minimum genetic distance approach was applied to infer pairwise relationships. The TN93 genetic distance threshold for network construction was selected within the commonly accepted range for recent HIV transmission (0.5–1.5%), balancing sensitivity to ongoing transmission with control of spurious genetic linkages. The final threshold of 1.41% was applied for subsequent network analysis of DRMs clustering. Genetic networks were visualized using R version 4.4.2.

Network centrality metrics were defined in accordance with established methodologies (de Arruda et al., 2014; Little et al., 2014). Degree centrality was defined as the number of direct connections incident on a node, with nodes having a degree ≥3 classified as highly connected. Betweenness centrality quantified the extent to which a node lay on the shortest paths between other nodes, reflecting its role as a network bridge. Closeness centrality was calculated as the reciprocal of the sum of the shortest path distances from a node to all other reachable nodes, indicating the efficiency of information flow within the network. For both betweenness and closeness centrality, values >0.75 (normalized to a 0–1 scale) were considered high. Nodes exhibiting high values across all three centrality measures were designated as “central nodes” within the network (Zeng et al., 2021).

All statistical analyses were performed using R software. Categorical variables across different antiretroviral therapy regimens were compared using the chi-squared test to assess differences in HIV drug resistance mutation profiles. Factors associated with DRMs clustering were evaluated using univariable and multivariable logistic regression analyses. Variables with p < 0.10 in the univariable analyses, together with those of established clinical relevance, were entered into the initial multivariable model. A stepwise selection procedure was applied to derive the final multivariable model. Potential multicollinearity among covariates was assessed using variance inflation factors (VIFs), with all VIF values <5 indicating the absence of significant multicollinearity. Model fit was evaluated using the Hosmer–Lemeshow goodness-of-fit test. All statistical tests were two-tailed, and a p value < 0.05 was considered statistically significant.

A total of 675 ART-treated individuals with virological failure were included in this study, of whom 398 (58.96%) were male and 277 (41.04%) were female. Most participants were aged 15–49 years, accounting for 77.19% (521/675) of the study population. The predominant route of HIV transmission was heterosexual contact, which accounted for 92.00% (621/675) of cases. Overall, 50.52% (341/675) of participants were married, and 94.37% (637/675) were of Uyghur ethnicity. Approximately 75.56% (510/675) of patients were receiving first-line ART, mainly TDF/AZT + 3TC + EFV/NVP–based regimens. More than half of the participants (55.11%, 372/675) had a viral load >10,000 copies/mL at the time of sampling (Supplementary Table S1). The median duration of ART was 58 months (interquartile range [IQR]: 38–83 months). CRF07_BC was the predominant HIV-1 genotype, accounting for 97.63% (659/675) of cases, whereas CRF01_AE was rare (0.74%, 5/675) (Supplementary Table S1 and Figure 1).

Among the 675 ART-treated individuals with virological failure, 407 (60.2%) were found to harbor drug resistance mutations (DRMs). The overall prevalence of ADR was 56.9% (384/675) during 2022–2023 (Figure 2a). Resistance to NRTIs, NNRTIs, and PIs was observed in 20.0% (135/675), 51.3% (346/675), and 7.4% (50/675) of participants, respectively. Notably, dual-class resistance to both NRTIs and NNRTIs was observed in 17.8% (120/675) of participants (Figure 2a). A total of seven NRTI-associated, five NNRTI-associated, and eight PI-associated DRMs were identified. The most common NRTI-associated mutation was M184V/I (27.3%), followed by S68G/N (7.4%) and K70R/E (5.7%). The most prevalent NNRTI-associated mutation was K103N/S (60.7%), followed by G190A/E/S (11.3%), E138A/K/Q/G (10.8%), and P225H (10.3%). All PI-associated mutations were observed at low frequencies, except for Q58E (10.07%) (Table 1 and Figure 2b).

Drug susceptibility was evaluated for 20 antiretroviral agents among sequences harboring any drug resistance mutations (DRMs), including NRTIs, NNRTIs, and PIs. High-level resistance to EFV and NVP was observed in 72.5 and 75.7% of sequences, respectively. Within the NRTI class, high-level resistance to FTC and 3TC was more prevalent than that to other NRTIs (Figure 2c). The primary NRTI resistance mutation M184V/I conferred high-level resistance to 3TC and FTC, whereas K65R was associated with intermediate-level resistance to abacavir (ABC) and tenofovir disoproxil fumarate (TDF). The most frequent NNRTI-associated resistance mutations included K103N/S, V106M/V, and G190A/E/S, which conferred high-level resistance to both NVP and EFV (Figure 2d).

DRMs were detected across all ART regimens examined. Among the 12 major NRTI-associated mutations, 11 were identified in both the TDF + 3TC + EFV/NVP and LPV/r + 3TC + TDF/AZT regimens, and 10 were observed in the AZT + 3TC + EFV/NVP regimen, whereas only three were detected in other regimens. The K65R mutation, which is associated with reduced susceptibility to TDF, was significantly more frequent among individuals receiving the TDF + 3TC + EFV/NVP regimen (p < 0.05). The K103N/S mutation showed the highest overall prevalence across all regimens. PI-associated DRMs were generally uncommon; however, Q58E was more frequently observed in the LPV/r + 3TC + TDF/AZT regimen, despite not typically conferring high-level resistance to LPV/r.

All 407 pol sequences harboring DRMs were included in the genetic network analysis using a TN93 genetic distance threshold of 1.41%. In total, 139 sequences (34.2%) formed 45 clusters, including 32 dyads and 13 clusters containing three or more sequences (Figure 3). After excluding sequences with only potential low-level resistance (score 10–14), 30 clusters containing clinically relevant DRMs (score ≥15) remained. Most clusters involved NNRTI-associated DRMs, with K103N/S, E138A/K/Q/G, V179D/E/F, and G190A/E/S being the most frequently observed. Among NRTI-associated DRMs, M184V/I and S68G/N predominated, whereas Q58E and L10F were the major PI-associated DRMs.

Centrality analysis was restricted to clusters containing at least three sequences. Within these clusters, 14.4% of nodes exhibited a degree centrality ≥3. A total of 18 central nodes were distributed across six genetic clusters, of which 10 harbored dual-class resistance mutations. Most central nodes (16/18) had viral loads >10,000 copies/mL (Figure 3a). Among NNRTI-associated DRMs, K103N/S was the most prevalent mutation (35 clusters), followed by E138A/K/Q/G and V179D/E/F (each in eight clusters), and G190A/E/S (seven clusters). Five central nodes harbored compound mutations, including K103N/S, P225H, and K238T/N (Figure 3b). M184V/I was the predominant NRTI-associated DRM (13 clusters), followed by S68G/N (seven clusters). Two major mutation clusters associated with central nodes involved L210W and M184V/I (Figure 3c). Among PI-associated DRMs, Q58E and L10F were detected within the network. L10F was observed only among individuals receiving first-line regimens, whereas Q58E was identified in patients failing both first- and second-line ART. This pattern suggests the possible presence of Q58E in clustered viruses prior to treatment initiation rather than acquisition solely under selective drug pressure (Figure 3d).

To identify factors associated with DRMs clustering, 14 independent variables were assessed using logistic regression analysis (Supplementary Table S2). Univariable analysis showed that ethnicity, viral load level, and current ART regimen were significantly associated with membership in a DRMs cluster (p < 0.05). In contrast, sex, age, marital status, educational level, transmission route, latest CD4⁺ T-cell count, duration of ART, presence of NRTI-, NNRTI-, or PI-associated DRMs, and viral subtype were not significantly associated with clustering (all p > 0.05; Supplementary Table S2). In the multivariable analysis, individuals with viral loads of 10,000–50,000 copies/mL (OR = 1.856, 95% CI:1.119–3.080) and>50,000 copies/mL (OR = 1.723, 95% CI: 1.035–2.870) had significantly higher odds of belonging to a DRMs cluster. In contrast, participants receiving the LPV/r + 3TC + AZT/TDF regimen were significantly less likely to belong to DRMs clusters (Figure 4).