Introduction

Antimicrobial resistance (AMR) is a global health problem1. Many efforts have been made to reduce the burden of AMR perils globally since 2013, yet threats from some species continue to rise regardless: drug-resistant Neisseria gonorrhoeae is one of five urgent threats2,3. Resistance to ceftriaxone, the last option for an empirical first-line antibiotic against Neisseria gonorrhoeae in most countries, has been reported and continues to emerge globally4. The gonococcal infection could become untreatable due to a high degree of AMR, which would increase serious complications: infertility, ectopic pregnancy, and increased transmission of HIV. The emergence of resistant pathogens to currently available antibiotics is very alarming; thus, the development of treatment options is imperative to tackle AMR.

The respiratory chain has recently garnered considerable scientific attention as a potential target for antibiotics. As a weapon to overcome AMR, compounds targeting the respiratory chain have been approved or entered clinical trials, for example, drugs against parasites, fungi, and particularly drug-resistant Mycobacterium tuberculosis5,6,7,8,9,10. However, most of them are competitive inhibitors of orthosteric sites. As respiratory enzymes are essential for life, their core structure is generally conserved across species. Structural similarity and substrate commonality with host proteins are risks for cross-reactivity, which could be a cause of side effects11. Therefore, an allosteric inhibitor is a more feasible choice as allosteric sites are evolutionarily less conserved in amino acid sequence than orthosteric sites, theoretically improving selectivity and reducing toxicity12,13. However, a systematic and strategic search for allosteric inhibitors has not been established yet, especially against membrane proteins; most respiratory enzymes are membrane proteins.

HCOs are terminal respiratory enzymes present in all three domains of life: bacteria, archaea, and eukaryotes. HCOs receive electrons from the respiratory chain and reduce molecular oxygen to water. This exergonic reaction is coupled with proton pumping across the membrane, which contributes to maintaining the proton motive force that is further used for ATP production14,15,16,17,18. HCOs are multi-subunit complexes, and their constitution varies among species; however, subunit I is a catalytic subunit common in all HCOs. It contains a low-spin heme and a binuclear center (BNC), the catalytic site formed by a high-spin heme and a copper ion. The low-spin heme first receives electrons and transfers them to the BNC for the reduction of oxygen16,17,18,19.

Eukaryotes originate from the symbiosis between Alphaproteobacteria and archaea20. Therefore, mitochondrial DNA-encoded subunits I to III of mitochondrial cytochrome c oxidase (mtCcO), analogous to eukaryotic HCOs, are descendants of the respiratory enzymes from bacteria20,21, and functionally important residues, and thus their core structure is conserved, although the remaining residues are not the same. In addition, mammalian mtCcO has 10 more subunits encoded by genomic DNA. The physiological roles of these subunits are not fully elucidated16. The core structures of fundamental proteins such as RNA polymerases or ribosomes are also similar among species; intriguingly, they have acquired additional subunits that modulate their function along molecular evolution22,23. Thus, we hypothesized that the surface of the core structure of HCOs, which are covered by additional helices in mammals, might contain allosteric sites, regulating their activity positively or negatively. The existence of an additional helix in mammalian mtCcO makes the pockets distinct from bacterial HCOs.

In this work, we identify an allosteric inhibitory site buried inside eukaryotic mtCcO. The steric conformation around the binding pocket of HCOs is highly conserved among bacteria and eukaryotes, yet the latter has an extra helix. The structural difference in the conserved allostery enables us to rationally identify bacterial HCO-specific inhibitors: an antibiotic compound against ceftriaxone-resistant Neisseria gonorrhoeae.

Results

An allosteric inhibitory site buried inside mtCcO

To test this hypothesis, first, we need to identify an allosteric inhibitory site. We started with mammalian mtCcO inhibitors obtained by random compound screening. We have previously found that an endogenous protein directly interacts with mtCcO and allosterically modulates mtCcO activity24,25. This finding led us to perform random compound screening that modulates mtCcO activity; we identified mtCcO inhibitors, chemically distinct from the known inhibitors, including carbon monoxide, nitric oxide (NO), or cyanides. We selected several allosteric inhibitors after studying their enzyme kinetics (Supplementary Fig. 1a–d). We then attempted to obtain complex crystal structures of mtCcO and our inhibitors. We focused on T113 hereafter as its binding site was buried within a mammalian-specific helix, COX7C (Fig. 1a, b). X-ray diffraction data with a resolution of 2.2 Å were obtained from a mtCcO crystal soaked with T113. We also determined the apo-structure of mtCcO under the same preparation conditions (Supplementary Table 1). The obtained complex structure showed a clear compound binding site that gave additional electron density compared with the protein, clearly found inside the mtCcO (Supplementary Fig. 2a). The Fo(T113)–Fo(DMSO) differential map confirmed that this electron density did not originate from water or lipid molecules, and showed the highest difference (Supplementary Fig. 2b). The binding pocket was different from the binding site for molecular oxygen or cytochrome c, the route for electron transfer pathway, proton pathway, or oxygen accessing channel16,19, suggesting that T113 is a genuine allosteric inhibitor.

Fig. 1: The allosteric site for T113 is buried inside eukaryotic mtCcO.
figure 1

a X-ray structure of mtCcO with T113. The electron density map (2Fo–Fc), contoured at 1σ, is shown in the left. The ribbon model of mtCcO with T113 in the sphere is shown on the right. T113 was covered by COX7C (red), hidden from the surface. b T113 was surrounded by 4 transmembrane helices (TM1–3, COX7C) and buried from the surface, viewed from the inter-membrane space. Protein molecular surface is shown as gray in the close-up view. Three helices of subunit I surrounding the allosteric site are shown as dark blue, the other helices of subunit I as pale blue, subunit COX7C in mtCcO is shown as red.

Full size image

Three out of four helices surrounding the allosteric site of mtCcO belong to subunit I, common in HCOs. Furthermore, we noticed that the steric conformation of the helices near the low-spin heme is well conserved in bacterial HCOs (Fig. 2a and Supplementary Fig. 3). These observations led us to assume that we could reasonably screen allosteric inhibitors from derivatives of our mtCcO inhibitors, which likely target the corresponding allosteric site of bacterial oxidases.

Fig. 2: Identification of bacterial HCO-specific allosteric inhibitors.
figure 2

a Allosteric sites of mitochondrial mtCcO (left) and E. coli bo3 UqO (center). Helices of subunit I are shown as blue, TM0 of subunit I as purple, subunit COX7C is shown as red, and the other helices as yellow. Merged views are shown (right) with E. coli bo3 UqO in gray. COX7C is close enough to cover the inhibitor. All the helices in subunit I are merged between mtCcO and E. coli bo3 UqO except TM0. b Screening of 434 chemicals at 50 µM against E. coli bo3 UqO. Data are presented as an average of a duplicate. c Dose-dependent inhibition of N4 specific on bo3 UqO enzymatic activity. Data are presented as an average value of technical replicates over two independent experiments. d Kinetic analysis of bo3 UqO with DMSO and 25 µM N4 molecule. Fitting lines are calculated by using Michaelis–Menten equation with the non-competitive inhibition model. Data are presented as an average value of technical replicate. Reproducibility was confirmed by two independent experiments. Source data are available as a Source Data file.

Full size image

Rational identification of bacterial HCO-specific inhibitors

For the preparation of a custom library, we applied in silico compound screening originated from two mtCcO inhibitors, including T113 and T151, which we obtained in the initial high throughput screening; T151 was also found to bind the allosteric site (Supplementary Fig. 1d). Structurally similar compounds to those were primarily collected by multiple ligand-based search algorithms from 80 million commercially available compounds. Then, our in-house algorithm integrated them, ranked them, and chose the first series of 285 compounds26,27. We added the second set of 149 compounds chosen by docking simulation, which screened the same 80 million compounds against the allosteric pocket of mtCcO. In total, we established a custom library consisting of 434 compounds that have a high probability of binding to the conserved allosteric site for HCOs. We tested this library against mtCcO, and as a result, 47 compounds inhibited mtCcO more than 40% at 50 µM (11.4%), higher than the usual hit rate from random screening, verifying that our custom library concentrated mtCcO inhibitors (Supplementary Fig. 4a).

We used E. coli bo3 ubiquinol oxidase (bo3 UqO) to test our hypothesis as a model bacterial HCO. Our custom library was screened against bo3 UqO, and we obtained 15 hit compounds that showed >40% inhibition for bo3 UqO at 50 µM (Fig. 2b). Among them, eight common inhibitors for both mtCcO and bo3 UqO, and more importantly, two specific inhibitors for bo3 UqO, were successfully acquired (Fig. 2c). As expected, one of these inhibitors, N4, fitted well with the allosteric inhibitory curve assessed by the Michaelis–Menten equation (Fig. 2d).

To obtain direct evidence that N4 binds to the corresponding allosteric site, we determined the structure of N4-bound bo3 UqO with a Fab fragment at 3.0 Å resolution using cryogenic electron microscopy (cryo-EM) (Fig. 3a, Supplementary Fig. 5e–h, and Supplementary Table 2). We also determined the apo-structure of bo3 UqO at 3.1 Å under the same preparation conditions (Supplementary Fig. 5a–d). Differential maps between holo- and apo-structures demonstrated explicit additional density, and this density was the top difference found (Supplementary Fig. 5i). Notably, the binding site was exposed to the surface and adjacent to transmembrane helix 1 (TM1), TM2, and TM3 of subunit I, which strongly corroborates our hypothesis (Fig. 3b, c). There were hydrogen bonds between Asp75, Arg71, and N4. The electrostatic potential surface of the binding sites showed that bo3 UqO has a more hydrophilic surface than mtCcO, suggesting it is unlikely that hydrophobic T113 binds to bo3 UqO (Fig. 3d). Meanwhile, relatively hydrophilic N4 is not a feasible binder for the mtCcO allosteric pocket, a more hydrophobic environment, explaining that N4 is a derivative of T113, yet they are mutually exclusive inhibitors for bo3 UqO or mtCcO, respectively. Mutants for amino acid residues around the inhibitor in the structure demonstrated a significant change in the inhibitory effect, confirming that N4 definitively binds bo3 UqO at the pocket (Supplementary Fig. 5j).

Fig. 3: Cryo-EM proof of the allosteric site of bo3 UqO.
figure 3

a Cryo-EM density map (left) and the ribbon model (right) of bo3 UqO in complex with N4 (red in the left, sphere in the right). The allosteric site of bo3 UqO is exposed. b N4 molecule in the allosteric site of bo3 UqO. The cryo-EM density map around N4 is shown in red. Molecular interactions between Asp75, Arg71, and N4 are shown as dotted lines. c N4 is accessible from the surface. The protein molecular surface at the plane of Fe in heme b is shown as gray. The whole N4 molecule was presented. d The electrostatic potential surface of the allosteric site of bo3 UqO (left) and mtCcO (right). e Growth inhibition by N4 for 24 h culture of wild-type E. coli (Wild), bo3 UqO-knockout strain (Δbo3), and bo3 UqO-dependent strain (Δbd). f A time- and concentration-dependent viability assay. A dotted line indicates inoculation (1 × 105 CFU/ml). Data are confirmed by biological triplicate (e) or duplicate (f). The color of the helices are as same as in Fig. 2. Source data are available as a Source Data file.

Full size image

A specific inhibitor of bacterial HCOs might have the potential as an antibiotic. To test this possibility, we assessed the inhibition of E. coli growth by N4. In E. coli, there are two branches for terminal oxidases in the respiratory chain, bo3 UqO and cytochrome bd oxidase (bd UqO). E. coli uses bo3 UqO in aerobic conditions, whereas it preferentially uses bd UqO in hypoxic conditions, as an adaptation to the environment28. In wild-type E. coli, there was no effect of N4 in growth; however, it significantly decreased E. coli growth in bo3 UqO-dependent strain (Fig. 3e). This growth inhibition was canceled in the bo3 UqO knockout strain, confirming that the growth inhibition by N4 was a result of HCO inhibition, not a nonspecific effect. To distinguish between the bactericidal and bacteriostatic effects of the compound, we performed a colony count assay and found that the effect of our inhibitor, N4, on E. coli is bacteriostatic (Fig. 3f).

Development of an antibiotic compound for Neisseria gonorrhoeae

Next, to expand our findings, we tested quinol-dependent NO reductase (bb3 qNOR), a distant family member of HCO from the pathogenic bacteria, Neisseria. In the genus Neisseria, the emergence and spreading of multi-drug resistant N. gonorrhoeae is a considerable global health concern29. qNOR reduces NO to nitrous oxide (N2O) and is a critical enzyme in denitrification, which oxidizes nitrogen compounds to produce energy in anoxic conditions30. Bacterial denitrification also plays a vital role in protecting them from exogenous NO produced by host immune cells, thereby being implicated in the pathogenicity of several bacterial species, including N. meningitidis and N. gonorrhoeae31,32. These suggest that qNOR in Neisseria can be a drug target. We used N. meningitidis qNOR because the amino acid sequences of N. meningitidis qNOR and N. gonorrhoeae qNOR are 98% identical, especially 100% in the transmembrane region, and importantly the structure of N. meningitidis qNOR has been solved33. We confirmed that the steric conformation of three helices (TM 3–5 in bb3 qNOR) around the low-spin heme is also conserved in bb3 qNOR (Fig. 4a)34. We employed the same approach to bb3 qNOR as for bo3 UqO. First, we screened our custom library against bb3 qNOR and obtained 52 hit compounds that showed >50% inhibition of bb3 qNOR at 50 µM (Fig. 4b). Among them, we found one bb3 qNOR-specific inhibitor, Q275. Q275 did not cross with either mtCcO or bo3 UqO (Fig. 4c). Q275 did not affect the cellular respiration of mammalian cells, nor cell viability of them (Fig. 4d, e). We confirmed that Q275 has an allosteric mode of inhibition (Fig. 4f). Further, amino acid-substituted mutants around the corresponding allosteric site of bb3 qNOR showed a considerable change in the inhibitory effect, confirming that Q275 binds to bb3 qNOR at the pocket (Fig. 4g).

Fig. 4: Identification of N. gonorrhoeae bb3 qNOR-specific inhibitor.
figure 4

a The steric conformation around the low-spin heme is conserved in bb3 qNOR (PDB: 6fwf). A red shadow indicates the allosteric site, distinct from the quinol binding site. b Screening of 434 chemicals at 50 µM on N. meningitidis bb3 qNOR (x-axis) and E. coli bo3 UqO (y-axis). Q275 is a specific bb3 qNOR inhibitor, shown in green. c Dose-dependent and specific inhibition of Q275 on bb3 qNOR activity. d Q275 has no effect on oxygen consumption rate in mammalian cells. N = 6 technical replicates for each group and reproduced in two independent experiments. e Q275 shows no cytotoxicity assessed in rat cardiomyocytes. N.C.; DMSO, P.C.; 1% Triton X-100. N = 6 technical replicates for DMSO, 3 for other groups. Reproducibility was confirmed by two independent experiments. f Kinetic analysis of bb3 qNOR with DMSO and 10 µM Q275 molecule. Fitting lines are calculated by using Michaelis–Menten equation with the non-competitive inhibition model. g The effect of amino acid substitution of bb3 qNOR. The positions of Val293, Val289, and Trp334 are shown in (a). h Minimum inhibitory concentration (MIC) of Q275 against the N. gonorrhoeae reference strain WHO F and a ceftriaxone-resistant FC428. Sodium nitrite (20 mM) was added to mimic a NO condition. The data represent duplicate experiments. i Colony count 24 h after Q275 treatment of the N. gonorrhoeae WHO F, FC428, and norB-deficient from both strains. Sodium nitrite (20 mM) was added to mimic a NO condition. A dotted line indicates inoculation (5 × 105 CFU/ml). c, f, g Data are presented as an average value of technical replicates over three (f, g) independent experiments. Source data are available as a Source Data file.

Full size image

Finally, we tested the antimicrobial effect of Q275 against a clinically isolated strain of N. gonorrhoeae that has resistance to ceftriaxone (FC428); the spread of the super-resistant N. gonorrhoeae has become a serious global health problem4. Notably, Q275 demonstrated antimicrobial effects against FC428 as well as the reference strain (WHO F) in NO-challenging (attack from immune cells) conditions, mimicking the infectious intravital environment (Fig. 4h). We established norB-deficient, encoding qNOR, strains in both WHO F and FC428 backgrounds. In NO-challenging condition (20 mM NaNO2), norB-deficient strains in both WT and FC428 did not grow, further supporting that growth inhibition is mediated via qNOR. To distinguish between bactericidal and bacteriostatic effects using these norB-deficient strains, we performed colony count after exposure to NaNO2 and found that targeting qNOR in N. gonorrhoeae is bacteriostatic (Fig. 4i). These results suggest that our approach developed a specific inhibitor against a pathogenic bacterial HCO on-demand with a narrow range of specificity, having the therapeutic potential to tackle AMR.

Allosteric inhibition in the substrate accessing channel

We next focused on the inhibitory mechanism. The binding site of T113 on mtCcO formed a narrow tunnel surrounded by four helices in total (Fig. 1b). The Fo(T113)–Fo(DMSO) differential map revealed several differences in the structure of mtCcO triggered by the binding of the inhibitor: the exit site of the proton pathway (Asp50, Asp51), a side chain of Ser382 of TM10 in subunit I, and hydroxyfarnesylethyl group of heme a (Supplementary Fig. 2b). These are the sites where structural change has been reported between the reduced and oxidized forms of mtCcO, raising the possibility that T113 might have reduced mtCcO35. To test this possibility, we analyzed the signature of mtCcO reduced through dithionite or by T113. Compared to a fully reduced state induced by dithionite, mtCcO mixed with T113 showed a minimum signature of reductive change in the Soret band and an increase in around 600 nm absorption spectra, which are hallmarks of reduced hemes, suggesting that T113 is not a reductive reagent (Supplementary Fig. 6). Therefore, the structural change we observed does not fully explain the mechanism of inhibition.

In the mechanism of allostery, an effector binding transmits the signal to the functional site, the orthosteric site, by a transition of conformational ensembles, which is often difficult to capture by structural analysis because of its snapshot nature or possible constraint in crystallization36. Therefore, to obtain a mechanistic insight, we applied molecular dynamics (MD) simulation of mtCcO with or without the inhibitor. Overall, MD simulations with the inhibitor (holo-MD) did not show significant structural deformation. Notably, the average trajectories with the inhibitor showed that TM2 of subunit I, one of the four helices forming the inhibitor binding pocket, bent to the direction of TM5/6, although TM3 did not change between the holo- and apo-MDs (Fig. 5a–c). This movement suggests that the channel for molecular oxygen, the hydrophobic cavity surrounded between TM2/4 and TM5/6, is constricted. The holo-MD demonstrated that the distance between Glu242 and Ile66, both of which face the oxygen channel and form a minimum cross-section part, became narrower in the presence of the inhibitor. In contrast, the distance between Ile312 and Val287, which face an exit path for the produced water, did not change (Fig. 5c, d). Furthermore, we performed additional simulations with ligand removal after ligand-bound MD (Re-apo MD). Re-apo MD demonstrated that TM2 relaxed to the initial apo structure, although the change in the oxygen channel did not relax during the Re-apo MD, suggesting that the effect on TM2 of the inhibitor that is next to TM2 is more direct than on the oxygen channel. These data suggested that T113 might interfere with oxygen access to the BNC, thereby inhibiting enzyme reaction.

Fig. 5: Allosteric obstruction of oxygen access to the BNC as a common mechanism of inhibition.
figure 5

a Axes of TM2 and 3 in a representative snapshot from a trajectory of apo-MD, holo-MD, and Re-apo MD. The axes were calculated by UCSF Chimera tools and represented by cylinders. Heme a and T113 molecules were shown as green sticks. b Distribution histograms of the angle for TM2 and 3. c The focused region in the MD simulation. Paths of the oxygen channel (red) and the water channel (green) were calculated by CAVER and shown as continued spheres. TM1-6 in subunit I is shown as dark blue, and the other helices of subunit I as pale blue, subunit COX7C in mtCcO is shown as red, and the other helices as yellow. d Distribution histograms of the distance between Glu242-Ile66 in the oxygen channel and the distance between Ile312-Val287 in the water channel. Arrows show the movement of TM2 by T113 in (a, c). Dot lines in b and d show the median value. e The 441.6-nm excited resonance Raman spectra of mtCcO hemes ((i), (iii)) without and ((ii), (iv)) with N62. mtCcO sample was laser irradiated twice with 10 min dark incubation in between. Spectra ((1); blue) and ((2); green) were obtained from 0 to 3 and 27–30 min in the first 30-min irradiation. After 10 min dark incubation, spectra ((3) black) and ((4) red) were obtained from 0 to 3 min and 27–30 min in the 2nd 30-min irradiation. The laser power was 1 mW for ((i), (iii)) and 0.1 mW for ((ii), (iv)). f A stopped-flow experiment showed that N62 inhibited the CO binding to CcO five times slower than the control. Absorbance changes at 440 nm are associated with the binding of CO to the reduced mtCcO. g Absorbance changes at 430 nm associated with binding of CO to the reduced wild-type bo3 UqO. Fitting lines are calculated by a one-phase decay model. The range of data is presented as a standard deviation of four replicates for CcO, and three replicates for bo3 UqO. Source data are available as a Source Data file.

Full size image

To test whether T113 has an allosteric effect in the oxygen channel, we performed resonance Raman spectroscopy, a sensitive method for detecting structural changes that cannot be assessed by X-ray crystal structural analysis. As T113 has autofluorescence, we screened its derivatives and chose N62 because of its higher affinity without autofluorescence (Supplementary Fig. 4b, c). We confirmed that N62 bound the same binding site in co-crystallography (Supplementary Fig. 4d) and gave a small difference around 440 nm in the heme absorption similar to that with T113 (Supplementary Fig. 4e). We used N62 for Raman spectroscopic analyses hereafter. Previously, it was reported that visible light induces photoreduction of mtCcO, where heme a is initially reduced, followed by heme a337. Figure 4e depicts the resonance Raman spectra of mtCcO with and without N62, focusing on the photoreduction of the hemes. The resonance Raman band at 1356/1372 cm−1 is assignable to the ν4 mode of the hemes (heme a and heme a3), an indicator of the redox status: 1356 cm−1 for the reduced state, 1372 cm−1 for the oxidized state38. In control conditions without N62, laser irradiation at a laser power of 1 mW demonstrated photoreduction of hemes. When the irradiation was stopped, available oxygen caused restoration of the oxidized hemes, reversing the redox marker bands. Reirradiation with the laser demonstrated a reduction of the hemes again (Fig. 5e(i)). By contrast, mtCcO with N62 displayed a different photoreduction behavior. With N62, a comparable reduction of the hemes was observed with a low laser power of only 0.1 mW. Notably, discontinuation of laser irradiation did not decrease the reduction marker band and reirradiation exhibited a further increase of it (Fig. 5e(ii)). These data indicate that N62 inhibits reoxidation of the photoreduced hemes by oxygen. To probe the oxygen binding more directly, we next analyzed the Fe-His stretching mode (νFe-His) of heme a3, the oxygen binding site, at 215 cm−1, representing the 5-coordinated, reduced state of heme a3. The 215 cm−1 bands will disappear when molecular oxygen binds to heme a3. mtCcO without N62 exhibited an increase/decrease cycle of the 215 cm−1 band (Fig. 5e(iii)); however, mtCcO with N62 demonstrated a continuous rise in the 215 cm−1 band (Fig. 5e(iv)), which provides an experimental evidence that N62 inhibits the binding of oxygen to heme a3. Inhibition of reoxidation of the heme a3 by oxygen lowers the threshold of the photoreduction, indicating that the structural change in mtCcO crystal was presumably caused by X-ray-induced reduction during the data acquisition.

To further strengthen the inhibitor’s effect on the oxygen channel, we performed a stopped-flow experiment in which we could directly assess the access of carbon monoxide (alternative to molecular oxygen) to the binuclear center39. As shown in Fig. 5F, N62 inhibited the CO binding to CcO five times slower than the control. Therefore, we concluded that allosteric inhibition in the oxygen-accessing channel plays a major role in their inhibition of mtCcO.

It is plausible that the allostery we found in mtCcO is also preserved between bacterial HCOs and their inhibitors. To this end, we created a single amino acid substitution in the amino acids that face the oxygen channel of bo3 UqO, which can be genetically manipulated. Among the mutants made, less bulky substitutions of Glu286 and Phe112, which correspond to Glu242 and Phe67 in mtCcO (Figs. 3c and 5c), reduced the inhibitory effect of N4 (Supplementary Fig. 6c). Further to obtain direct evidence, we performed a stopped-flow experiment with bo3 UqO. As N4 has absorbance around 400–450 nm, we used N62 for this experiment. N62 is a derivative of T113 and a common inhibitor for both mtCcO and bo3 UqO. First, we confirmed that the inhibitory effect of N62 on bo3 UqO was also reduced by the mutants in the oxygen channel as same as N4 (Supplementary Fig. 6d). A stopped-flow experiment with bo3 UqO demonstrated that the inhibitor also slowed the CO binding in bo3 UqO, as it did on mtCcO (Fig. 5g). Collectively with these multimodal analyses, we conclude that our HCO inhibitors allosterically obstruct molecular oxygen/NO entry to the BNC by a conformational change of a transmembrane helix of subunit I, thereby inhibiting HCO function.

Discussion

Drugs with a unique mode of action are critical for expanding our antimicrobial options to overcome AMR, especially in cases where the pathogens, including N. gonorrhoeae, have acquired resistance to all the currently available antibiotics, spreading globally, and would become untreatable40. The respiratory chain has been a prospective target for antibiotics development; allosteric inhibitors are more desirable than orthostatic ones to minimize the risk of side effects, considering the importance of the respiratory chain in life. Here, we identified conserved allostery in HCOs and an additional helix that is only found in eukaryotic mtCcO. This structural difference in the allostery allowed us to isolate the specific inhibitors for two different bacterial HCOs, including an antibiotic against a ceftriaxone-resistant N. gonorrhoeae, which is one of five urgent threats in the 2019 report from the Centers for Disease Control and Prevention3.

This study is a proof of concept, and the compound is still in the early stage of drug development; however, our findings will pave the way for the development of antibiotics with a different mechanism of action. The comparison of the published structures of HCOs revealed that the steric conformation around the low-spin heme is very much conserved in all HCOs from bacteria, yeasts, plants, and mammals (Supplementary Fig. 3)41, suggesting that the allosteric site identified by us is likely to be conserved. Hence, our approach can generate specific inhibitors, potential antibiotics, for each bacterial HCO on-demand with a narrow range of specificity. The development of narrow-spectrum agents is in line with the current requirement for minimizing the effect on the host microbiome and preventing widespread resistance42. Targeting HCOs may not be simple, as pathogenic bacteria often have multiple terminal oxidases in their respiratory chain. Therefore, a desirable HCO inhibitor as an antibiotic is the one that specifically targets a particular infectious stage where pathogens critically require the HCO to adapt to the growing environment43, as we have shown qNOR as a therapeutic target for N. gonorrhoeae in this study, or can be used with other antibiotics as a combination therapy. Our results suggest that the effects of UqO inhibitor on E. coli and qNOR inhibitor on N. gonorrhoeae were both bacteriostatic. Although bactericidal action sounds preferable, the superiority of bactericidal action over bacteriostatic has rarely been documented44. Further research and development of HCO inhibitors are necessary to reach the clinical arena.

One of the four helices in mtCcO surrounding the allosteric pocket is the genome-encoded subunit COX7C, covering the surface as if it hid the allosteric site. Our phylogenic analysis indicated that the allosteric site of the ancient HCOs was exposed, and it has been sealed in eukaryotic, mitochondrial HCOs during molecular evolution (Supplementary Fig. 7). The eukaryotic-specific subunit might protect mtCcO from access to inhibitors, or there might be an endogenous inhibitor that negatively regulates mtCcO activity at a specific time. The evolutional role of the acquired subunit requires further research.

Mechanisms of inhibition

Combined with stopped-flow experiments, resonance Raman spectroscopy, and mutant analysis, we concluded that our HCO inhibitors allosterically obstruct molecular oxygen/NO entry to the BNC by a conformational change of a transmembrane helix of subunit I, thereby inhibiting HCO function. However, notably in the case of bo3 UqO, the site N4 binds is the quinol binding site. The quinone-bo3 UqO structure confirms that N4 occupies the space where the substrate binds45. Asp75 and Arg71 are the same amino acids N4 used for molecular interaction as quinone does45, suggesting that N4 inhibits bo3 UqO by obstructing the substrate binding. These findings suggest that the space surrounded by TM0 and TM1–3 of bo3 UqO works as both a substrate binding site and an allosteric inhibition site we proposed. TM0 is only present in quinol oxidases, including bo3 UqO. It effectively stabilizes hydrophobic ubiquinol in the transmembrane region so that bo3 UqO can use as a substrate; the existence of TM0 makes the ubiquinol oxidase family unique. TM0 is not found in other types of HCOs, in which the allosteric site we proposed is distinct from the substrate binding site, as shown in qNOR which does not have TM0.

Regarding mtCcO, we demonstrated that the inhibitory mechanism on oxygen entry plays a significant role in T113 and N62 inhibition for mtCcO; however other inhibitory mechanisms might also be involved, as in the case we discussed for N4 on bo3 UqO. We found the structural change in Asp50/51 and Ser382 in the mtCcO-T113 crystal structure. We reasoned that it was caused by photoreduction during sample preparation and radiation-induced reduction for the following reasons. (1) change in Asp50/51 and Ser382 is found in the reduced form of CcO structure; however, the simple addition of T113 did not cause the reduction of CcO (Supplementary Fig. 6a, b). (2) laser irradiation during Raman data acquisition caused CcO reduction. Figure 5e suggests that the mtCcO inhibitor lowered the threshold of photoreduction. (3) MD simulation with the inhibitor did not cause the structural change in Asp50/51 (Supplementary Fig. 6e). These observations, however, did not eliminate the possibility that the binding of T113 induces the structural change found in Asp50/51 and Ser382. The change in these residues might affect proton pumping or electron transfer; Asp50/51 and Ser382 are essential residues for proton pumping, especially in mtCcO, forming the H channel as suggested by Yoshikawa et al.16, although H channel is only found in mtCcO. Rich and colleagues indicated that H channel works as dielectric well17. Furthermore, the Sharma group recently reported that conformational change in the Ser382 carrying domain affects electron transfer46. Thus, perturbation in the region may cause inhibition of proton pumping or electron transfer. Further study is warranted.

Our MD simulation of 100 ns gave us an important clue to the allosteric mechanism of the mtCcO inhibitor; however, Re-apo MD did not show that the oxygen channel relaxes apo. Longer MD may be needed to clarify molecular mechanistic details.

Future perspective for finding allosteric modulators

Our approach can be applied to finding allosteric modulators in other therapeutic targets. Enzymes generally acquire additional domains or subunits along molecular evolution22,23. As the respiratory chain in bioenergetics is fundamental and essential for life, respiratory enzymes other than HCOs are also conserved among species, and their core structures too. The size of the molecule of these respiratory enzymes is bigger in eukaryotes than in their bacterial counterparts. They could likely contain allostery inside the protein at the boundary of the structures between eukaryotes and bacteria, leading to the development of antibiotics, as the respiratory chain is a proven target for antibiotics. Furthermore, any fundamental molecule essential for life and conserved among species could be a potential target. Also, the additional peptides might contain a positive allosteric site at the border of their core structure; a positive allosteric modulator for the loss-of-function human disease could be a therapeutic direction.

Generally, the search for an allosteric site is challenging, requiring considerable experimental work for each target protein and difficult to apply to others. The concept of buried conserved allostery will help develop a systematic approach, which has been critically desired. The number of protein structures has been considerably increased by the emergence of cryo-EM and recent advances in structural prediction such as Alphafold2 and RoseTTAfold47,48. Comparing protein structures among species, not only static structures but also ensembles generated by MD will accelerate finding buried conserved allosteric sites. Thus, in conclusion, this study will open fresh avenues in protein science and therapeutic development, especially for antibiotics with different mechanisms of action.

Methods

Preparation of resting oxidized bovine heart mtCcO crystals

Bovine heart cytochrome c oxidase solution and crystals were prepared as described in a previous study49. Before freezing, the crystals were treated with 2 mM compound or DMSO in the final medium.

X-ray diffraction experiment and structural determination

X-ray experiments were carried out at SPring-8 beamlines BL26B1/B2. The diffraction data processing and scaling were carried out using XDS50. The initial phase was calculated by MOLREP51 using a model of PDB code 5B1A after removing nonprotein molecules. To improve the density and remove the model bias, the maximum entropy method in Phenix52 was performed. For the calculation of differential electron density map between with and without compound, Fo–Fo map calculation in CNS53 was used with a ligand-omitted model. Rebuilding was performed by using COOT54. Models were refined with REFMAC555 and phenix.refine56 with the atomic parameters revised by Dr. Tomitake Tsukihara. To remove the model bias, all electron density maps used in rebuilding were calculated with the compound-omitted model. The electron density maps were confirmed with replicate samples. Refinement statics are provided in Table S1.

Enzyme assay of mtCcO

Purified mtCcO and compounds were incubated on 96-well plates at 30 °C for 30 min in assay buffer (pH 7.4 50 mM potassium phosphate, 0.1% bovine serum albumin, 0.025% 14:0 Lyso PG). To start the enzymatic reaction, reduced cytochrome c (1.5 µM) was added to the mixture, then the plates were read at 550 nm using a plate reader. The slope for 60 s was calculated and the change in mtCcO activity was exported by comparing it to the DMSO-treated sample as a negative control.

Absorption and resonance Raman spectroscopy

mtCcO solutions for spectroscopy were prepared to 15 µM mtCcO with and without 50 µM N62 in 60 mM sodium phosphate buffer and 0.2% n-decyl-β-d-maltoside (pH 6.8) and their spectra were measured at 4 °C unless otherwise noted. Absorption spectra were measured in a 3 mm path cuvette with a spectrophotometer (lamda650, PerkinElmer). The spectra were presented with 1-cm equivalent absorbance values. Oxidized mtCcO was measured under air and dithionite-reduced mtCcO was measured under N2. Resonance Raman spectra were measured in a spinning cell at 2000 rpm with the excitation wavelength at 441.6 nm (a HeCd laser, Kimmons). Raman scattered light was detected by a Raman spectrometer (500M, SPEX) attached with CCD (Symphony II, Horiba) with a 90° scattering geometry. Resting oxidized mtCcO solution for Raman measurements was prepared and transferred into the spinning cell sealed with a rubber septum in a glove box (MM3-H60S, MIWA) at an oxygen level below 100 ppm. The used buffer and rubber septum were degassed and incubated in the glove box overnight. The sample solution thus prepared contained only a small amount of residual oxygen, which ensured the observation of a five-coordinate, reduced heme a3 fraction in the photoreduction measurement, while avoiding immediate full oxidation of all heme a3.

Molecular dynamics simulations

The model system was constructed from the reduced bovine mtCcO crystal structure (PDB ID code 3AG2) to perform all-atom MD simulations. All crystal water molecules present in the PDB ID code 3AG2 atomic structure were removed. The missing hydrogen atoms of the systems were added with the tleap program in AmberTools 2020 (ref: https://ambermd.org/CiteAmber.php). The Amber ff14SB force field for protein, water, and ions57, and Amber lipid14 force field for lipids58 were adopted in the MD simulations. The entire 13-subunit monomeric mtCcO was immersed in a lipid bilayer formed by 50% phosphatidylcholine and 50% phosphatidylethanolamine using CHARMM-GUI59. The protein with lipid bilayers was then solvated with the TIP3P water model in rectangular boxes. The boxes were prepared with margins of 20 Å apart from each protein to avoid artificial effects by its periodically repeated images. The parameters for metal centers and amino acids were obtained from an earlier study60. Briefly, CuA and heme a were oxidized and heme a3 and CuB were reduced. The protonation states of residues were determined by the electrostatic continuum method. The protonation states of the important residues are as follows: all the propiones of hemes a and a3 were deprotonated, Arg438 and 439 were de-protonated, Asp364 and Glu242 were de-protonated, His290 and His291 were protonated in the δ-position. To neutralize each system with counterions, some water molecules were randomly selected and replaced with a set of Na+ or Cl− ions. The particle-mesh Ewald (PME) method was employed for Coulomb interactions61. All the MD simulations were performed using the GROMACS 2019 software62. The entire membrane–protein–solvent system consisted of 302,004 atoms for apo-MD, and 302,029 atoms for holo-MD. The MD simulations were performed under constant-NPT (300 K and 1 atm) ensemble using 2 fs time step.

The following treatments (a–c) were considered to equilibrate each system: (a) Energy minimization was performed for 10,000 steps. (b) A 100 ps MD equilibration was performed with the V-rescale thermostat63 at 300 K under 1 atm, and (c) a 100 ps MD equilibration was performed with the Berendsen coupling64 at 300 K under 1 atm. Each snapshot was recorded every 1 ps65. Finally, a 100 ns MD simulation was performed under constant-NPT (300 K and 1 atm) ensemble as production runs and the last 80 ns were used as analyses. To obtain statistically reliable trajectories, we performed three runs with different initial velocities for the apo or holo-structures.

Preparation of a custom compound library

For the first group, we calculated the Tanimoto coefficient between four mtCcO inhibitors and each of the 80 million commercially available compounds on the basis of 2D fingerprints (MACCS keys, ECFP4, FCFP4, and GpiDAPH3) and 3D shape metrics (ComboScore) with the software Pipeline Pilot (BIOVIA, Dassault Systèmes), MOE (Chemical Computing Group), and ROCS (OpenEye Scientific Software). Compounds showing high Tanimoto coefficient values were clustered by ECFP4 fingerprint, and then compounds in the center of each cluster were selected for enzyme assay. For the second group, mtCcO/T113 complex structure was prepared appropriately using a preparation wizard in Schrödinger suite 2016-1 (Schrödinger, LLC), and then a receptor grid generation was performed for the docking protocol. 80 million commercially available compounds were subjected to molecular docking against mtCcO using Glide 7.0. Compounds selected based on the Glide score, which shares the same space as T113 in the complex, were clustered by ECFP4 fingerprint, and then compounds in the center of each cluster were selected for enzyme assay. Finally, 434 compounds were selected and purchased.

Preparation of bo 3 UqO

The plasmid coding bo3 UqO with a carboxyl-terminus His-tag on subunit II was a gift from Gennis laboratory66. bo3 UqO was purified based on the method as described previously67. Briefly, the plasmid was transformed into the C43 (DE3) Δcyo E. coli strain, and cells were grown in an M63 medium. bo3 UqO was solubilized by 1% n-dodecyl-β-d-maltoside (C12M) and purified with TALON resin, Super-Q resin, and size-exclusion chromatography, and then stored in a buffer containing 25 mM Tris–HCl, 200 mM NaCl, and 0.02% C12M (pH 7.5) at a concentration of 10–15 mg/mL.

Establishment, screening, and preparation of monoclonal antibody

The purified E. coli bo3 UqO was reconstituted into liposomes (egg yolk phosphatidylcholine: E. coli polar lipid = 3:1). Female BALB/c mice were immunized five times with 0.1 mg doses of the reconstituted bo3 UqO with 50 µg of lipopolysaccharide, at intervals of 7 days. Single-cell suspensions were prepared from the spleens of the immunized mice, and the cells were fused with P3U1 myeloma cells, using the conventional polyethylene glycol (PEG) method68. Screening of antibodies was performed by three methods, liposome enzyme-linked immunosorbent assay (L-ELISA), denatured ELISA (D-ELISA), and size0exclusion chromatography (SEC)69. For L-ELISA, the purified bo3 UqO was reconstituted into liposomes containing biotinyl-PE, and was immobilized on Immobilizer Streptavidin plates (Nunc). High-affinity antibodies that formed stable complexes with the purified bo3 UqO were selected by SEC on a Superdex 200 5/150 column (GE Healthcare). Antibodies that recognized the native conformation of E. coli bo3 UqO were assayed by D-ELISA with SDS-denatured E. coli bo3 UqO. Three selected clones were isolated by the limiting dilution-culture method, and monoclonal hybridoma cell lines producing anti-bo3 UqO antibodies were established. The Fab fragment of an anti-bo3 UqO monoclonal antibody was prepared as described69. The purified bo3 UqO proteins were mixed with the Fab fragment for 24 h at 4 °C, and the bo3 UqO–Fab complexes were purified by SEC, followed by incubation with N4. We tested three clones by cryo-EM and selected the one for structural determination.

Enzyme assay of bo 3 UqO

Purified UqO and compounds were incubated on 96-well plates at 25 °C for 60 min in assay buffer (pH 7.4, 50 mM potassium phosphate, 1 mM ethylenediamine-N’,N’,N’,N’-tetraacetic acid (EDTA), 0.1% bovine serum albumin, 0.15% C12M). To start the enzymatic reaction, a reduced 40 µM of UQ-1 was added to the mixture, then the plates were read at 278 nm using a plate reader. The slope was calculated and the enzymatic activity change was exported by comparing it to the DMSO-treated sample as a negative control.

Cryo-EM data acquisition

A 3 µL sample was applied to a glow-discharged Quantifoil grid (R1.2/R1.3 300 mesh, copper), blotted in 100% humidity at 8 °C, and plunged into liquid ethane using Vitrobot MkIV. The Cryo-EM analysis was initially performed by using a Talos Arctica microscope at 200 kV with a Falcon3 detector at KEK. Data collection was performed by using a Glacios microscope with a Gatan K2 Summit detector in the counting mode at RIKEN SPring-8. Movies were acquired at ×45,000 magnification with an accumulated dose of 50.0 electrons per Å2 over 40 frames. The pixel size was 0.889 Å. The data were automatically acquired by the beam-image shift method using the SerialEM software70.

Cryo-EM data processing and model building

Cryo-EM data processing was performed using RELION 3.171. Raw movie stacks were motion-corrected using MotionCor272 of RELION’s own implementation. The CTF parameters were determined using the CTFFIND4 program73. Data processing workflows for UqO (12,388 movie stacks from Glacios) and UqO/compound complex (7173 movie stacks from Glacios) are summarized in Fig. S3A–H, respectively. The final resolution was estimated by gold-standard Fourier shell correlation (FSC) between the two independently refined half maps (FSC  =  0.143). For model building against the EM map, an initial model was generated by homology modeling with the MOE program suite (MOLSIS Inc.) and docked into the final map using Chimera74. The model was then manually refined using the COOT program54. Real-space refinement in the Phenix program was performed with the atomic parameters of hemes revised from the default52. A summary of the model parameters and the cryo-EM map statistics are provided in Table S2.

E. coli growth inhibition assay

E. coli was grown overnight in LB medium and diluted to OD600 of 0.0132. 5 µL of cells was added to each well of the 96-well plate and treated by a compound at concentrations of 128 µg/mL or two-fold serial dilutions in 100 µL LB medium. Plates were then incubated at 37 °C for 24 h without shaking and read at 600 nm using a plate reader to quantify the cell growth. C43 (DE3) stain and C43 (DE3) Δcyo were provided by Dr. Yoshio Nakatani and used as wild-type strain and bo3 UqO-knockout strain, respectively. C43 (DE3) ΔcydΔappx was provided by Dr. Robert Gennis and used for bo3 UqO-dependent strain. As N4 at 32 µg/ml demonstrated growth inhibition from the above assay, we conducted a time- and concentration-dependent viability assay at 32, 128 µg/ml of N4. An inoculum was made of around 1 × 105 CFU/ml and then incubated at 37 °C without shaking. The surviving colonies were counted at fixed times (12, 24, 48, and 72 h).

Preparation of bb 3 qNOR

Neisseria meningitidis bb3 qNOR was provided by Shiro laboratory. Expression and purification were performed as previously reported with modification34.

N. gonorrhoeae culture and minimum inhibitory concentration (MIC)

The N. gonorrhoeae reference strain WHO F75 and a ceftriaxone-resistant FC42876 were used for antimicrobial susceptibility testing. N. gonorrhoeae was cultured in GW medium under 5% CO2 at 37 °C without shaking77. We outsourced the medium preparation to Research Institute for the Functional Peptides Co., Ltd. The micro broth dilution MIC method was used to quantitatively measure the in vitro antibacterial activity of Q275 against the bacterial strains. The lowest concentration of antibiotic that prevented the growth was interpreted as the MIC. Sodium nitrite (20 mM) was added to the medium to mimic a NO-challenging condition (attack from immune cells). For the colony count assay, we prepared an inoculum of around 5 × 105 CFU/ml and then incubated it at 37 °C for 24 h without shaking. The surviving colonies were counted.

Preparation of norB-deficient N. gonorrhoeae

Neisseria gonorrhoeae mutants were constructed as described previously78. Briefly, in order to construct the N. gonorrhoeae norB-deficient mutant, a 4.3-kb DNA fragment from N. gonorrhoeae FA1090 chromosomal DNA containing the norB gene was amplified with the primers norB-1 and norB-2 by PrimeSTAR Max DNA polymerase (Takara Bio) and cloned into the SmaI site of the pMW119 vector (4.2 kb) (Nippon gene) to construct pHT1729 (8.5 kb). The 6.5-kb DNA region of pHT1729 was amplified with the primers norB-3 and norB-4 by PrimeSTAR Max DNA polymerase and was ligated with a kanamycin-resistance gene (kan) to construct pHT1739. A 3-kb DNA fragment containing the norB allele, in which norB structural gene was replaced with kan gene, was amplified with the primers norB-1 and norB-2 from pHT1739, and transformed into N. gonorrhoeae strains and kanamycin-resistant clones were selected, resulting in norB-deficient mutants. Primers used are,

norB-1TCGAGCTCGGTACCCGATGTAGAACTCTTTATCCACTTTCGGCAG

norB-2CTCTAGAGGATCCCCAGGCGGGCAGCCGCCGTTTCCAACGGTTTG

norB-3GGGAAAACCCTGGCGGTTTTAGCCTGAAAATGGAAACCG

norB-4CATAGCTGTTTCCTGTTTGAGAGCTCCTTTTAATAAATC

Stopped-flow spectroscopy

To prepare the fully reduced enzyme, an enzyme solution (1 µM of mtCcO or 0.5 µM bo3 UqO), mixed with DMSO or 200 µM N62, was alternately incubated under vacuumed pressure and N2 atmosphere, then reduced by dithionite. A CO solution was prepared by alternately incubating under vacuumed pressure and a 20% CO atmosphere for mtCcO or 100% CO for bo3 UqO. To investigate the binding of CO to the reduced enzyme, an equal volume of the reduced enzyme solution and the CO solution were mixed in a stopped-flow spectrometer at 4 °C.

Oxygen consumption rate

An XFe96 Extracellular Flux Analyzer (Agilent) was used to measure the oxygen consumption rate (OCR). Mouse C2C12 cells were seeded in an XFe96 culture microplate (1.0 × 104 cells per well) a day before measurement. The OCR of each well was measured in unbuffered DMEM containing 2 mM l-glutamine, and 1 mM sodium pyruvate under basal conditions and in response to 1 µM oligomycin A (Oligo A), 2 µM fluorocarbonyl cyanide phenylhydrazone (FCCP), and 0.5 µM rotenone + 0.5 µM antimycin A (Anti/Rote). OCR for ATP production was calculated as basal OCR−Oligo A OCR using Mito Stress Test Report Generator.

Cytotoxicity

A microplate LDH assay (Cytotoxicity LDH Assay kit, Dojindo, CK12) was used to determine the cytotoxicity of Q275. The mammalian cells used were rat neonatal cardiomyocytes, cultured in DMEM supplemented with 10% fetal bovine serum. Cells were seeded into a 96-well, flat bottom, tissue culture-treated plate and incubated at 37 °C with 5% CO2. After 24 hours, the medium was aspirated and replaced with fresh medium containing test compounds. After 24 h of incubation, LDH activity in the supernatant was assessed according to the manufacturers’ guidance.

Phylogenic analysis

To acquire the amino acid sequence from the database, the sequence whose crystal structure was elucidated and registered in the PDB was preferentially selected (PDB ID was shown). For other amino acid sequences, representative ones in the taxon were selected from HCO enzyme types A, B, C, and NOR according to the previous report41. These indicated as follows: Bos taurus (5B1A_1), Homo sapiens (5Z62_1), Saccharomyces cerevisiae (6T15_11), Pseudomonas stutzeri (3MK7_1), Rhodobacter capsulatus (6XKX_4), Thermus thermophilus (A2: 2YEV_1, ba3: 1EHK_1), Escherichia coli (1FFT_1), Bacillus subtilis (bo3: 6KOB_1, aa3: P24010), Mus musculus (BBA31677), Aquifex (A1: O67935, A2: O67937), Rhodothermus marinus (CAC08532), Paracoccus denitrificans (P08305), Yersinia pestis (WP042593520), Pseudomonas aeruginosa (MXH36788), Snodgrassella alvi (WP_100091491), Rickettsiales bacterium (MBB19300), Lacrimispora indolis (VDG67555), Halobacterium (WP_010902450), Natronomonas pharaonis (CAA71525), Geobacillus stearothermophilus (3AYG), and Pseudomonas aeruginosa (3O0R).

Multiple amino acid sequence alignment was performed using Clustal Omega v1.2.1 after trimming the N-terminus to align the overall length of the alignment79. A default parameter was used, and the conservation of regions essential for the proton pathway (proton channel) in HCOs, such as metal ligands and the catalytic tyrosine residues, was visually confirmed. Neighbor-joining (NJ) trees were inferred with the BioNJ algorithm80 using PhyML v2016011581 with the following parameters: 1000 bootstraps, 100 seeds, and correction for multiple substitutions. A majority consensus tree was created and visualized in Dendroscope v3.7582.

Statistics and reproducibility

There were no statistics applied except for the analysis of SD for Fig. 5f, g. Experiments were repeated at least two times and confirmed its reproducibility.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.