Discovery and In Silico Characterization of Anatolian Water Buffalo Rumen-Derived Bacterial Thermostable Xylanases: A Sequence-Based Metagenomic Approach

This study involved shotgun sequencing of rumen metagenomes from three Anatolian water buffalos, an exploration of the relationship between microbial flora and xylanases, and in silico analyses of thermostable xylanases, focusing on their sequence, structure, and dynamic properties. For this purpose, the rumen metagenome of three Anatolian water buffalos was sequenced and bioinformatically analyzed to determine microbial diversity and full-length xylanases. Analyses of BLAST, biophysicochemical characteristics, phylogenetic tree, and multiple sequence alignment were performed with Blastp, ProtParam, MEGA11 software, and Clustal Omega, respectively. Three-dimensional homology models of three xylanases (AWBRMetXyn5, AWBRMetXyn10, and AWBRMetXyn19) were constructed by SWISS-MODEL and validated by ProSA, ProCheck, and Verify3D. Also, their 3D models were structurally analyzed by PyMOL, BANΔIT, thermostability predictor, What If, and Protein Interaction Calculator (PIC) software. Protein–ligand interactions were examined by docking and MD simulation. Shotgun sequence and Blastp analyses showed that Clostridium (Clostridiales bacterial order), Ruminococcus (Oscillospiraceae bacterial family), Prevotella (Bacteroidales bacterial order), and Butyrivibrio (Lachnospiraceae bacterial family) were found as dominant potential xylanase-producer genera in three rumen samples. Furthermore, the biophysicochemical analysis indicated that three xylanases exhibited an aliphatic index above 80, an instability index below 40, and melting temperatures (T_m) surpassing 65 °C. Phylogenetic analysis placed three xylanases within the GH10 family, clustering them with thermophilic xylanases, while homology modeling identified the optimal template as a xylanase from a thermophilic bacterium. The structural analysis indicated that three xylanases possessed the number of salt bridges, hydrophobic interactions, and T_m score higher than 50, 165, and 70 °C, respectively; however, the reference thermophilic XynAS9 had 43, 145, and 54.41 °C, respectively. BANΔIT analysis revealed that three xylanases exhibited lower B′-factor values in the β3-α1 loop/short-helix at the N-terminal site compared to the reference thermophilic XynAS9. In contrast, six residues (G79, M123, D150, T199, A329, and G377) possessed higher B′-factor values in AWBRMetXyn5 and their aligned positions in AWBRMetXyn10 and AWBRMetXyn19, relative to XynAS9 including Gln, Glu, Ile, Lys, Ser, and Val at these positions, respectively. MD simulation results showed that the β9-η5 loop including catalytic nucleophile glutamic acid in the RMSF plot of three xylanases had a higher fluctuation than the aligned region in XynAS9. The distance analysis from the MD simulation showed that the nucleophile residue in AWBRMetXyn5 and AWBRMetXyn10 remained closer to the ligand throughout the simulation compared with XynAS9 and AWBRMetXyn19. The most notable difference between AWBRMetXyn5 and AWBRMetXyn10 was the increased amino acid fluctuations in two specific regions, the η3 short-helix and the η3-α3 loop, despite a minimal sequence difference of only 1.24%, which included three key amino acid variations (N345, N396, and T397 in AWBRMetXyn5; D345, K396, and A397 in AWBRMetXyn10). Thus, this study provided computational insights into xylanase function and thermostability, which could inform future protein engineering efforts. Additionally, three xylanases, especially AWBRMetXyn5, are promising candidates for various high-temperature industrial applications. In a forthcoming study, three xylanases will be experimentally characterized and considered for potential industrial applications. In addition, the amino acid substitutions (G79Q, M123E, D150I, T199K, A329S, and G377V) and the residues in the β3-α1 loop will be targeted for thermostability improvement of AWBRMetXyn5. The amino acids (N345, N396, and T397) and the residues on the β9-η5 loop, η3 short-helix, and η3-α3 loop will also be focused on development of the catalytic efficiency.