Role of PheE15 Gate in Ligand Entry and Nitric Oxide Detoxification Function of Mycobacterium tuberculosis Truncated Hemoglobin N

The truncated hemoglobin N, HbN, of Mycobacterium tuberculosis is endowed with a potent nitric oxide dioxygenase (NOD) activity that allows it to relieve nitrosative stress and enhance in vivo survival of its host. Despite its small size, the protein matrix of HbN hosts a two-branched tunnel, consisting of orthogonal short and long channels, that connects the heme active site to the protein surface. A novel dual-path mechanism has been suggested to drive migration of O2 and NO to the distal heme cavity. While oxygen migrates mainly by the short path, a ligand-induced conformational change regulates opening of the long tunnel branch for NO, via a phenylalanine (PheE15) residue that acts as a gate. Site-directed mutagenesis and molecular simulations have been used to examine the gating role played by PheE15 in modulating the NOD function of HbN. Mutants carrying replacement of PheE15 with alanine, isoleucine, tyrosine and tryptophan have similar O2/CO association kinetics, but display significant reduction in their NOD function. Molecular simulations substantiated that mutation at the PheE15 gate confers significant changes in the long tunnel, and therefore may affect the migration of ligands. These results support the pivotal role of PheE15 gate in modulating the diffusion of NO via the long tunnel branch in the oxygenated protein, and hence the NOD function of HbN.


Introduction
Mycobacterium tuberculosis (Mtb) poses a serious threat to the public health worldwide, infecting nearly one third of the global population. The remarkable adaptability of tubercle bacillus to cope with hazardous level of reactive nitrogen/oxygen species within the intracellular environment contributes to its pathogenicity. An enhanced level of nitric oxide (NO) and reactive nitrogen species produced within activated macrophages during infection act as a vital part of host defense, limit the intracellular survival of Mtb, and contributes in restricting the bacteria to latency. Nevertheless, Mtb has evolved efficient resistance mechanisms by which toxic effects of NO and nitrosative stress can be evaded. One of the unique defense mechanisms by which Mtb protects itself from the toxicity of NO relies on the oxygenated form of truncated hemoglobin N (HbN), which catalyzes the rapid oxidation of NO to harmless nitrate [1][2][3]. Compared to horse heart myoglobin, the nitric oxide dioxygenase (NOD) reaction catalyzed by Mtb HbN is ,15-fold faster, suggesting that it may be crucial in relieving nitrosative stress [4].
Despite having single domain architecture, the NO-scavenging ability of Mtb HbN is comparable to flavoHbs that are integrated with a reductase domain and known to have a high NOD activity. It is thus important to understand what structural and dynamical features contribute to the efficiency of its enhanced NOscavenging function, and therefore ensure survival of the bacillus under nitrosative stress. X-ray crystallographic studies revealed that Mtb HbN hosts a protein matrix tunnel composed by two orthogonal branches [5,6]. In addition, computational simulations performed by some of the authors suggested that Mtb HbN has evolved a novel dual-path mechanism to drive migration of O 2 and NO to the distal heme cavity [7,8]. According to such a mechanism (Fig. 1), access of O 2 to the heme cavity primarily involves migration through the tunnel short branch (,10 Å long, shaped by residues in helices G and H). Binding to the heme then regulates opening of the tunnel long branch (,20 Å long, mainly defined by helices B and E) through a ligand-induced conformational change of PheE15 residue, which would act as a gate. It has been recently shown that the opening of PheE15 in the oxygenated protein is also affected by the N-terminal Pre-A motif [9].
Since the NOD function of HbN depends on the diffusion of NO to the O 2 -bound heme through the long tunnel branch, the PheE15 gate emerges as a fundamental residue in determining the overall efficiency of the NO scavenging. Accordingly, the NOD function of HbN must result from a balanced tuning of the opening/closing events of the gate. Moreover, the functional implication of PheE15 in assisting the NOD activity is supported by the preservation of this residue in mycobacterial HbNs, while it is replaced by other residues in truncated hemoglobins O and P [10,11]. However, to the best of our knowledge, no experimental data have yet been reported to examine the gating role of PheE15 and its influence on the NOD activity conducted by Mtb HbN. In this context, this study has been undertaken to probe the role of PheE15 in protein function. To this end, several PheE15 gate mutants have been tested experimentally for their NOD function. In addition, molecular dynamics (MD) simulations have been performed to analyze the structural changes in the topology of long tunnel and the alterations in the protein dynamics, paying attention to the ligand migration properties through the tunnel. Our results confirm the critical role played by E15 in ligand migration along the long channel.

Strains, Plasmids and Culture Conditions
Escherichia coli strains, JM109 and BL21DE3 were used for the cloning and expression of recombinant genes. Bacterial cultures were grown in Luria-Bertani (LB) or Terrific Broth (containing 24 g of Yeast Extract, 12 g of Bacto-Tryptone, 12.54 g of K 2 HPO 4 , 2.31 g of KH 2 PO 4 ) medium at 37uC at 180 r.p.m. When required, ampicillin and kanamycin (Sigma) were added at a concentration of 100 and 30 mg/ml, respectively. Plasmids, pBluescript (Stratagene) and pET28C (Novagen) were used for cloning and expression of recombinant genes as described earlier [3,12]. The oligonucleotides were custom synthesized by Integrated DNA Technologies Inc. NO (98.5%) was obtained from Sigma Aldrich and saturated NO was prepared as mentioned previously [13]. Heme content of the cell was measured as noted in previous studies [13].

Site-directed Mutagenesis and Construction of PheE15 Gate Mutants of HbN
Recombinant plasmid pPRN [3] was used as a source of HbN gene for the site directed mutagenesis. PheE15 mutants to Ala, Tyr, Ile or Trp were generated using a PCR approach. PCR amplified genes were cloned at NdeI-BamHI site of pET28c and expressed under T7 promoter as described previously [3]. Recombinant HbN and its mutant proteins were purified from the cell lysate of E. coli using metal affinity chromatography following standard procedures. Authenticity of mutants was confirmed after nucleotide sequencing. Figure 1. Representation of the long and short branches of the tunnel system. The graphical display is based on the X-ray crystallographic structure of Mtb HbN (PDB entry 1IDR), and the access routes of O 2 and NO in the dual-path ligand-modulated mechanism proposed for this protein are indicated. The gating residue PheE15 (residue 62) is shown in the two conformations found in the X-ray structure as sticks. doi:10.1371/journal.pone.0049291.g001

Measurements of Heme Content and NOD Activity
Total heme content was determined following the procedure described earlier [14]. Heme concentration was calculated from the absorption difference at 556 and 539 nm for the sodium dithionite-reduced and ferricyanide-oxidized sample. NOD activity of cells or purified protein was monitored polarographically as described previously [3,15]. NO consumption buffer assay contained 60 mM K 2 HPO 4 , 33 mM KH 2 PO 4 , 7.6 mM (NH 4 ) 2 SO 4 , 1.7 mM sodium citrate, 10 mM glucose and 200 mg/ml chloramphenicol. NO uptake rate of O 2 -bound HbN and its mutants was measured from the slope of curving traces recorded in the presence of specified concentration of NO following established protocols [3,15].

Measurement of Ligand Binding
O 2 and CO binding was checked from the absorption spectra of O 2 -and CO-bound species [9]. CO difference spectra were recorded between 350 to 600 nm after bubbling CO into the protein sample cuvette and recording the difference spectra against sodium-dithionite reduced protein. Oxygen equilibrium curves of HbN mutants were checked following the published procedure [16] to check their p 50 value. The association rate for CO binding to HbN mutants was determined by flash photolysis. A concentrated stock solution of deoxyHbN was diluted anaerobically (,100 mM) into a cuvette (1 mm path length) containing CO (1 mM). The fully liganded sample of ferrous HbN was photodissociated by 0.3 ms excitation pulse from a dye laser. The bimolecular rebinding time courses were collected as described elsewhere [17]. A minimum of five traces were collected and averaged for each experiment.

NO-oxidation by an Oxygenated Adduct of HbN and its Mutants
Wild type and mutant proteins were fully oxygenated by exposing the deoxygenated protein samples to air and checking their absorption spectra, which gave a specific Soret peak at 415 and two peaks, a and b, at 570 and 540 nm, very similar to oxy form of hemoglobin. With a gas tight Hamilton syringe, NO (5 mM) was sequentially added to the oxygenated protein (40 mM), and absorption spectra were recorded after each addition to follow the conversion into the oxidized form. The NO-induced oxidation of mutants was compared with the profile determined for the wild type protein.

Molecular dynamics simulations
The dynamical behaviour of the oxygenated form of HbN mutants was examined by means of extended MD simulations and compared to the results reported in previous studies for the wild type protein [7,8]. The X-ray structure of wild type Mtb HbN (PDB entry 1IDR, chain A, solved at 1.9 Å resolution) was used as  starting point for simulations. Mutants were generated by replacing PheE15 by Ala, Ile and Tyr in the X-ray structure of the wild type protein. In all cases simulations were performed using the same protocol adopted in our previous studies [7,8].
Briefly, the enzyme was immersed in a pre-equilibrated octahedral box of TIP3P [18] water molecules. The final systems contained the protein and around 8600 water molecules (ca. 28,270 atoms).
The system was simulated in the NPT (1 atm.; 298 K) ensemble using SHAKE [19] to keep bonds involving hydrogen atoms at their equilibrium length, periodic boundary conditions, Ewald sums for treating long range electrostatic interactions [20], and a 1 fs time step for the integration of Newton's equations. All simulations were performed with the parmm99SB force field [21] and employing heme parameters developed in previous works [7,8].
MD simulations were performed with the PMEMD module of the AMBER10 program [22]. The geometry of the models was relaxed by energy minimization carried out in three steps where hydrogen atoms, water molecules and finally the whole system were minimized. Equilibration was performed in successive 50 ps runs where the temperature was gradually increased from 100 K to 298 K in four steps at constant volume, followed by an additional step run at constant pressure for 100 ps. Then, a series of 100-150 ns MD simulations (at 298 K and 1 atm) were run.
The analysis of the trajectories was performed using frames collected every 1 ps during the production runs. Furthermore, 75 MD simulations (25 per mutant) were run to explore the pathways for ligand access (free NO in solution) to the heme cavity in PheE15Tyr, PheE15 Ile and PheE15Ala. To this end, five structures of the protein were taken from the last 25 ns of the trajectories run for the oxygenated mutants. These snapshots were used as starting points for unrestrained simulations run in presence of NO, which was placed at random positions around the protein (5 distinct random positions per protein snapshot). The systems were thermalized following the protocol mentioned above, and MD simulations were run up to 20 ns using the same simulation conditions.

Essential Dynamics
The dynamical behavior of HbN and its mutants was explored by means of essential dynamics [23,24]. Residues 1-15 were excluded as this region is very flexible and would mask the

Ligand Migration Profiles
The effect of PheE15 mutation on the migration of ligands was examined using two techniques. The preferred docking sites and migration pathways were identified using MDpocket [25]. Then, the migration free energy profiles were obtained using Multiple Steered Molecular Dynamics (MSMD) [26].
MDpocket is a pocket detection program that uses a fast geometrical algorithm based on a Voronoi tessellation centered on the atoms and the associated alpha spheres, which are clustered and filtered giving origin to pockets and channels [27]. These pockets are then used to identify docking sites and migration paths along the protein matrix. Analyses were performed using 10000 snapshots taken equally spread over the last 50 ns of the trajectories. The minimum and maximum alpha sphere radius was 2.8 Å and 5.5 Å , respectively. The identified cavities were superposed in time and space and a density map was generated from this superposition. Stable cavities are identified as highdensity 3D isocontours, while low-density isocontours denote transient or nearly non-existent cavities in the MD simulation.
MSMD simulations were run to evaluate the free energy profiles of ligand migration through the tunnels using Jarzynski's equality [26], which estimates the free energy from an ensemble of irreversible works along the same reaction coordinate. A steering potential forces the motion of the probe with constant velocity along the reaction coordinate. The reaction coordinate was the iron-ligand distance, the force constant was 200 kcal mol 21 Å 21 and the pulling velocity 0.025 Å ps 21 . The free energy profile of ligand migration along a tunnel was obtained following the computational scheme reported in previous studies [8]. At least twenty SMD simulations were run pushing the ligand through the tunnel from the solvent towards the iron. The starting snapshot for each SMD was taken from the final structure of an equilibrated MD simulation with the ligand placed at a fixed distance from the iron. Typically those distances are chosen from the preferred docking sites found in MDpocket analysis in conjunction with short MD simulations run to examine the motion of the free diatomic ligand along the tunnel.

Results and Discussion
Four gate mutants, where PheE15 in wild type (wt) HbN was replaced by Ala (PheE15Ala), Ile (PheE15Ile), Tyr (PheE15Tyr) and Trp (PheE15Trp), were created by site-directed mutagenesis. Mutations were chosen to span a wide range of sizes, varying from the small methyl group in Ala to the large indole ring in Trp, the replacement of the planar benzene by the branched chain of Ile, and the conservative mutation of PheE15 by Tyr. It was expected that the distinct chemical nature of the side chains would translate into differences in the ligand migration through the long tunnel, which in turn should lead to differences in the NOD activity measured for the mutants. Thus, in the absence of relevant structural alterations in the tunnel due to the mutations at position E15, which might be relevant for the bulky Trp, it was expected that replacement of PheE15 to Ala should open permanently the long tunnel, whereas mutation to Trp should occlude the access of ligands. Likewise, the branched side chain of Ile was expected to limit the accessibility of diatomic ligands through the tunnel. Finally, the conservative PheE15Tyr mutation was a priori expected to have little effect on the migration properties.

Effect of Mutations at PheE15 Gate of HbN on O 2 /CO Binding
Absorption spectra of O 2 and CO bound forms of mutants were indistinguishable from that of wt HbN (Fig. 2), suggesting that all these mutants bind O 2 . The ability to bind O 2 was further assessed by measuring the p 50 values (Table 1). For the wt protein, the p 50 value was 0.019 (mm Hg), which compares well with previous data [1]. When the PheE15 gate mutants were compared with wt HbN, no significant difference in their p 50 profile was observed suggesting that the O 2 binding properties of mutants are not impaired. To substantiate these results, the rate of CO association was also determined to characterize the kinetics of ligand association ( Table 1). The k on value determined for the wt protein

Mutations at PheE15 Gate of HbN Alter its NOD Activity
Even though mutants exhibit very similar O 2 binding properties, relevant differences were observed in their NO metabolizing activities ( Table 2). The NOD activity of PheE15Ile was reduced to around 55% of the activity measured for the wt enzyme, whereas a larger reduction (around 70%) in NOD activity was observed for the PheE15Trp mutant. Unexpectedly, mutation of PheE15 to Tyr also yielded a significant reduction (around 65%) in the NOD activity, suggesting that the apparently conservative replacement of benzene by phenol has a drastic influence on the ligand migration properties of the mutant. Finally, the PheE15Ala mutant also exhibited a significant reduction (around 49%) in the NOD activity.

Oxidation of NO by Oxygen Adduct of PheE15 Gate Mutants of HbN
The NO oxidation profile of wt HbN and its mutants was determined by titrating the oxy forms of protein with NO in a time course manner (Fig. 3). The addition of NO (5 mM) resulted in the appearance of a partially oxidized spectrum of HbN and repeated addition of NO solution to this sample resulted in a fully oxidized spectrum shifting Soret peak of oxyHbN from 415 to 405 nm. In contrast, similar additions of NO to the oxygenated form of mutants PheE15Trp, PheE15Ile and PheE15Tyr did not change the spectra to the oxidized form immediately and changed slowly after 30-35 min of exposure, indicating very slow NO oxidation. Finally, NO oxidation by PheE15Ala mutant displayed a spectral profile similar to PheE15Ile and PheE15Tyr, which is intermediate between those observed for the wt protein and the PheE15Trp mutant.

MD Simulations
The preceding data indicate that the different PheE15 mutations do not affect the binding of O 2 to the heme. However, since the mutants exhibit a distinctive reduction in the NOD activity, the PheE15 residue has to play a key role in mediating the access of NO to the oxygenated protein. On the basis of these findings, MD simulations were run with a twofold purpose: to examine the structural integrity of the overall protein fold, and to identify local changes in the long branch tunnel that might affect the ligand migration. To this end, a series of 100-150 ns MD simulations were run for the heme-bound O 2 forms of PheE15Ala, PheE15Ile and PheE15Tyr mutants, and the results were 7. Representative trajectories followed by free NO ligand through the protein matrix. MD simulation of free NO located at random positions from the protein were followed to investigate the pathways leading to the heme cavity for oxygenated forms of (A) PheE15Ala, (B) PheE15Ie and (C) PheE15Tyr. The distinct pathways are indicating by showing the position of NO (represented as sticks) along the trajectory using different colors: long branch (blue), short branch (brown), EH (green), and other (yellow). Note that helix G is displayed as highly transparent cartoon for the sake of clarity. doi:10.1371/journal.pone.0049291.g007  [7][8][9]. Choice of the simulation time, however, was also dictated by consistency in the analysis of the structural and dynamical properties of the simulated systems, either related to the overall features of the protein backbone or to the ligand migration through the protein matrix. Such a convengence was achieved for all the simulated systems but for the PheE15Trp mutant (see below). For the particular case of the PheE15Tyr mutant, two MD simulations differing in the starting orientation of TyrE15 were run. Thus, on the basis of the conformational preferences found for PheE15 in both the X-ray structure [5] and previous MD studies [7,8], the phenol ring was oriented in the closed or open conformations, which prevent or facilitate ligand migration through the tunnel, respectively. When TyrE15 was in the closed orientation, the side chain remained stable, and no transitions between open and closed states were found along the whole trajectory (Fig. S1). There was only a conformational change leading to the transient formation of a hydrogen bond between the hydroxyl group of TyrE15 and the carbonyl group of Ile115. When the MD simulation started from the open conformation, the side chain of TyrE15 suddenly changed to the closed conformation after the first 12 ns and then remained stable until the end of the trajectory (Fig. S1), even though few attempts to form transient hydrogen-bond interactions with Ile115 can be observed. Since the two simulations exhibited the same behavior, hereafter discussion of PheE15Tyr mutant will be limited to the trajectory starting from the closed conformation.
Finally, the MD simulation run for the PheE15RTrp mutant was extended up to 200 ns. However, the results will not be presented here as the analysis of the trajectory points out that the structural changes induced by the mutation both in the tunnel and in the global protein structure are still not fully converged. It seems that a proper description of the structural rearrangements triggered by this mutation would require longer simulations, making it necessary to be cautious for not overinterpreting the structural changes due to the PheE15RTrp mutation, which cannot be easily accommodated in the tunnel (see Fig. S2 for details).

Structural Analysis
For all the simulations, inspection of both the time evolution of the potential energy (Fig. S3) and the rmsd ( Fig. 4; see also Fig. S4) determined for the backbone atoms of the protein core (excluding the first 15 residues) supports the integrity of the simulated proteins. Thus, the rmsd of the protein core (ranging from 1.2 to 2.1 Å for the different mutants) remains very stable after the first few nanoseconds. In contrast, a much larger rmsd profile showing notable fluctuations along the trajectory is obtained when the whole protein backbone is included in the analysis. This finding indicates that the pre-A segment is very flexible in all the mutants and can adopt a large number of conformational states, as noted in the formation of distinct structural arrangements for residues 1-15 (data not shown). Therefore, the mutation does not alter the large conformational flexibility found for the pre-A segment in the wt HbN [9,28].
The geometrical arrangement of residues TyrB10 and GlnE11 in the distal cavity is a structural feature of particular relevance, as it has been proposed that opening of the gate is modulated by the oxygen-sensing properties of the TyrB10-GlnE11 pair [8,29]. In the oxygenated wt HbN these residues form a network of hydrogen bonds, where the hydroxyl group of TyrB10 is hydrogen-bonded to the heme-bound O 2 , and also accepts a hydrogen bond from the side chain of GlnE11. Due to these interactions, the GlnE11 side chain lies closer to the benzene ring of PheE15 than in the deoxygenated protein, and the enhanced steric clash favors opening of the PheE15 gate upon O 2 binding. Accordingly, it might be argued that disruption of the hydrogen bonds formed by TyrB10 and GlnE11 could explain the reduction in NOD activity of the mutants. However, the analysis of the trajectories completely rules out this possibility (Fig. S5). Thus, in all the mutants the hydroxyl group of TyrB10 is hydrogen-bonded to the hemebound O 2 (average distances about 2.85 Å ), and the side chain amide nitrogen of GlnE11 is hydrogen-bonded to the TyrB10 hydroxyl group (average distances about 3.0 Å ). These interactions reproduce the hydrogen bonds found in the X-ray structure, as the corresponding distances (averaged for subunits A and B) are 3.15 and 2.95 Å . Therefore, the hydrogen-bond network found in wt HbN is not affected by the PheE15 mutations, and hence one should expect that the steric pressure exerted by GlnE11 is retained in the mutated proteins.

Dynamical Analysis of the Protein Backbone
Binding of O 2 to the heme also changes the dynamical motion of the protein backbone [7,8]. Thus, essential dynamics analysis of oxygenated HbN reveals that the major motions in the deoxygenated protein affect helices C, G and H, while the largest contribution to protein flexibility comes from helices B and E in the oxygenated protein. Since helices B and E define the walls of the tunnel long branch, the increased motion of these helices should facilitate the transition between open and closed states of the gate, thus influencing the ligand migration through the tunnel. This dynamical alteration agrees with the large scale conformational change observed experimentally upon binding of NO to the ferric form of wt HbN [30], and with the occurrence of distinct conformational relaxation processes found in the kinetics of CO recombination to the protein encapsulated in gels [31].
In order to investigate whether mutation of the PheE15 gate influences the dynamics of the protein backbone and eventually affects the migration of NO in the oxygenated protein, we determined the essential dynamics of the mutants and compared them with the wt protein. Diagonalization of the positional covariance matrix for the backbone atoms points out that few motions account for a significant fraction of the protein dynamics. Nearly 50% and 70% of the backbone conformational flexibility is accounted for by the first 4 and 10 principal components (Table  S1). The two first essential motions, which mainly involves motions of helices B, E and H, loop F and the hinge region around helix C, account for 22-37% of the structural variance in the backbone of the mutants, in agreement with the value (36%) determined for wt HbN.
The similarity between the structural fluctuations of the protein backbone in wt HbN and its mutants was measured by means of the similarity index j AB (see Text S1), which takes into account the nature of the essential motions and their contribution to the structural variance of the protein [32]. When the 10 most relevant motions are considered, the similarity index varies in the range 0.60-0.65, which is slightly lower than the self-similarities obtained for wt HbN and mutants (Table 3). These results point out that mutations preserve to at large extent the dynamical behavior of the protein. For the sake of comparison, the similarity indexes determined for the 10 most essential motions between wt HbN and the mutants TyrB10Phe and GlnE11Ala only amount to 0.47 and 0.38, respectively [27], indicating that the hydrogen-bond network formed by the TyrB10-GlnE11 pair has a larger impact on the protein dynamics than mutation of PheE15.

Ligand Migration
Since the preceding results did not reveal significant changes neither in the hydrogen-bond network of TyrB10-GlnE11 pair nor in the dynamical behavior of the protein skeleton, the reduction in NOD activity determined for PheE15 mutants might reflect a reduced accessibility of ligands for migrating through the long tunnel. To corroborate this hypothesis, we determined the shape of internal tunnels by using MDpocket, which provides a grid that encloses the protein, where each grid point is assigned an occupancy value that denotes the accessibility of the volume associated to that point [23].
For the wt HbN the tunnel long branch appears as a continuous cavity leading from the protein surface to the heme cavity, as noted by the continuous progression of the accessible volume isocontour shown in Fig. 5. In contrast, the representation of the same isocontour for the mutants reveals a discontinuous channel, even in the case of the mutation to Ala.
For mutants PheE15Ile and PheE15Tyr, such a disruption affects the region that surrounds the mutated residue, suggesting that ligand migration through the tunnel long branch is impeded by the steric hindrance imposed by the side chains of the mutated residues (Fig. 5C,D). This is confirmed by the free energy profiles determined by MSMD calculations for these mutants, where the limiting step for ligand migration is associated to surpassing the gate (Fig. 6). In fact, the free energy profiles are intermediate between those obtained for the wt HbN with the PheE15 gate in either closed or open conformational states. Compared to the wt protein in the open state, mutation of the gate slightly destabilizes the minimum located at around 11 Å (i.e., the docking site before the gate), and increases the barrier required to pass above the side chain of the mutated residue in order to access the heme cavity. Therefore, a small diatomic ligand such as NO is expected to be trapped effectively at the highly hydrophobic entrance of the long branch, thus enhancing the local concentration, but access to the heme cavity is mainly limited by the hindrance due to the side chain of the mutated gate.
In contrast with the preceding findings, the isocontour is continuous around the mutated gate in the PheE15Ala mutant, as expected from the small side chain of Ala. However, compared to the wt protein (Fig. 5A), there is a significant occlusion at the entrance of the channel (Fig. 5B), which stems from a slight readjustment of helices B and G that reduces the width of the channel entry by 0.4-0.6 Å (as measured from distances between Ca atoms of residues located at the helical ends). Accordingly, the probability of the ligand to be trapped at the entrance of the tunnel long branch is lower compared to wt HbN. In turn, this finding raises the question about the existence of alternative pathways that can justify the remaining NOD activity (around 51% compared to wt protein) retained by PheE15Ala mutant. To this end, we examined the trajectory followed by a free NO in a water box containing the PheE15Ala mutant and determined the routes leading to the heme cavity in 25 independent MD simulations. The results (see Table 4) showed that 18 out of 25 trajectories were successful in allowing the ligand (NO) to achieve the heme cavity within the simulation time. Nevertheless, only 5 trajectories showed that NO migrates via the tunnel long branch, and in 7 simulations the entry pathway involved the EH tunnel (defined by residues in helices E and H) reported by Daigle et al. [33]. Remarkably, in 3 trajectories NO was able to reach the heme cavity through the short branch, and in other 3 trajectories NO accessed the heme cavity via a distinct channel (see Fig. 7A). Thus, even though the reduction in NOD activity can be related to the occlusion of ligand access at the beginning of the tunnel long branch, this effect is counterbalanced by the existence of three alternative entry pathways, which arise from slight structural alterations in the protein skeleton. This finding reveals the delicate structural balance imposed by the PheE15 gate, which not only regulates ligand migration, but also contributes to avoid the collapse of helices B and E, thus preserving the structural integrity and ligand accessibility along the tunnel long branch.
For the sake of completeness, this analysis was also performed for PheE15Ile and PheE15Tyr mutants. For these mutants the access routes leading to the heme cavity were significantly different. Thus, the trajectories primarily involved migration through the tunnel long branch (Table 4), which should be considered the main pathway for ligand migration in the oxygenated protein. In this pathway, the ligand remained docked in a region before the gate for significant periods of time (as noted in the high density of NO molecules located before the mutated gate residue; see Fig. 7B,C) until it was able to surpass the barrier due to the side chain of Ile/Tyr and access the heme cavity. There were some attemtps to get the heme cavity through the short tunnel, but they were unsuccessful. Finally, only in very few cases the ligand entered the protein through the tunnel EH, leading to the docking site located before the mutated gate residue. Therefore, these findings are in contrast with previous modeling studies where the short tunnel was reported to be the main route for NO diffusion [34]. Furthermore, present results indicate that the main route in PheE15Ile and PheE15Tyr is the tunnel long branch, which supports the role of the tunnel long branch in the dual-path migration mechanism [7,8] and particularly the functional relevance of the PheE15 gate.

Conclusion
The experimental data collected for the mutants indicate that the PheE15 mutation does not affect the binding of O 2 to the heme, as noted in the similar absorption spectra of O 2 and CO bound forms, as well as in the similar p 50 and k on values determined for both wt and mutated proteins. These findings suggest that the PheE15 mutation has little impact on the tunnel short branch, which is proposed to be the main pathway for migration of O 2 to the heme cavity. However, since the mutants exhibit a distinctive reduction in the NOD activity, PheE15 residue has to play a key role in mediating the access of NO to the oxygenated protein. In agreement with the dual-path migration mechanism [7,8], the main pathway for NO migration is the tunnel long branch. This finding points out the delicate structural balance imposed by the PheE15 gate, which not only regulates ligand migration, but also contributes to avoid the collapse of helices B and E, thus preserving the ligand accessibility along the tunnel long branch. Overall, the results presented herein demonstrate the pivotal role of PheE15 in modulating the NOD function of Mtb HbN, thus confirming the suggestion by Milani and coworkers [5] about the gating role of PheE15 in HbN. Table S1 Structural variance (%) of the protein backbone. The contribution of the first essential motions to the structural variance is gindicated for the mutated forms of HbN. The cumulative value is given in parenthesis.