A comprehensive summary of the patient cohort enrolled in our study is presented in Supplementary Table 1 . With patient age spanning from 23 to 71 years, the broad age distribution underscores the extensive impact of Hp infection across the entire population. The diverse age range vividly reflects the ubiquitous nature of Hp’s influence. Notably, Hp isolates derived from these patients underwent rigorous drug testing to ascertain their susceptibility to six antibiotics, namely Amoxicillin, Clarithromycin, Levofloxacin, Metronidazole, Tetracycline, and Furazolidone. Interestingly, 86.5% of strains (45/52) are resistant to Metronidazole.

The statistic for enrolled patients and analyzed strains are summarized in Supplementary Table 2 . The patient cohort included 21/52 (40.4%) male patients and 31/52 (59.6%) female patients. The median age is 50. The rate of resistance are 17.3% for Amoxicillin, 44.2% for Clarithromycin, 38.5% for Levofloxacin and 86.5% for Metronidazole. All strains are sensitive to Furazolidone. Only one strain is resistant to Tetracycline.

To investigate how the genotype information of Hp can be used to predict antibiotic resistance, we established a workflow to measure drug sensitivity phenotype using disk diffusion test and genotype using whole genome sequencing ( Figure 1 ). The relation between genotype and phenotype was modeled using univariate statistical test, unsupervised machine learning (phylogenetic analysis) and supervised machine learning (support vector machine).

The pivotal aspect of our investigation involved subjecting Hp isolates obtained from the patients to cutting-edge next-generation sequencing. This technological approach facilitated the comprehensive acquisition of whole genome data. To ensure the reliability and quality of our findings, patients lacking accurate clinical data were thoughtfully excluded, as were patients yielding insufficient SNP data due to data quality concerns. Seven patients were excluded because lack of clinical information. Five patients were removed from downstream analysis as the total number of SNP is less than 3000. This stringent curation led to the inclusion of a robust cohort comprising 52 patients. In total, our study unveiled a median of 66,198 SNPs identified among these patients. The significant SNPs for Amoxicillin, Clarithromycin, Levofloxacin, Metronidazole and Tetracycline were summarized in ( Supplementary Tables 3–7 ).

Leveraging this pool of extracted SNPs, we conducted a series of rigorous downstream analyses. These analyses encompassed single SNP statistical tests, the construction of phylogenetic trees, and the implementation of supervised machine learning using Support Vector Machine. This multi-faceted approach allowed us to draw comprehensive insights into the genetic intricacies underpinning antibiotic resistance in Hp strains within our patient cohort.

Three modeling approaches were undertaken to address the relationship between antibiotic resistance and genotype, including univariate statistical test, unsupervised machine learning and supervised machine learning. Univariate statistical test serves as a foundational method to identify basic relationships and differences between resistant and sensitive strains. This approach is valuable for its simplicity and clarity, providing an initial understanding of the data. For the unsupervised machine learning component, phylogenetic analysis was employed. Phylogenetic trees facilitates the visualization and clustering of different strains, offering insights into how genetic variations might influence antibiotic resistance. For supervised machine learning, the Support Vector Machine (SVM) algorithm was selected for its effectiveness in classification tasks, especially in high-dimensional spaces. This approach is particularly beneficial when labeled data is available and predictive modeling is desired.

To evaluate the model performance, we used sensitivity and specificity throughout our study. Sensitivity and specificity are two key metrics used to evaluate the performance of a diagnostic test or model. Briefly, in our study sensitivity measures the proportion of actual positives (resistant strains) that are correctly identified by the model. Specificity measures the proportion of actual negatives (sensitive strains) that are correctly identified by the model.

We identified 2483 significant SNPs for Amoxicillin ( Supplementary Table 3 ), including one significant site on 23S rRNA (A1834G). Using individual SNP for Amoxicillin resistance prediction, we evaluated the model with sensitivity and specificity. In this case the presence of a SNP suggested that the sequenced strain was an Amoxicillin resistant one. This assumption allowed us to determine the sensitivity and specificity of the individual models based on one single SNP site. We plotted the sensitivity and specificity of all resulted models ( Figure 2A ). The sensitivity ranged from 0 to 1, with a mean sensitivity of 0.45. The specificity also ranged from 0 to 1, with a mean specificity of 0.46. It was difficult to obtain a good tradeoff between model sensitivity and specificity ( Figure 2B ). For A1834G of 23S rRNA, the model achieved 100% specificity and only 22.2% sensitivity.

To enhance model performance, we took all 2483 significant sites as input for the construction of phylogenetic tree. The phylogenetic tree indicated that sensitive and resistant strains were well separated ( Figure 2C ).

In our investigation, we identified a total of 266 significant SNPs associated with Clarithromycin resistance ( Supplementary Table 4 ), among which a noteworthy site was found on the 23S rRNA gene (A2147G). To predict Clarithromycin resistance, we utilized individual SNPs as predictive markers and assessed the model’s performance. The model’s effectiveness was evaluated by plotting the sensitivity and specificity of the resulting models ( Figure 3A ). Sensitivity values ranged from 0 to 1, with a mean sensitivity of 0.4, while specificity values ranged from 0 to 1, with a mean specificity of 0.58. Notably, except for A2147G, achieving a balanced tradeoff between model sensitivity and specificity proved to be challenging ( Figure 3B ). However, for the specific case of A2147G on the 23S rRNA gene, the model demonstrated 87% sensitivity and 96% specificity.

To further enhance the predictive capabilities of the models, we incorporated all 266 significant SNP sites as input for the construction of a phylogenetic tree. The outcome of this phylogenetic analysis indicated that strains with different sensitivities to Clarithromycin could be distinctly separated ( Figure 3C ). This result highlighted that through phylogenetic analysis, different strains of Helicobacter pylori were distinctly categorized based on their individual reactions to the antibiotic Clarithromycin, as evidenced by their separation in the phylogenetic tree.

We identified a total of 200 significant SNPs associated with Levofloxacin resistance ( Supplementary Table 5 ), among which a significant site was found on the 23S rRNA gene (A2147G). To predict Levofloxacin resistance, we employed individual SNPs as predictive markers and evaluated the model performance. The performance metrics of the resultant models were assessed by plotting sensitivity and specificity ( Figure 4A ). Sensitivity values ranged from 0 to 1, with a mean of 0.4, while specificity values ranged from 0 to 1, with a mean of 0.58. Interestingly, there appeared to be a negative correlation between sensitivity and specificity ( Figure 4B ). Notably, the model for A2147G of the 23S rRNA gene achieved 65% sensitivity and 71% specificity.

To maximize the separation in resistant and sensitive strains, we utilized all 200 significant SNP sites as input to construct a phylogenetic tree. The resulting phylogenetic tree exhibited distinct separation between Levofloxacin-sensitive and resistant strains ( Figure 4C ). This phylogenetic analysis provided additional support for the differentiation between bacterial strains based on their Levofloxacin resistance profiles.

We successfully identified a total of 379 significant SNPs associated with Metronidazole resistance ( Supplementary Table 6 ), among which an important site on the 23S rRNA (G1517A) stood out. Employing individual SNPs for Metronidazole resistance prediction, we evaluated the performance metrics of our models. The sensitivity and specificity of all resulting models were plotted and analyzed ( Figure 5A ). The spectrum of sensitivity scores ranged between 0 and 1, with an average sensitivity value of 0.41. Correspondingly, specificity scores ranged between 0 and 1, yielding an average specificity of 0.52. Remarkably, one single SNP model emerged as an outlier, achieving nearly 100% sensitivity and an impressive 100% specificity ( Figure 5B ). However, when focusing on the known G1517A variant of the 23S rRNA, the model exhibited only 2% sensitivity and 50% specificity.

For unsupervised clustering of all strains based on SNP information, we leveraged all 379 significant sites as input for the construction of a phylogenetic tree. The resulting tree compellingly demonstrated that Metronidazole-sensitive and resistant strains could be distinctly separated ( Figure 5C ).

Furthermore, we conducted a comprehensive investigation into the potential impact of two different gastric mucosal protective agents (Group A and B) on treatment outcomes for Hp infections. More information about treatment scheme is documented in method. Hp is eradicated in 83.3% of strains treated with A and 85.7% of strains treated with B (p = 1). For the eight strains not successfully eradicated by Amoxicillin and Clarithromycin, four strains are resistant to both Amoxicillin and Clarithromycin and the other four stains are resistant to either Amoxicillin (n =2) or Clarithromycin (n =2).

Using the known A1834G of 23S rRNA, the univariate model achieved 100% specificity and only 22.2% sensitivity. The mean sensitivity and specificity for single SNP model are only 0.45 and 0.46. All single SNP models with 100% specificity have sensitivity less than 50%. Therefore, we aimed to investigate the application of the supervised machine learning algorithm SVM in developing a predictive model for Amoxicillin resistance. Here we employed all significantly different SNP sites between resistant and sensitive strains to train SVM model. The feature selection procedure reduced the number of dimensions for the classification problem. Through this approach, we successfully generated an SVM model with 100% sensitivity and 100% specificity in predicting sensitive and resistant Hp strains ( Figure 6 ). For the four strains (two lowly sensitive strains and two intermediate strains regarding Amoxicillin response) not used in model training, the model predicted all four strains as Amoxicillin sensitive.

Although this model allows us to obtain a good tradeoff between sensitivity and specificity, it is subjected to potential bias as we have much more variables than observations. To maximize the generalizability of the model, we implemented two validation procedures. In the first procedure, the original dataset was randomly split into training set (60%) and test set (40%). We obtained a model with 66.6% sensitivity and 100% specificity in the test set. In the second procedure, we employed 5-fold cross-validation. All five resulted models displayed 100% specificity, while the maximal sensitivity obtained is 50% in three models.

Using the known A2147G on the 23S rRNA gene, the model to predict Clarithromycin resistance demonstrated 87% sensitivity and 96% specificity. We further asked whether we could improve the model performance if we incorporate all significantly different SNP between Clarithromycin sensitive and resistant strains. Consequently, our investigation delved into the application of a supervised machine learning algorithm, specifically Support Vector Machine (SVM), in the development of a predictive model for Clarithromycin resistance. Remarkably, our efforts yielded a model that achieved 100% sensitivity and 100% specificity in predicting Clarithromycin-sensitive and resistant Hp strains ( Figure 7 ).

Similar to the case of Amoxicillin resistance prediction with SVM, the resulted model strikes a perfect balance of sensitivity and specificity that is subjected to overfitting. The validation was implemented by two procedures. In the split procedure (60% training and 40% test), the model reached 100% sensitivity and 100% specificity in the test data. Using 5-fold cross-validation, all five resulted models exhibited 100% sensitivity and 100% specificity.

In total we collected sample from 52 patients with Hp infections. All enrolled patients have no intake of antibiotics within one month, no consumption of PPI or Chinese traditional medicine within two weeks. UBT test was employed to confirm the positivity of Hp infection. Biopsy samples were preserved in cultivation media and transported to laboratory at 4 degree. The studies involving human participants were reviewed and approved by the seventh Medical Center of PLA General Hospital Ethics Committee (No.2017-74) and other departments of gastroenterology from different hospitals applied and followed this content’s introduction. The patients/participants provided their written informed consent to participate in this study.

Patients were divided into two groups (A and B). The patients in the group A were treated with 1000 mg tid of hydrotalcite, 20 mg bid of rabeprazole, 1000 mg bid of amoxicillin and 500 mg bid of clarithromycin; The patients in group B were treated with colloidal bismuth pectin 300 mg bid + rabeprazole 20 mg bid + amoxicillin 1000 mg bid + clarithromycin 500 mg bid for 10 days. At least 28 days after the end of treatment, all patients received ¹³C urea breath test to evaluate and compare the Hp eradication rate between two groups.

Biopsy was inoculated onto Columbia blood agar plates supplemented with designated antibiotics. Inoculation was performed by direct contact of mucosa side and agar. The contact is gentle and even to ensure successful inoculation. Preservation media was applied to evenly cover the whole plate. The biopsy was removed and used in to confirm the presence of Hp. Agar plates were incubated in 5% O₂, 10% CO₂, 85% N₂ at 37 degree for 3 or 4 days.

Hp colony was validated with morphology (Gram negative) and biochemical assays. Single Hp colony was transferred to a new Columbia agar plate and inoculated evenly. After 3-4 days of culture in incubator with 5% O₂, 10% CO₂, 85% N₂ at 37 degree, Hp colonies were cryo-preserved in -80 degree for next generation sequencing.

Determination of antibiotic sensitivity was performed with disk diffusion test according to standardized Kirby-Bauer procedure. White paper disks containing antibiotics were arranged onto the agar plates. The size of circular zones of poor bacteria growth surrounding paper disks was measured and used to define antibiotic sensitivity using customized thresholds ( Supplementary Figure 1 ).

DNA extraction was performed with Magen HiPure Bacterial DNA kits. Final DNA output was dissolved in 30 μl TE buffer and analyzed with Qubit (Thermofisher) to determine concentration.

10 ng DNA was used for library construction. DNA was fragmented with sonification to obtain fragments of 300 bp. Beads-based selection was performed to select DNA fragments after end repairing and adapter ligation. Library amplification was performed with Super Canace High Fidelity enzyme. Amplified libraries and size selected and purified. Agilent 2100 was used for library QC. The resulted libraries were normalized by molar concentration and sequenced with paired-end 150 bp mode on illumina platform.

Our mutation calling process commenced with rigorous quality control measures aimed at eliminating adapters and low-quality reads, ensuring the generation of clean reads. These clean reads were then meticulously aligned to the Helicobacter pylori reference genome, specifically the Helicobacter pylori 26695 genome (accessible at https://www.ncbi.nlm.nih.gov/nuccore/AE000511) using BWA (v0.7.17) with parameters “mem -t 2 -M -Y”. The resulting alignment data, in BAM format, was sorted and indexed with samtools (v1.9), which was the basis for mutation calling through the use of bcftools (v1.9) with parameters “call -vmO z”. To enhance the accuracy of our findings, only SNP sites demonstrating at least 2x coverage across all samples were retained for subsequent analysis. We performed “bcftools filter -O v -s LOWQUAL -e ‘\’’QUAL<10 || FMT/DP <5’\’’ –SnpGap 5 –set-GTs” to filter out SNPs of low quality. The resulting SNPs were annotated with Annovar (v191024).

The determination of significant SNPs was achieved via rigorous statistical assessments. Our analysis incorporated two independent statistical tests, namely the Chi-squared test and Fisher’s exact test. Employing a stringent approach, only SNPs exhibiting p-values below 0.05 in both tests were considered significant and subsequently utilized in our downstream analytical endeavors.

To gauge the efficacy of our approach for each antibiotic, we formulated confusion matrices specific to each significant SNP. The calculation of Sensitivity involved determining the ratio between true positives and model-predicted positives, while Specificity was calculated as the ratio between true negatives and model-predicted negatives. The confusion matrix is comprised of true positive (TP), false positive (FP), false negative (FN) and true negative (TN). To be precise, the formulars are “sensitivity = TP/(TP+FN)” and “specificity = TN/(TN+FP)”. These calculations were executed in the R programming environment.

Phylogenetic trees were constructed using SNP data identified in sample subsets of varying sizes, encompassing 5, 10, 20, 30, 40, and 50 samples, yielding 56671, 22601, 18749, 13938, 10397, and 3470 SNPs, respectively. Furthermore, we also plotted phylogenetic trees using only the significant SNPs identified through the aforementioned statistical tests, providing a comprehensive visual representation of genetic variation. The brief step for phylogenetic analysis is as follows. The SNPs in VCF file were converted to FASTA format and imported to MEGA (v10.2), then maximum likelihood tree was constructed with default parameters. The resulting tree file (nwk fomat) was then imported to an online tool called iTOL (https://itol.embl.de) for visualization.

To train classifier to predict antibiotic sensitivity using SNP information, input data was prepared using significant SNP for antibiotic under investigation. To facilitate the SVM modeling process, SNP data was systematically converted into numeric values, serving as the input for our SVM model. Briefly, input data was a patient-SNP matrix consisting of 0 or 1. SNP value of “1” suggests a polymorphism (non-wild type status), while SNP value of “0” indicates same genotype as the reference strain.

For Amoxicillin antibiotic resistance, we had 39 sensitive strains, 9 resistant strains, 2 intermediate strains and 2 lowly-sensitive strains. For Clarithromycin, we had 27 sensitive strains, 23 resistant strains and 2 intermediate strains. As the label “intermediate” and “lowly-sensitive” had very few observations, we aimed to build the two class SVM classifier using only sensitive and resistant strains. As for Clarithromycine, we had 27 sensitive strains, 23 resistant strains and 2 intermediate strains. Only sensitive and resistant strains were used for SVM.

The implementation of Support Vector Machine (SVM) modeling was facilitated through the e1017 package in the R programming environment. Employing a “radial” kernel, with a cost set at “1” for optimal SVM performance, we ensured model reproducibility by utilizing set.seed before the SVM modeling process. Sensitivity and specificity metrics were calculated using the above mentioned formulars: “sensitivity = TP/(TP+FN)” and “specificity = TN/(TN+FP)”.

In the model validation process for the Support Vector Machine (SVM) analysis, two distinct validation procedures were implemented to ensure the robustness and reliability of the model. The first procedure involved partitioning the original dataset into two subsets: 60% of the data was used as the training set, where the SVM model was trained to understand the patterns and relationships in the data. The remaining 40% constituted the test set, which was utilized to evaluate the model’s performance, specifically its ability to accurately predict outcomes on new, unseen data. This approach of splitting the dataset provides a straightforward way to check the model’s efficacy and generalizability.

The second procedure employed was 5-fold cross-validation, a more rigorous validation technique. In this method, the dataset was divided into five equal parts. In each of the five iterations of the process, a different fold was used as the validation set, while the remaining four folds collectively served as the training set. This cycle ensured that each part of the dataset was used both for training and validation. Cross-validation is particularly valuable as it mitigates the risk of overfitting, ensuring the model’s performance is not overly tailored to a specific subset of data. By averaging the results from all five folds, a more comprehensive and reliable assessment of the model’s performance is obtained, enhancing confidence in its predictive accuracy. Together, these two validation procedures provide a thorough examination of the SVM model’s capabilities, contributing to its credibility and applicability in practical scenarios.