Site-selective dynamics of ligand-free and ligand-bound azidolysozyme

Proteins are essential for function and sustaining the life of organisms. Experimentation and computational studies have clarified that protein function involves both structure and dynamics.^1–3 However, characterizing structural and functional dynamics of proteins at the same time under physiological conditions in the condensed phase, which is prerequisite for understanding cellular processes at a molecular level, remains challenging.¹ Vibrational spectroscopy, in particular, two-dimensional infrared (2D IR) spectroscopy, has been shown to be a powerful tool for studying the structural dynamics of various biological systems.⁴ One of the particular challenges is to obtain structural and environmental information in a site-specific manner. To address this, significant effort has been focused on the development and application of various infrared (IR) reporters^5,6 that absorb in the frequency range of 1700–2800 cm⁻¹ to discriminate the signal from the strong protein background.^7,8 Such IR probes have provided valuable information about the structure and dynamics of complex systems. For example, nitrile probes have helped to clarify the role of electrostatic fields in enzymatic reactions^9,10 or to elucidate the mode of drug binding to proteins.^11,12 Isotope edited carbonyl spectroscopy was used to characterize the mechanism of protein folding and amyloid formation^13,14 or the structure and function of membrane proteins.^15,16 Additional molecular groups, such as thiocyanate,¹⁷ cyanamide,¹⁸ sulfhydryl vibrations of cysteines,¹⁹ deuterated carbons,²⁰ carbonyl vibrations of metal-carbonyls,^21–23 cyanophenylalanine,²⁴ and azidohomoalanine (Aha),²⁵ have also been explored.

In this work, AlaN₃, an analog of azidohomoalanine that has been shown to sensitively report on local structural changes while still being minimally invasive,^26,27 is used as the probe. This modification can be incorporated into proteins through established expression techniques.²⁸ The asymmetric stretch frequency of –N₃ is at $∼2100$ cm⁻¹ and has a comparatively high extinction coefficient of 300–400 M⁻¹ cm⁻¹, which makes it a suitable spectroscopic reporter.²⁵ Aha has been used for biomolecular recognition after incorporation into the peptide directly^25,26 or in the vicinity of the binding area of a PDZ2 domain²⁹ to detect the water-specific response of azide vibrations when attached to small organic molecules³⁰ or to probe the frequency shift and fluctuation due to its sensitivity to the local electrostatic environments and dynamics.^27,31 Such studies confirm that AlaN₃ and/or Aha are environment-sensitive IR probes and suitable modifications for site-specific investigations of protein dynamics. It should be noted that actual experiments are carried out in D₂O rather than in H₂O to minimize overlap with the bend–libration combination band of H₂O that absorbs in the vicinity of the azide asymmetric stretch vibration.

With its picosecond time resolution, IR spectroscopy provides direct information about the structural dynamics around a probe molecule with high temporal resolution.^4,32 Moreover, introducing IR probes with isolated vibrational frequencies overcomes the problem of spectral congestion that complicates discrimination and analysis of desired vibrational bands. With that, the inter- and intramolecular coupling between degrees of freedom or the local structure or dynamics of biological systems can be specifically probed and characterized. Such an approach relies on the sensitivity of the probe to report on changes in the vibrational frequencies induced by alterations in the local electrostatic interactions in the vicinity of the probe.²⁴ 

IR spectroscopy is a potentially advantageous technique to characterize ligand binding to proteins.^33,34 Its success depends, in part, on the notion that when a ligand binds to a protein, the frequency of an infrared active vibration shifts due to the different electric field in solution—often water—and in the protein binding site due to the Stark effect. Such an approach often requires the ligand to be modified, e.g., through the addition of a suitable label, such as –CN as in benzonitrile. The relationship between shifting the IR response of the label and the influence on the binding affinity has been explicitly assessed through computer simulations for benzonitrile in the active site of wild type (WT) and mutant lysozyme.³⁴ 

Alternatively, the protein can be selectively modified by attaching a spectroscopic label at strategic positions so that the binding process and functional dynamics can be interrogated with the functional and unmodified ligand. This has the potential advantage that interactions between the ligand and the surrounding protein are unaltered. These interactions contribute the majority of the enthalpic part to the binding free energy and, therefore, directly affect the affinity of the ligand and its rate of unbinding. The dynamics of WT lysozyme without and with labeled alanine (AlaN₃) has been recently found to provide position-specific information about the spectroscopy and dynamics of the modification site.²⁷ The structure of the protein with the labeled Ala residues is shown in Fig. 1 together with the binding site lined by residues Leu84, Val87, Leu91, Leu99, Met102, Val111, Ala112, Phe114, Ser117, Leu118, Leu121, Leu133, and Phe153. Following Ref. 34, cyano-benzene (PhCN) was used as the ligand to allow direct comparison. The WT structure was used here (a) to compare directly with earlier results²⁷ and (b) because PhCN has a comparatively small binding free energy toward the WT protein (ΔG_bind = −0.5 kcal/mol), which suggests that the interaction between the ligand and the protein is weak.³⁴ For the L99A mutant protein, ΔG_bind = −3.9 kcal/mol for PhCN,³⁴ which compares with an experimentally determined value of $∼−3.5$ kcal/mol for iodobenzene from isothermal titration calorimetry.³⁵ 

In this work, changes in 1D- and 2D-IR signatures of the azido group attached to all alanine (Ala) residues of lysozyme upon binding of cyano-benzene (PhCN) are determined. In addition, the changes of the environmental dynamics around all AlaN₃ are quantified for ligand-free vs ligand-bound lysozyme. Such differences are experimentally observable and yield valuable insight into the energetics and dynamics of protein–ligand binding. This paper is structured as follows: First, the methods are presented. This is followed by results and discussion of the structural dynamics and the spectroscopy and frequency fluctuation correlation functions (FFCFs). Finally, conclusions are drawn.

MD simulations for the WT and all azide-labeled (AlaN₃) modified lysozymes were carried out using an adapted version of the Chemistry at Harvard Molecular Mechanics (CHARMM) program³⁶ with an interface to perform simulations with the reproducing kernel Hilbert space (RKHS) potential energy surface (PES).³⁷ The initial structure is that of T4 lysozyme with benzene in the cavity (PDB 1L83).³⁸ Benzene was replaced by cyano-benzene (PhCN) maintaining carbon atom positions.³⁴ The protein was solvated in explicit TIP3P water³⁹ using a cubic box of size (78)³ ų. First, all systems were minimized, followed by heating and equilibration. Next, 2 ns NVE production simulations were carried out with and without the PhCN ligand present in the active site for all 16 protein variants with Ala replaced by AlaN₃. Bonds involving H-atoms were constrained using the SHAKE⁴⁰ algorithm, and all nonbonded interactions were evaluated with shifted interactions using a cutoff of 14 Å and switched at 10 Å.⁴¹ Snapshots for analysis were recorded every 5 fs.

The total energy function for the system consisted of the CHARMM force field⁴² for the protein, the TIP3P³⁹ parameterization for water, and a reproducing kernel (RKHS)^43,44 based representation for the spectroscopic probe (–N₃).²⁷ For the –N₃ label, a full-dimensional, accurate potential energy surface (PES) calculated at the pair natural orbital based coupled cluster [PNO-LCCSD(T)-F12/aVTZ]^45,46 level is available. The RKHS-representation exactly reproduces the reference energies and provides analytical derivatives, which allows us to carry out high-accuracy, energy-conserving MD simulations.^37,44 This is of particular relevance for applications to vibrational spectroscopy. Details of the representation and validation of the full-dimensional PES are provided in Ref. 27. For the PhCN ligand, parameterization is based on Swissparam,⁴⁷ while charges are fitted to the electrostatic potential MP2/aug-cc-pVDZ level of theory.³⁴ 

The effect of ligand binding on the overall flexibility of the modified protein can be assessed from considering the root mean squared fluctuation (RMSF) of the C_α atoms. Depending on the position at which the –N₃ label is located, the changes in RMSF range from insignificant (Ala82N₃ or Ala160N₃) to major (Ala73N₃ or Ala112N₃) (see Figs. 2 and S1). Changes in flexibility can either be local (Ala82N₃ or Ala160N₃) or global (Ala73N₃ or Ala112N₃) similar to the RMSF variation. In addition, the average RMSF for the ligand-bound protein can be larger or smaller than that for the ligand-free protein. This suggests that insertion of the ligand can either rigidify the protein or increase its flexibility.

It is also of interest to compare all ligand-free –N₃-labeled RMSFs (blue traces in Figs. 2 and S1), which provides information about the influence of attaching the spectroscopic label on the protein dynamics. As an example, flexibility at position Asp20 is, in general, high. The same applies to the region Ile100 to Met120. In general, attaching –N₃ to alanine residues along lysozyme leads to differential changes in the flexibility depending on the location of the modification site. Hence, although the azide is considered a “minimally invasive” probe, depending on the modification site the local or even global dynamics can still be affected.

For a more global assessment of the changes in dynamics and potential couplings between different parts of the protein induced by attaching the –N₃ label and/or introducing the PhCN ligand, dynamical cross-correlation maps^48,49 (DCCMs) were calculated from the trajectories using the Bio3D package.⁵⁰ Dynamic cross-correlation matrices are based on the expression

where r_i and r_j are the C_α atom positions of the respective ith and jth amino acids and Δr_i corresponds to the displacement of the ith C_α atom from its trajectory-averaged position. DCCMs report on the correlated and anticorrelated motions within a protein, and difference DCCMs (ΔDCCM) provide a global view of the positionally resolved differences in the dynamics. In the following, only absolute values for C_ij and differences between them that are larger than 0.5 are reported. The DCCMs are symmetrical about the diagonal, and for clarity, positive correlations (for DCCM) or positive differences in C_ij (for ΔDCCM) are displayed in the lower right triangle and negative values or differences in C_ij are displayed in the upper left triangle.

The DCCM for lysozyme with Ala129N₃ ligand-free, ligand-bound, and the difference between the two is shown in Fig. 3. These maps reveal ligand-induced differences in the correlated and anticorrelated motions with appreciable amplitudes [see features A–D in Fig. 3(c)]. For ligand-free lysozyme, there are pronounced couplings between residues [130, 147] and [20, 25]/[32, 37] for anticorrelated motions and residues [68, 80] and [103, 112] for correlated motions. As demonstrated in Fig. 3(b), upon binding the ligand to lysozyme, the DCCM shows different coupled residues compared to the ligand-free protein. As an example, residues [35, 45] and [55, 68] are affected more for anticorrelated motions, whereas for correlated motions the coupling is between residues [5, 15] and [55, 65]. Note that these effects may not be visible in the ΔDCCM as the magnitude of the difference between the two systems may be below the threshold of 0.5 in the ΔC_ij.

In the difference map [Fig. 3(c)], feature A indicates coupling between residues [135, 145] and [20, 25]/[30, 42], whereas feature B refers to coupled residues [65, 75] and [58, 65]. Furthermore, feature C demonstrates prominent variations between residues [84, 95] and [117, 125], whereas for feature D residues [129, 140] and [140, 147] are strongly correlated. These findings suggest that residues couple both locally (feature B/D) and through space (feature A/C). It should also be pointed out that residues involved in features A–C are among those with higher RMSF (see Fig. S1). Interestingly, the region around residue Ala146 with larger differences ΔC_ij displays correlated dynamics with spatially close residues around residue Asp20 (white licorice in Fig. 1). On the other hand, the pronounced differences in the RMSF of Ala129N₃ (see Fig. S1) for residues [42, 57] do not show up in ΔDCCM because their C_ij coefficients are below the threshold of 0.5.

Difference DCCMs (ΔDCCM) between WT and ligand-bound (WT-PhCN) lysozyme are shown in Fig. 4(a) together with the difference map between WT and azido modified lysozyme at position 134 (Ala134N₃) [see Fig. 4(b)]. With the PhCN ligand bound to the protein, the ΔDCCM compared with that for the ligand-free protein is sparsely populated [see Fig. 4(a)]. This indicates that the conformational dynamics of the two systems is similar. Contrary to that, a larger number of differences in the dynamics between WT and Ala134N₃ arise, as shown in Fig. 4(b). Finally, the ΔDCCM between the ligand-free (Ala134N₃) and ligand-bound (Ala134N₃-PhCN) labeled lysozyme at Ala134 shown in Fig. 4(c) demonstrates that the “contrast” further increases. The major difference in the conformational dynamics between the ligand-free and ligand-bound protein arises for coupled residues [57, 65] with [105, 120] (feature A), [103, 112] with [112, 122] (feature B), [60, 68] with [47, 52]/[55, 62] (feature C), and [57, 65] with [17, 25]/[32, 42] (feature D). Interestingly, as mentioned before for residue 129, residues involved in features A and B are also among those with higher RMSF (see Fig. S1). However, the DCCM and ΔDCCM provide considerably more information about changes in the average dynamics. Additional difference DCCMs (ΔDCCM) for the remaining AlaN₃ residues are shown in Figs. S2–S15.

Correlated motions in hen egg white lysozyme (HEWL) labeled with ruthenium dicarbonyl attached to His15 were also investigated from a normal mode analysis of the C_α atoms and based on an elastic network model.²¹ This analysis also reported both correlated and anticorrelated motions for low-frequency normal modes as is found in the present case. Furthermore, such analyses help us to more quantitatively characterize the slaving motion of the protein to the surrounding solvent.²¹ 

Using instantaneous normal mode (INM) analysis,^27,37,51 the frequency trajectory ω(t) of the asymmetric stretch vibration of the –N₃ label was determined. Based on this, the 1D infrared spectra corresponding to the azide asymmetric stretch vibration for each of the 16 AlaN₃ residues was computed for the ligand-free and ligand-bound protein [see Figs. 5(a), 5(b), and S16]. A direct comparison of the maximum position of the infrared line shape in Fig. 5(c) shows that for three N₃-modified alanine residues (Ala41, Ala98, and Ala130), the difference in the absorption frequency is insignificant. For positions Ala63, Ala73, and Ala160, the differences are 8, 3, and 2 cm⁻¹, respectively, whereas for the other residues the change is within 1 cm⁻¹.

Such frequency changes can be measured with state-of-the art experiments,²⁴ and their magnitude is also consistent with previous simulations of the vibrational Stark effect for the –CN probe in PhCN with red shifts of up to 3.5 cm⁻¹ in going from the WT to the L99A and L99G mutants of T4-lysozyme.³⁴ Similarly, the 1D and 2D infrared spectroscopy of –CO as the label for the insulin monomer and dimer found⁵² that the relative shifts of the spectroscopic response was correctly described, whereas the absolute frequencies may differ by some 10 cm⁻¹. In a recent study, such an approach found a splitting of 13 cm⁻¹, compared with 25 cm⁻¹ from experiment, for the outer and central –CO labels in cationic trialanine in water.⁵³ Hence, MD simulations together with instantaneous normal modes are a meaningful approach to determine relative frequency shifts, whereas capturing absolute frequencies in such simulations requires slight reparameterization of the underlying force field, e.g., through morphing techniques.^54,55

The magnitude of frequency shifts found from the present simulations is also comparable with −3 cm⁻¹ reported from experiments of the nitrile stretch in ligand IDD743 bound to WT vs V47N mutant hALR2³³ or a +6 cm⁻¹ blue shift of the –CO vibrational frequency due to the binding of 19-NT to the Asp40Asn mutant of the protein ketosteroid isomerase compared to the WT.⁵⁶ Thus, differences of $∼1$ cm⁻¹ for the frequency of the reporter in different chemical environments can be experimentally detected.²⁴ 

From the frequency trajectories, the frequency fluctuation correlation function (FFCF) can be determined, which contains valuable information on relaxation time scales corresponding to the solvent dynamics around the solute. The FFCFs are fit to an empirical expression

which allows analytical integration to obtain the line shape function⁵⁷ using an automated curve fitting tool from the SciPy library.⁵⁸ As was found for the RMSFs and 1D IR spectra, the FFCFs from the simulations with and without the ligand bound to the protein can be very similar or differ appreciably [see Fig. 6(a)]. The slow decay time, τ₂, of the –N₃ asymmetric stretch mode of the label is typically shorter for the ligand-bound protein compared to that without PhCN [see Fig. 6(b)], although exceptions exist. For Ala97N₃, Ala112N₃, and Ala134N₃, the slow relaxation time τ₂ is faster by 75% up to a factor of $∼2.5$⁠, and for Ala146N₃, the slow time scale, τ₂, differs by a factor of $∼3$ between the ligand-free (τ₂ = 5.13 ps) and ligand-bound (τ₂ = 1.61 ps) lysozyme. For the other alanine residues, the τ₂ times between ligand-free and ligand-bound lysozyme are similar. As an exception, for Ala129N₃, the decay is slowed down by $∼50$% for PhCN-bound lysozyme. All FFCFs without (Fig. S17) and with (Fig. S18) the ligand are given in the supplementary material and the fitting parameters are reported in Table I.

As a last feature of the FFCF, it is found that the static component Δ₀ can differ appreciably between ligand-free and ligand-bound lysozyme [see Fig. 6(c)]. The static offset Δ₀ is an experimental observable and characterizes the structural heterogeneity around the modification site. There are only four alanine residues for which the static offset is similar (Ala41, Ala49, Ala82, and Ala93) for ligand-bound and ligand-free lysozyme. For all others, the differences range from 15% to a factor of $∼3$⁠. As an example, for Ala73N₃ the difference for $Δ02$ between bound and ligand-free lysozyme is a factor of $∼2.5$ (⁠$Δ02=0.18$ vs 0.07 ps⁻² or Δ₀ = 0.42 vs Δ₀ = 0.26 ps⁻¹), and for Ala146N₃, they differ by a factor of $∼3.5$ (⁠$Δ02=0.52$ vs 0.15 ps⁻², i.e., Δ₀ = 0.72 vs Δ₀ = 0.39 ps⁻¹). Thus, the environmental dynamics around the spectroscopic label can be sufficiently perturbed by binding of a ligand in the protein active site to be reported directly as an experimentally accessible quantity with typical errors²⁴ between 0.1 and 0.3 cm⁻¹ (⁠$∼0.05$ ps⁻¹). Hence, the differences found from the simulations are well outside the expected error bars from experiment.

Nonvanishing static components of the FFCF were also reported from experiments. For trialanine (Ala)₃, a value of Δ₀ = 5 cm⁻¹ was reported,⁵⁹ compared with Δ₀ = 4.6 cm⁻¹ from MD simulations (0.94 ps⁻¹ vs 0.86 ps⁻¹) with multipolar force fields.⁵³ Similarly, CN⁻ in water features a nonvanishing tilt angle by τ = 10 ps⁶⁰ with Δ₀ ∼ 0.1 ps⁻¹ ∼ 0.5 cm⁻¹.⁶¹ Finally, 2D IR experiments for p-cyanophenylalanine bound to six distinct sites in a Src homology 3 domain reported static components ranging from Δ₀ = 1.0 to 3.7 cm⁻¹ (corresponding to 0.19–0.70 ps⁻¹).²⁴ 

Up to this point, the structural dynamics was characterized on the time scale required to converge the infrared spectroscopy within a given conformational substate. It was previously shown that IR spectra and FFCFs converge on the 1–2 ns time scale.⁵² During this time, the average C_α root mean squared deviation (RMSD) for the modified lysozymes compared with the WT reference structure is 1.5 Å with a minimal value of 0.85 Å (for Ala63N₃) and a maximum value of 2.6 Å (for Ala42N₃). However, attaching an azide label to alanine residues may incur structural changes of the protein on longer time scales. For example, it is known that lysozyme samples open and closed conformations,⁶² which occurs on the microsecond time scale.⁶³ 

In order to probe and quantitatively assess dynamics on longer time scales, additional simulations were carried out for ligand-free WT and modified lysozyme with buried and solvent exposed labels. The radial distribution functions g(r) for water around the middle nitrogen atom (NB; see Fig. S19) of the azide label show that the attached N₃ is either solvent exposed (residues 41, 42, 49, 63, 73, 74, 82, 93, 97, 112, 134, and 160) or buried (residues 98, 129, 130, and 146). First, two independent 10 ns simulations were run for Ala82N₃ (solvent exposed), Ala98N₃ (buried), and the WT protein (for comparison) using the RKHS representation for the –N₃ label [see Fig. S19 for g(r)]. Because the RKHS PES is computationally more expensive than evaluating an empirical energy function, a second set of simulations was run with an empirical parameterization for the N–N bond and the N–N–N angle. The N–N bond potential was harmonic with r_e = 1.14 Å and k_e = 877.4 kcal/mol, while for the N–N–N angular potential the parameters were θ_e = 180° and k_e = 46.7 kcal mol⁻¹ rad⁻²; otherwise, the energy function remained unchanged. For these simulations, the OpenMM⁶⁴ implementation of CHARMM was used. With this setup, one 100 ns simulation was carried out for each of the three systems.

For all these simulations, the root mean squared deviation (RMSD) and the RMSF with respect to the initial structure after heating were determined (see Figs. S20 and S21). The top panel in Fig. S20 reports the RMSF for WT (black), Ala82N₃ (blue), and Ala98N₃ (red) from the 10 ns production run with the N₃-label described by the RKHS PES. The results show that the residues around [70–80] and [105–116] have higher fluctuation for Ala98N₃ compared to the other two systems. Moreover, Ala82N₃ generally has lower RMSFs for all residues, except for Thr142. The bottom panel of Fig. S20 shows the corresponding RMSDs. Ala82N₃ stabilizes at $∼1.2$ Å, whereas for Ala98N₃ it increases to $∼2.8$ Å and WT lysozyme at ∼2.2 Å. This suggests that Ala82N₃ is structurally more stable than the other two systems on the 10 ns time scale.

Similarly, for the 100 ns simulations using the empirical energy function throughout, the bottom panel of Fig. S21 shows that the RMSDs for WT and Ala98N₃ stabilize at $∼5.0$ and $∼3.5$ Å, respectively, whereas Ala82N₃ stabilizes at $∼2.5$ Å. The RMSF results for 100 ns confirm that Ala82N₃ still fluctuates least compared with the other two systems. For some of the residues, such as [32–45], [65–83], [102–130], and [139–143], the WT lysozyme has larger fluctuations compared to the other two systems. For comparison, it is of interest to note that the RMSD from a 50 µs simulation for the M6I mutant of lysozyme⁶³ reached between 4 and 5 Å, and the experimentally reported⁶² RMSD between the open and closed structure for the I3P mutant is $∼5$ Å, all of which is consistent with the present findings.

Furthermore, the DCCMs for the WT and the two modified lysozymes were determined for the two sets of simulations (see Figs. S22–S24). For the WT protein, the DCCM from 10 to 100 ns simulations are comparable with some difference between residues [30–42] and [17–25] together with a few additional smaller features (see Fig. S22). For Ala82N_3, the DCCM from the 10 ns simulation using the RKHS PES for the label shows reduced anticorrelations but is overall similar to that from the 100 ns simulation using the conventional FF. On the other hand, for Ala98N₃, a pronounced coupling between residues [61–75] and [23–26] emerges (see Fig. S25), which was not found from simulations on the 10 ns time scale. In addition, the anticorrelations for some of the residues are less pronounced on the 100 ns time scale.

Finally, it is also of interest to determine the number of water oxygen atoms as a function of simulation time within a given cutoff (here 4 Å was chosen) of the middle nitrogen atom NB of the azide label for Ala82N₃ and Ala98N₃ (see Fig. S26). This confirms that throughout the 2 ns, the azide label is hydrated for Ala82N₃ and often “dry” for Ala98N₃. This is also consistent with the radial distribution functions (see Fig. S19). The findings from the simulations with the empirical force field for the entire system on the 100 ns time scale are similar: NB of Ala82N₃ is hydrated throughout the simulation, whereas for Ala98N₃ the middle nitrogen atom of the label is usually “dry” or is occasionally partially hydrated (see Fig. S27).

In summary, this work demonstrates that the 1D and 2D IR spectroscopy of the azide bound to alanine residues within one conformational substate for WT lysozyme provides valuable site-specific and temporal information about ligand binding of PhCN to the active site of WT lysozyme. The static component Δ₀ of the FFCF, which is an experimentally accessible observable, shows pronounced differences between the ligand-bound and ligand-free protein and can serve as a useful indicator for ligand binding. Changes in the maximum of the infrared absorbance are of the order of one to several cm⁻¹, which can be detected with state-of-the-art experiments.²⁴ The contrast between ligand-free and ligand-bound lysozyme increases when the azido-label is present, as demonstrated for Ala134N₃. The structural dynamics of the modified lysozymes depends on the position of the alanine residue along the chain. Given that even within the same conformational substate the spectroscopy and dynamics depend on the position of the modification, it is expected that when sampling a wider distribution of conformations, differences in the dynamics and spectroscopy persist. However, quantifying this is outside the scope of this work. Finally, this work also lays the foundations to investigate the binding of biologically relevant substrates, such as N-acetyl-D-glucosamine (NAG) that is hydrolyzed at the β-(1,4)-glycosidic bond by hen egg white lysozyme in Gram-positive bacteria.⁶⁵ 

See the supplementary material for Figs. S1–S27.

The authors gratefully acknowledge financial support from the Swiss National Science Foundation through Grant No. 200021-117810 and the NCCR-MUST.

The authors declare no conflict of interest.

The data that support the findings of this study are available from the corresponding author upon reasonable request.