BACKGROUND
numerous reports have shown that enzyme promiscuity, defined as the capacity of an enzyme to perform activities other than the function for which they have evolved–using the same active site– is an extremely common event  <cit> .
the biological implications of these ‘secondary’ activities have been broadly discussed at different levels.
the redundancy of enzymatic functions has been hypothesized to lead to a ‘plastic’ metabolic network, important for organismal evolution  <cit> .
at the protein level, these secondary activities have been hypothesized to serve as raw material for the evolution of new activities  <cit> .
thus, promiscuity could represent an evolutionary advantage that could be selected for as part of a mechanism to acquire novel enzyme functions  <cit> .
indeed, not only positive selection but also neutral evolution, which leads to accumulation of non-conserved mutations usually away from catalytic active sites, have been proposed to lead to promiscuous enzymes  <cit> .
methods for functional classification of protein sequence data based on molecular evolution theory assume that members of protein families will diverge from the consensus sequence as functional shifts take place.
however, the sequence differences associated with such functional shifts, mainly at early evolutionary stages when enzyme promiscuity and sign epistasis occur , are difficult to detect.
moreover, the contribution of point mutations may vary broadly both in terms of the derived trade-off between promiscuous and primary activities  <cit>  and their closeness to the active site of enzymes  <cit> .
therefore, our inability to predict enzyme promiscuity, whose evolutionary nature remains to be fully understood, hampers current enzyme classification systems.
the abovementioned conundrum relates to the issue of ‘hidden’ or adjacent information in biology, a recurring theme in many fields with heterogeneous large datasets.
such heterogeneity has been probabilistically tackled using hidden markov models   <cit> .
hmm profiles, a broadly used method for sequence classification and protein functional annotation , can be defined as a motif definition with some probabilities involved.
indeed, a recent study partly based on the usage of hmm profiles shows how these tools can yield important information on the origin of enzymatic functions  <cit> .
here we hypothesized that hmm profiles provide an efficient approach to identify subtle functional shifts involving enzyme promiscuity.
to test this we focus on the 8-isomerase hisa or n’--5-aminoimidazole-4-carboxamide ribonucleotide  isomerase  enzyme family .
hisa is involved in l-histidine biosynthesis, and due to its internal symmetry, it has been largely used to test hypotheses related to the evolution of the 8-fold from smaller  subunits .
moreover, the hisa homologs from the class actinobacteria have been renamed as pria, from phosphoribosyl isomerase a  <cit> , due to their broader specificity.
indeed, in addition to their ancestral hisa activity, pria enzymes also have a proficient n’- anthranilate  isomerase or trpf activity  .
thus, where a trpf gene is lacking, as in the actinobacteria or gram positive eubacteria with high  content, pria participates in the biosynthesis of both l-histidine and l-tryptophan.fig.
<dig> enzymatic activities of the hisa / pria enzyme superfamily.
hisa is a ubiquitous enzyme that catalyzes the conversion of n’--5-aminoimidazole-4-carboxamide ribonucleotide  into n’--5-aminoimidazole-4-carboxamide ribonucleotide .
pria, exclusively found in the actinobacteridae class of the high g + c gram positive bacteria, evolved from hisa to catalyze, in addition to the hisa activity, the conversion of phosphoribosyl anthranilate  into 1--1-deoxyribulose 5-phosphate .
this latter activity is the same as that catalyzed by the ubiquitous trpf enzyme.
subhisa has recently evolved from pria , and it is found exclusively in certain, but not all, corynebacterium species.
as with hisa, subhisa is mono-functional enzyme
in addition to pria, we have recently reported a closely related group of homologs specific to a certain sub-clade of the genus corynebacterium, which was renamed as subhisa  <cit> .
the subhisa enzyme sub-family evolved from pria in a process resembling sub-functionalization, and it lacks any trpf activity, rendering a mono-substrate profar specific hisa-like enzyme.
the hisa, pria and subhisa enzyme sub-families share 35 % sequence identity between them, and several of their catalytic amino acid residues have been identified .
notably, at least one residue that is conserved at the sub-family level, gly <dig>  ser <dig> and thr <dig>  respectively, was identified early on and further characterized  <cit> .
more recently, additional pria  <cit>  and subhisa  <cit>  specific residues, including arg <dig> and asn <dig>  respectively, have been discovered as sub-family specific residues.
these residues are part of the n-terminal phosphate-binding site  motif, where substrate specificity takes place, in equivalent positions within the 8-fold adopted by these enzymes.
despite these well-established sequence, structural, and functional differences, semi-automatic annotation methods based in multiple sequence alignments  and structural similarities, such as those used by the scop and cath databases, continue to group some hisa, pria and subhisa sequences within one single enzyme family, and under the same enzyme classification  number.
for instance, at time of submission of this paper, cath groups all these enzymes under the enzyme superfamily,  <dig> .
<dig>  or aldolase class i; and the pria enzyme from mycobacterium tuberculosis, despite being comprehensively characterized  <cit> , is still annotated as a hisa enzyme, i.e.
only under the ec number  <dig> .
<dig> , and not by both this number and  <dig> .
<dig> , which refers to trpf or pra isomerase activity.
in this paper we report that with the use of hmm profiles, and in accordance with their substrate specificities, it is possible to successfully classify sequences of unknown function as hisa or pria.
we identify evolutionary intermediates between hisa and pria enzyme sub-families within a large metagenome.
selected enzymes, including single amino-acid mutants constructed following evolutionary analyses, as well as a  <dig> truncated enzyme, which may be useful to explore hypotheses related to the evolution of the 8-fold, were synthesized and biochemically characterized.
our results demonstrate that substrate specificity within enzyme sub-families, and therefore enzyme promiscuity, can be identified with hmm profiles.
methods
construction of pria hidden markov model profiles
the sequence alignment and modeling system program   <cit>  was used to predict the remote intermediates.
hmm profiles were built with the w <dig>  script  <cit> .
to score the training set of sequences with the hmm the hmmscore program was used.
alignments of the hits were performed with the program align2model  <cit> , and editions and redundancy analyses were performed with belvu alignment viewer  <cit> .
the best ten hits were aligned to the hmm as above, and the hits of the original entries were removed from the initial training set to avoid re-sampling in the following iteration.
the resulting multiple sequence alignment  was edited, trimmed to the aligning region, and made non redundant at 80 % identity.
construction of uniprot database and profile validation
a uniprot database  consisting of all sequences in the size range of  <dig> to  <dig> residues was created.
from each scoring run, the sequences with the top  <dig> e-values were retrieved along with their scores and e-values for further analysis.
the data of all rounds were summarized in a matrix describing the e-values of each sequence obtained per iteration.
therefore, the entry mij describes the e-value.
the ith is the sequence obtained when evaluated with the profile hmm from the jth iteration.
construction of metagenomic database
to generate an initial scoring set, sequences were retrieved from the camera database  <cit>  using the blast wizard tool in the webpage http://camera.calit <dig> net/ .
pria and subhisa sequences from actinobacteria, and hisa sequences from other eubacteria and archeae, were blasted against all metagenomic orf peptides  with default parameters.
pria query was p <dig> from streptomyces coelicolor , plus the enzymes annotated as hisa  p9wmm <dig> from mycobacterium tuberculosis , q8g4s <dig> from bifidobacterium longum, q4jw <dig> from corynebacterium jeikenium , q0rfx <dig> from frankia alni, q5yyp <dig> from nocardia farcinica, a8lx <dig> from salinispora arenicola, a6wcu <dig> from kinecoccocus radiotolerans, a1r <dig> from arthrobacter aurescens, a7bd <dig> from actinomyces odontolyticus, and a5csk <dig> from clavibacter michiganensis.
one subhisa sequence, o <dig> from corynebacterium glutamicum, was also used for blast.
hisa queries included one archea: p <dig> from thermus thermofilus and six other eubacteria, two proteobacteria: p <dig> from escherichia coli and q7n8d <dig> from photorhabdus luminiscens; one cyanobacteria: b0c <dig> from acarychloris marina; and two firmicutes: o <dig> from bacillus subtilis and q2rgw <dig> from moorella thermoacetica.
the top  <dig> hits from each run were gathered and repeated sequences above 80 % sequence identity were removed to produce a non-redundant sequence set.
phylogenetic analysis of evolutionary intermediates
all msa were built using muscle within the software seaview  <cit> .
in order to define the best-fitting model for our data, the program prottest  <cit>  was used.
the output of this program was used as input for mrbayes, which was run for two million generations  <cit> .
the sequences used, alignment and tree files were deposited in the treebase repository   <cit> .
construction of truncated enzyme variants
a library of truncated enzyme variants was constructed by systematically removing one by one the amino acids located at the c-terminus end of the protein, covering a total region of  <dig> residues.
the deletion process was performed by the application of a codon-based oligonucleotide method that we developed, referred to as ‘truncagenesis’.
briefly, the method relies on synthesizing a priming segment of the target oligonucleotide and removing a small fraction of the support, carrying the growing oligonucleotide, every three nucleotides during the synthesis process.
the removed fractions are accumulated in a second synthesis column where they are finished with addition of an antisense codon and an appropriate restriction site for cloning purposes, when truncation is targeted at the c-terminus of a protein.
due to technical limitations in the synthesis of long oligonucleotides, coverage of the  <dig> amino acid targeted region was accomplished by the synthesis of five sets of truncated primers, with each set containing  <dig> truncated primers.
sequences of the five sets and the primers contained in set <dig> are shown in additional file 1: table s <dig>  truncated variants were amplified by pcr using the forward primer ndefw and the corresponding reverse set.
the pcr products, as well as the cloning vector, were double digested with the restriction enzymes ndei and hindiii, purified, and ligated.
approximately  <dig>  colonies were obtained, which represents library coverage of 186x.
plasmids from  <dig> colonies were isolated and sequenced, validating that the explored region was randomly shortened, as expected.
cloning, site-directed mutagenesis and functional analysis of hisa homologs
sequences retrieved from camera, termed cam <dig> and cam <dig>  and hisa from acidimicrobium ferrooxidans genes, were commercially synthesized .
codons were optimized for e. coli heterologous expression.
pria from b. longum was cloned from genomic dna generously provided by frabizio arigoni .
the full sequence is shown in additional file 1: table s <dig>  all genes, including cam-derived truncated library, were cloned into pask plasmid  <cit>  using the ndei and hindiii restriction sites.
the pask plasmid is a version of the commercially available vector pask_iba3plus , from which the ndei site was deleted by site-directed mutagenesis using the oligonucleotides paskndefor  and pasknderev .
we have previously used pask for complementation assays involving enzymes with very low promiscuous activities.
in addition, a his-tag from vector pet15b was cloned using the ecori and hindiii restriction sites.
the pask final version contains an ndei site, following the n-terminal his-tag and lacks the c-terminal strep-tag.
cam mutants’ a81s and a81g were constructed using a site-directed mutagenesis commercial kit  and pask_cam <dig> as template.
sequences of oligonucleotides used are: a81s_f  and a81s_reverse  a81g_forward  and a81_r .
the hisa enzyme from a. ferrooxidans was cloned into the pet22b  vector using ndei and hindiii restriction sites, and pria from b. longum, cam1_a81g and a81s variants were cloned in the pet28a  vector using the ndei and xhoi restriction sites.
in vivo activity complementation assays, as well as in vitro steady-state enzyme kinetics, were performed as previously  <cit> , other than pask derivatives were used and m <dig> minimal media was supplemented with anhydrous tetracycline at 20 ng/ml .
construction of 3d structural models
after failed attempts to obtain x-ray crystallographic structures of the proteins investigated, homology models using 2vep structure from streptomyces coelicolor  <cit>  as the initial template were built for cam <dig> , cam <dig>  and cam1_ <dig> .
ten thousand decoys were produced with rosetta  <dig>   <cit>  and clustered for each independent target.
the model with the lowest energy of the widest and most populated cluster was chosen as the final structure for each sequence.
protons were placed afterwards, such that an optimal hydrogen-bonding network was achieved, using the whatif package  <cit> .
two additional independent models were produced in silico, following the same methodology as before, in which the last  <dig> residues of the cam <dig> sequence were replaced with a randomly generated sequence  <cit> .
RESULTS
construction of hmm profiles for classification of hisa, pria and subhisa sequences
our main goal was to obtain a hmm profile that could be used to classify the hisa/pria enzyme family into the three known functional groups, each representing one enzyme sub-family: hisa , pria  or subhisa .
moreover, in addition to classifying enzyme sub-families, the use of an hmm profile would allow detection of evolutionary intermediates between hisa and pria enzymes, where divergence is far too large to allow pinpointing key differences at the sequence level.
in contrast, the differences between pria and subhisa sub-families, which show less sequence divergence despite a big functional leap, facilitate the analysis of this process by simpler sequence comparisons .
the functional assignment for each enzyme sequence family, i.e.
hisa, pria and subhisa, was done accordingly to the amino acid residues previously identified as specific for each sub-family, as describe in the background section and references therein .
the seed sequence used for the initial hmm profile was that of the well-characterized pria enzyme from s. coelicolor, referred to as pria_scoe  <cit> .
hisa sequences, from bacillus subtilis and relatives, were also used for construction of hmm profiles.
probably due to the larger heterogeneity of hisa enzymes in the databases, results obtained using pria gave better resolution early on during these analyses .
thus, the pria hmm profile was further developed, and subsequent pria models incorporated new sequences in each iteration, with redundancy being eliminated at 80 % in each step.
for instance, “q5yyp5” from nocardia farcinica and “a1t8w5” from mycobacterium vanbaaleni, which were present in the model of iteration  <dig>  were removed, as their similarity is above this threshold.
the sequences that were included at each iteration, as well as those that were removed are shown in table  <dig>  the top ten hits obtained after all iterations were aligned to the model as part of the msa that was used to build the resulting hmm.
therefore, every iteration involved:  construction of a hmm profile, which was used to score a training set;  subtraction of the top scoring sequences from the training set;  alignment of these sequences with the hmm and the previous sequences; and  editing of the msa to build the model for the next iteration.
this process yielded one hmm per iteration until it was stopped at the seventh iteration, as no further actinobacterial sequences were retrieved.table  <dig> sequences included at each hmm profile
uniprot entry1
streptomyces coelicolor
thermobifida fusca
kinecococcus radiotolerans
acidothermus cellulolyticus
nocardia farcinica
mycobacterium vanbaaleni
saccharopolyspora erythrea
mycobacterium abscessus
mycobacterium leprae
arthrobacter
sp.
rhodococcus erythroplis
salinispora arenicola
nocardioides
sp.
corynebacterium jeikeium
frankia
sp.
marine actinobacteria
micrococcus luteus
actinomyces naeslundi
clavibacter michiganensis
janibacter
sp.
leifsonia xyli
bifidobacterium longum
propionibacterium acnes
actynomyces odontolyticus
photorhabdus luminescens
desulforudis audaxviator
moorella thermoacetica
pelotomaculum thermoprpionicum
desulfotomaculum reducens
heliobacterium modesticaldum
geobacillus thermodentrificans
bacillus subtilis
1accession number into uniprot.
the pria homologs from the class actinobacteria are highlighted in bold.
hisa homologues from firmicutes, proteobacteria and cyanobacteria are not shadowed.
the sequences included in each hmm profile are shown as x. iter#  refers to the number of iteration made to construct the hmm.
rf stands for refined model obtained after iteration 6
a total of seven hmm profiles were obtained and evaluated.
optimal profiles were defined as pria sequences with low variance between the assigned e-values and with significantly lower e-values for pria than for hisa.
based on these criteria, profiles  <dig> and  <dig> were the only ones that could differentiate between the hisa and pria sub-families.
for the sake of clarity, only results from profiles  <dig> and  <dig> are shown in fig.
<dig>  however, none of the profiles could separate subhisa  from the pria group .
in addition, using profiles  <dig> and  <dig>  some pria sequences coming from deep-rooted actinobacteria, as well as some hisa sequences  coming from organisms belonging to the firmicutes or gram-positive with low  content bacteria , could not be classified.
the region that includes these outlying sequences represents the transition between these enzyme groups, labeled as transition zone  in fig.
<dig> fig.
<dig> performance of hmm profiles against uniprot and camera databases.
hisa, pria, subhisa and ‘transition zone’ sequences are shown as red, blue, orange and green dots, respectively.
uniprot codes corresponding to sequences removed to construct the hmm refined profile are indicated.
the functionally characterized evolutionary intermediates from camera database are marked as cam <dig> and cam2
to improve the hmm profile, a refined profile, constructed from profile  <dig>  was obtained.
this was done after arbitrarily removing the outlying sequences “a4j708” and “q2rgw1” of the firmicutes desulfotomaculum reducens and moorella thermoacetica, as well as “q7n8d1” and “a7bd07” from the enteric bacterium photorhabdus luminescens and the actinobacterium actinomyces odontolyticus, respectively .
these outlying sequences, in particular the latter, may provide interesting candidates for future functional characterization with very unique evolutionary histories, but for the purpose of hmm construction they were removed, as previously advised  <cit> .
this step leads to a more homogeneous sequence sub-sample and increases the resolution of the new model.
indeed, this refined pria model no longer misclassifies firmicutes sequences as pria enzymes, and instead correctly classifies them as hisa enzymes.
use of the refined profile allowed us to classify uniprot hisa and pria sequences, from which  <dig> selected non-redundant sequences were retrieved.
these sequences were successfully grouped, with a difference of nineteen orders of magnitude, into hisa  and pria .
in addition, the pria sequences were slightly separated into two sub-groups.
the first sub-group corresponds to dual-substrate pria enzymes from actinobacterial organisms , while the second group corresponds to subhisa enzymes from the previously identified sub-clade of the genus corynebacterium .
the distribution of these enzymes was 50 % hisa,  <dig>  % pria,  <dig>  % subhisa.
however, only  <dig>  % was assigned to the tz sub-set .
identification of metagenomic evolutionary intermediates between hisa and pria
once our hmm refined profile was validated using uniprot, it was used to search for promiscuous evolutionary intermediates in the large and diverse metagenomic database camera.
at the time when these sequence searches were done, this database had been shown to contain higher diversity, useful for sequence bioprospecting  <cit> .
to define the e-value threshold between pria and subhisa groups, the highest e-value of the pria group and the smallest e-value of the subhisa group obtained from uniprot were selected.
the former was subtracted from the latter and the result was divided by two, giving a value that was added to the highest e-value found for pria.
this is equivalent to subtracting the same value from the lowest e-value found in hisa.
the result of this calculation gives the threshold e-value for classification to functional groups.
thus, sequences with lower e-values were considered pria, and sequences with higher e-values were considered subhisa.
the same procedure was applied to calculate the e-value threshold between subhisa or hisa and the functionally ambiguous enzymes present in the tz.
use of the refined profile allowed us to classify the previously built metagenomic database, consisting of  <dig> non-redundant sequences.
these sequences were grouped into the three known sub-families with a distribution of  <dig>  % hisa,  <dig>  % pria and  <dig>  % subhisa.
moreover, the tz group now included up to 40 % of the sequences, increasing the likelihood of finding evolutionary intermediates.
unlike the e-value distribution of the uniprot sequences, the distribution of e-values of camera sequences shows a continuum that covers the e-value distance between hisa and pria, i.e.
the entire tz .
as previously stated, and also shown in fig.
<dig> and additional file 1: figure s <dig>  hisa and pria sequences have a glycine or a serine residue, respectively, in position  <dig>  however, none of the presumed subhisa sequences  has the threonine in this position, which is the distinctive feature of this sub-family  <cit> .
this may be due to the fact that camera is based on marine samples, and species belonging to the genus corynebacterium may not be present, at least not prominently, in marine environments.
this observation makes it unlikely that these sequences are in fact subhisa sequences, and may well be divergent pria homologs from uncultivated distantly related taxa.
lacking of further data providing insights into the nature of these sequences, we found it unwise to further characterize these enzymes functionally.fig.
<dig> phylogenetic analysis of hisa / pria evolutionary intermediates.
a bayesian phylogenic tree of hisa homologs is shown.
minimum branch support is  <dig> .
hisa, pria and intermediate sequences, as cam <dig>  are grouped in different clades, labeled in red, blue and green, respectively.
the key enzymes hisa_afer  and pria_blon  are also shown.
the evolutionary intermediate, reflected at the codon level in the n-pbs motif, is shown.
functionally analyzed proteins, and their first n-pbs amino acid codon usage, are marked with an asterisk
the tz group includes an additional  <dig> sequences , which have either glycine or serine in position  <dig>  making them closer to either hisa or pria, respectively.
interestingly, exceptions to this observation were found in two sequences that have an alanine in this position, which we called cam <dig> and cam <dig>  alanine is one of the possible transitional states from glycine  to serine  at the codon level.
the mutational path ggc  → gcc  → ucc  modifies the n-terminal pbs, from glyglygly to alaglygly and finally into serglygly .
as the sequences of cam <dig> and cam <dig> suggest that they are evolutionary intermediates, these putative enzymes were investigated further by means of sequence, structural and functional analyses.
at the sequence level, cam <dig> and cam <dig> were found to be almost identical, other than their c-terminal ends.
as shown in the msa of additional file 1: figure s <dig>  cam <dig> has  <dig> amino acids, which is the average length of both hisa and pria enzymes.
in contrast, cam <dig> has only  <dig> amino acids, of which all except the last  <dig> residues are identical to cam <dig> .
thus, a phylogenetic analysis using only the complete sequence of cam <dig>  and not cam <dig>  was performed.
this analysis also included selected sequences obtained from public databases, as retrieved at time of submission, using standard sequence similarity searches.
the resulting phylogenetic tree is shown in fig.
<dig>
the hisa clade  encompasses hisa enzymes from proteobacteria, firmicutes and cyanaobacteria, and it serves as the root.
in the transition between hisa and pria, we found the hisa enzyme from acidimicrobium ferrooxidans , an organism that belongs to the high  content gram-positive bacteria, but a different order than actinobacteria.
moreover, the outgroup of the pria clade  is indeed cam <dig>  as expected, the enzyme “q8g4s5” from bifidobacterium longum , an outgroup of the actinobacteria belonging to the order bifidobacteriales, as well as the actinobacterial deep-rooted genus actinomyces, appear as sister taxa of cam <dig>
functional analysis of evolutionary intermediates
we then decided to experimentally investigate the specificities, towards pra  and profar , of hisa_afer and pria_blon, as well as cam <dig> and the truncated version cam <dig> .
all genes were cloned into suitable plasmids for complementation, protein over-expression and purification procedures.
however, despite our multiple attempts to express and purify cam <dig>  this hypothetical protein could not be overexpressed to a level where purification could be attained.
nevertheless, in vivo complementation of a hisa minus e. coli mutant, termed hfrg <dig>  <cit> , showed that hisa_afer, pria_blon, and the putative evolutionary intermediates cam <dig> and cam <dig>  all have profar isomerase or hisa activity.
however, cells complemented with cam <dig> showed significantly reduced growth.
moreover, based on similar in vivo complementation experiments, but using a trpf minus e. coli mutant termed fbg-wf  <cit> , pra isomerase activity could only be found for pria_blon and cam <dig>  but not for cam <dig> or hisa_afer .table  <dig>
in vivo and in vitro characterization of selected hisa/pria homologs and mutants
in vitro activitya
k
k
cat s−1
k
cat/k
m s−1 μm−1
k
k
cat s−1
k
cat/k
m s−1 μm−1
aeach data point comes from at least three independent determinations using freshly purified enzyme.
n.d., activity not detected
the aforementioned in vivo results were confirmed by in vitro steady-state enzyme kinetics, which led to kcat and km parameters for all enzymes analyzed, other than for cam <dig>  which could not be purified.
a decreased profar isomerase activity for hisa_afer was detected, which may be due to optimal ph  <dig>  growth conditions of the organism from which this enzyme was obtained  <cit> .
no differences in activity, however, could be detected when the enzyme assays were done at different ph values, ranging from  <dig>  till  <dig>  .
it was also noted that cam <dig> has a decreased pra isomerase activity, consistent with the proposed nature of a promiscuous evolutionary intermediate .
indeed, inspection of the kinetic parameters of all these enzymes shows that cam <dig> has a very poor km value for pra , at least one order of magnitude higher than that found for pria_blon .
this observation implies that cam1’s low trpf or pra isomerase activity relates to its inability to properly bind the substrate pra.fig.
<dig> comparison of steady-state enzyme kinetics.
a graphical comparison of the catalytic efficiencies  of wild type and mutant enzymes that were characterized is shown.
values for profar  and pra  isomerase activities, expressed as log <dig>  are compared.
data from hisa_afer , pria_blon , cam <dig> , cam1_a81s  and cam1_a81g , as well as from pria_scoe  and the ser81thr mutant of pria from s. coelicolor, labeled as pria_scoe  , is included
altogether, these results imply  the existence of hisa enzymes able to bind pra promiscuously, or early pria enzymes with dual-substrate specificity, in the root of actinobacteria;  that cam <dig> has the potential for being a hisa to pria evolutionary intermediate, as predicted by our hmm refined profile and suggested by its nascent pra isomerase activity; and  despite cam <dig> having profar isomerase activity in vivo, as suggested by its truncated nature and the difficulties encountered during its characterization, this protein may be an artifact derived from sequence misannotation, and not a naturally occurring 6-barrel .
mutagenesis and structural analysis of cam <dig> and cam <dig> variants
serine  <dig>  located at the n-terminal pbs motif of pria, is part of the active site and its specific role in the binding of pra , and not profar , has been suggested by enzyme kinetics and x-ray crystallographic structural analysis  <cit> .
cam <dig>  as previously mentioned, has an alanine in this position.
site-directed mutagenesis of this residue into glycine  and serine  was done.
as can be observed in table  <dig> and fig.
<dig>  the mutants cam1_a81g and cam1_a81s retain both activities.
cam1_a81g presents similar profar and pra isomerase activities to those found for cam <dig>  whereas cam1_a81s activities were only reduced 8-fold, i.e.
not even an order of magnitude.
this observation shows that this position in cam <dig> can be mutated without major consequences upon catalytic efficiency, which contrasts with the dramatic effect that mutation of ser <dig> in pria, even into a highly similar amino acid such as threonine, has upon its isomerase activities .
crystallographic structural characterization was not possible; therefore, we built homology-based structural predictions for cam <dig> and cam <dig>  as expected, cam <dig> is predicted to fold as a complete  <dig> barrel with two unstructured tails of ~ <dig> amino acids at the n and c terminal regions.
in contrast, cam <dig> is predicted to lack the last two units of the  <dig> barrel, leading to a putative 6-barrel.
the model of cam <dig> shows two unstructured tails, that of the n-terminal region identical to cam <dig>  and a  <dig> amino acid c-terminal tail different from cam <dig> .
given the different enzymatic activities found in cam <dig> and cam <dig>  and starting from the sequence of cam <dig>  we designed a library of truncated variants.
in this library, cam <dig> was truncated by one amino acid at a time, until reaching  <dig> amino acids, with the aim of isolating 6-barrels with activity.fig.
<dig> structural analysis of cam <dig> and cam <dig>  a. structural homology model of cam <dig>  corresponding to a full 8-barrel, and cam <dig>  the  <dig> subunit missing in cam <dig> is shown in gray.
n and c-terminal tails are shown in magenta, loops, alpha helix and beta sheet are marked in green, red and yellow, respectively.
b. structural superposition of cam <dig> , cam2_ <dig>  cam2_r <dig> and cam2_r <dig>  where n and c-terminal tails are shown in blue, red, green and orange, respectively.
c. sequence alignment of the 12-amino acid c-terminal variable region of cam <dig>  cam2_ <dig> , cam2_r <dig> and cam2_r <dig>  the pairwise rmsd between cam <dig> and cam <dig>  cam <dig> and cam2_r <dig>  cam <dig> and cam2_r <dig>  and cam <dig> and cam2_ <dig> was  <dig>  Å,  <dig>  Å,  <dig>  Å,  <dig>  Å, respectively
a library of  <dig>  clones was obtained by pcr assembly and cloning into pask vector, as described in the materials and methods.
the cam <dig> truncated library was screened for profar isomerase activity by complementation using the e. coli hisa− strain hfrg <dig>  colonies that could grow in minimal medium were pooled together and used for screening of pra isomerase activity in the e. coli trpf− strain fbg-wf, as before.
approximately  <dig> positive clones were isolated and their sequences obtained and analyzed.
the results show that only cam <dig> variants with lengths between  <dig> and  <dig> amino acids could be recovered.
according to our structural model, the last  <dig> amino acids correspond to a tail with no secondary structure .
thus, all cam <dig> truncated variants that are active, most likely fold as a complete  <dig> barrels.
we sequenced cam <dig> variants immediately after screening for the profar isomerase activity.
not even variants with the same length to cam <dig>  but with a c-terminal region identical to cam <dig>  or even slighter longer , had profar isomerase activity.
moreover, deletion of the non-conserved c-terminal end of cam <dig>  rendering a protein with  <dig> amino acids , abolished the incipient profar isomerase activity found in cam <dig> .
all together, these analyses support the idea that cam <dig> is an artifact, probably coming from misannotation, and not a naturally occurring protein.
nevertheless, cam <dig> could provide an interesting protein model to investigate the evolution of the  <dig> fold from smaller subunits.
to further investigate cam <dig>  we modeled the structure of cam2_ <dig>  as can be observed in fig.
5b, the structure of the c-terminal regions of both cam <dig> and cams_ <dig> have high similarity.
although these observations are based on computational models, it was interesting to find differences between these two enzymes at the n-terminal region, which corresponds to the first  <dig> amino acids of this tail.
cam <dig> can condense a small alpha helix at this end, while cam <dig> and cam2_ <dig> are predicted to be unable to do so.
to test if the 12-residue sequence found at the c-terminal of cam <dig> is responsible for the structural variation at its n-terminal region, two additional independent structural models were produced .
the last  <dig> residues of the sequence of cam <dig> were replaced with a randomly generated sequence.
the first variant had an amino acid distribution that mimicked that of cam <dig>  while the second variant had an amino acid distribution that mimicked that observed in nature.
these models were termed cam2_r <dig> and cam2_r <dig>  respectively, and no alpha helix condensation could be found at their n-terminal tails .
discussion
our report represents an effort to classify enzyme sub-families with different substrate specificities, aiming to predict enzyme promiscuity in non-annotated databases.
for this purpose, we used hmm profiles and  <dig> isomerases involved in l-histidine and l-tryptophan biosynthesis, two ancestral pathways, for proof-of-concept.
the use of hmm profiles allowed us to distinguish between the closely related enzyme sub-families hisa  and pria , but not the evolutionarily more closely related subhisa .
after validation of our pria hmm profile we classified a subset of hisa homologs contained in the camera metagenomic database.
interestingly, discrete groups were no longer obtained; on the contrary, a continuum of e-values was found .
this seems to correlate with a higher degree of diversity contained within non-cultivated microorganisms, as previously widely noted  <cit> .
by means of msa, we detected two evolutionarily intermediate sequences, cam <dig> and cam <dig>  with a unique pbs sequence motif.
the mutational path from glycine to serine requires a two-nucleotide substitution, with alanine being one of the plausible transition amino acids.
although such amino acid transition has been previously observed in the laboratory after directed evolution experiments  <cit> , to the extent of our knowledge, it is rare to find homologous enzymes that enable the study of a natural mutational path.
thus, the existence of such a scenario, represented here by evolution of pria from hisa via cam, support the idea that cam <dig> and cam <dig> can be assumed to be enzymes representing possible evolutionary trajectories.
along these lines, it would be interesting to take into account residues unrelated to the n-pbs but that directly mediate substrate specificity, such as the loop-embedded arg <dig> and trp <dig> involved in the metamorphosis of the active site of pria  <cit> .
dedicated phylogenetic analysis of hisa, cam and pria enzymes, placed cam <dig> as the outgroup of pria sub-family, between the genera acidimicrobium and bifidobacterium.
in accordance with the bacterial tree of life  acidimicrobium species, a group of early divergent gram positive with high  content bacteria, is the outgroup of bifidobacterium.
moreover, the genus bifidobacterium is the outgroup of the order actinomycetales, which includes species where pria has been found as a response to the lack of a trpf gene  <cit> .
in agreement with this, hisa_afer was shown to have specificity towards profar, while pria_blon can accept both pra and profar substrates.
these observations suggest that evolution of pria from hisa could have occurred in the divergence events that gave rise to the gram positive with high  content bacterial group, dated to  <dig> billion years ago  <cit> .
mutation of residue  <dig>  could be tolerated by cam <dig>  which contrasts with the effect of mutating this residue in pria from s. coelicolor   <cit> .
the reason why cam <dig> could withstand mutations in this key catalytic residue is not understood.
it implies, however, an increased robustness, a characteristic that could be intrinsic to evolutionary intermediates, in this case between the hisa and pria sub-families.
given recent reports, it is tempting to speculate that this observation has to do with the role of epistasis during early enzyme evolution, which involves promiscuous enzyme states .
this could be exploited as an attribute for detection and annotation of enzyme promiscuity from sequence data using hmm profiles.
given the nature of neutral mutations, the development of systematic approaches to detect enzyme promiscuity at the sequence level has been hindered  <cit> .
based on molecular signatures obtained by means of graph kernel support vector machines, it has been suggested previously that it is possible to make predictions with an accuracy of 85 %  <cit> .
more recently, random forests have been used to predict enzyme substrate-specific residues, which are intimately related to enzyme promiscuity  <cit> .
however, although these methods claim to predict function at the fourth level of the enzyme classification number, a situation that remains to be experimentally validated, incorrect annotation in current databases, incorrect training of the method, domain fusions, and bias in the knowledge of functional residues, among other factors, are likely to hamper their efficiency.
our report expands on this repertoire of approaches, with the advantage of being computationally simpler, and possibly also more sensitive as it is based in hmms.
two plausible scenarios have been hypothesized for the evolution of the 8-fold.
first, based on the symmetry seen in hisa and its remote homolog hisf, whose substrates/products are also symmetrical, duplication and subsequent fusion of  <dig> half barrels has been put forward  <cit> .
this hypothesis has been supported experimentally by reconstitution of complete  <dig> barrels from two  <dig> halves  <cit> .
second, an ‘asymmetrical’ pathway, which considers two, three or more duplication and fusion events, proposes that the  <dig> barrel could have originated from  <dig> and  <dig> subunits  <cit> .
indeed, it has been shown experimentally that  <dig> barrels are viable, and thus they may also provide intermediaries in the course of evolution and appearance of the  <dig>  <cit> .
isolation of a functional hisa-like  <dig> barrel in this study suggests that the latter ‘non-symmetrical’ hypothesis may also be viable for symmetrical proteins.
unfortunately, despite several efforts, we could not obtain crystallographic structural data for either cam <dig> or cam <dig>  this observation, together with our inability to isolate  <dig> active proteins after screening of a cam <dig> truncated library, suggests that cam <dig> is not a naturally occurring protein, but rather a fictitious protein arising from misannotation of metagenomic dna.
nevertheless, our homology-based structural modeling supports the idea that cam <dig> maintains some activity because its amino acid sequence at the c-terminal end allows a stable secondary structure.
this predicted structural feature may stabilize the  <dig> fold, and thus its active site, analogous to cases in which remote functional mutations have been found  <cit> .
in the absence of the c-terminal pbs, needed to bind the symmetrical profar substrate, an interesting possibility along these lines would be that cam <dig> undergoes dimerization and domain swapping reconstitute its active site.
CONCLUSIONS
as expected, reconstruction of the evolution of pria from hisa, using hmm profiles and phylogenetics, suggest functional shifts involving intragenic epistatic mutations.
these mutations seem to be compatible with a stability-activity tradeoff that allows a broader exploration of sequence space by nature.
indeed, as the evolutionary intermediates identified here behave as promiscuous enzymes, a key evolutionary raw material, intragenic epistasis is confirmed as a mechanism driving functional shifts.
the use of hmm provides a convenient approach for gaining insights into these evolutionary processes, as was shown here for the evolution of pria from hisa, but not for subhisa from pria, as the latter are evolutionarily closely related and thus highly similar.
isolation of a  <dig> protein with hisa activity provides further evidence for the hypothesis that extant  <dig> may have evolved not only from symmetrical subunits.
additional file
additional file 1: figure s <dig>  sequence alignment of functionally characterized proteins.
figure s <dig>
e. coli trpf
−
and hisa
− complementation assays with selected proteins.
figure s <dig>  catalytic activity ph dependence of hisa_afer.
table s <dig>  oligonucleotides used to construct cam <dig> truncated variants.
table s <dig>  protein sequence of synthetic genes.
lianet noda-garcía and ana l. juárez-vázquez contributed equally to this work.
competing interests
the authors declare that they have no competing interests.
authors’ contributions
lng and aljv performed experiments, analyzed results, and helped draft the manuscript.
mcaa and fbg designed and performed the hmm profile analyses.
lng, aljv, and fbg designed and performed the sequence and phylogenetic analyses.
mct, lng and fbg designed the structural modeling analyses, which were specifically performed by mct.
eavc and gmm helped with functional analyses of enzymes.
pg, lng and fbg designed and analyzed the library of truncated proteins, which was specifically constructed by pg.
fbg conceived the study, participated in its design and coordination, and helped to draft the manuscript.
all authors read and approved the final manuscript.

 numerous reports have shown that enzyme promiscuity, defined as the capacity of an enzyme to perform activities other than the function for which they have evolved–using the same active site– is an extremely common event. the biological implications of these ‘secondary’ activities have been broadly discussed at different levels. the redundancy of enzymatic functions has been hypothesized to lead to a ‘plastic’ metabolic network, important for organismal evolution. at the protein level, these secondary activities have been hypothesized to serve as raw material for the evolution of new activities. thus, promiscuity could represent an evolutionary advantage that could be selected for as part of a mechanism to acquire novel enzyme functions. indeed, not only positive selection but also neutral evolution, which leads to accumulation of non-conserved mutations usually away from catalytic active sites, have been proposed to lead to promiscuous enzymes. methods for functional classification of protein sequence data based on molecular evolution theory assume that members of protein families will diverge from the consensus sequence as functional shifts take place. however, the sequence differences associated with such functional shifts, mainly at early evolutionary stages when enzyme promiscuity and sign epistasis occur, are difficult to detect. moreover, the contribution of point mutations may vary broadly both in terms of the derived trade-off between promiscuous and primary activities and their closeness to the active site of enzymes. therefore, our inability to predict enzyme promiscuity, whose evolutionary nature remains to be fully understood, hampers current enzyme classification systems. the abovementioned conundrum relates to the issue of ‘hidden’ or adjacent information in biology, a recurring theme in many fields with heterogeneous large datasets. such heterogeneity has been probabilistically tackled using hidden markov models. hmm profiles, a broadly used method for sequence classification and protein functional annotation, can be defined as a motif definition with some probabilities involved. indeed, a recent study partly based on the usage of hmm profiles shows how these tools can yield important information on the origin of enzymatic functions. here we hypothesized that hmm profiles provide an efficient approach to identify subtle functional shifts involving enzyme promiscuity. to test this we focus on the 8-isomerase hisa or n’--5-aminoimidazole-4-carboxamide ribonucleotide isomerase enzyme family. hisa is involved in l-histidine biosynthesis, and due to its internal symmetry, it has been largely used to test hypotheses related to the evolution of the 8-fold from smaller subunits. moreover, the hisa homologs from the class actinobacteria have been renamed as pria, from phosphoribosyl isomerase a, due to their broader specificity. indeed, in addition to their ancestral hisa activity, pria enzymes also have a proficient n’- anthranilate isomerase or trpf activity. thus, where a trpf gene is lacking, as in the actinobacteria or gram positive eubacteria with high content, pria participates in the biosynthesis of both l-histidine and l-tryptophan.fig. enzymatic activities of the hisa / pria enzyme superfamily. hisa is a ubiquitous enzyme that catalyzes the conversion of n’--5-aminoimidazole-4-carboxamide ribonucleotide into n’--5-aminoimidazole-4-carboxamide ribonucleotide. pria, exclusively found in the actinobacteridae class of the high g + c gram positive bacteria, evolved from hisa to catalyze, in addition to the hisa activity, the conversion of phosphoribosyl anthranilate into 1--1-deoxyribulose 5-phosphate. this latter activity is the same as that catalyzed by the ubiquitous trpf enzyme. subhisa has recently evolved from pria, and it is found exclusively in certain, but not all, corynebacterium species. as with hisa, subhisa is mono-functional enzyme in addition to pria, we have recently reported a closely related group of homologs specific to a certain sub-clade of the genus corynebacterium, which was renamed as subhisa. the subhisa enzyme sub-family evolved from pria in a process resembling sub-functionalization, and it lacks any trpf activity, rendering a mono-substrate profar specific hisa-like enzyme. the hisa, pria and subhisa enzyme sub-families share 35 % sequence identity between them, and several of their catalytic amino acid residues have been identified. notably, at least one residue that is conserved at the sub-family level, gly ser and thr respectively, was identified early on and further characterized. more recently, additional pria and subhisa specific residues, including arg and asn respectively, have been discovered as sub-family specific residues. these residues are part of the n-terminal phosphate-binding site motif, where substrate specificity takes place, in equivalent positions within the 8-fold adopted by these enzymes. despite these well-established sequence, structural, and functional differences, semi-automatic annotation methods based in multiple sequence alignments and structural similarities, such as those used by the scop and cath databases, continue to group some hisa, pria and subhisa sequences within one single enzyme family, and under the same enzyme classification number. for instance, at time of submission of this paper, cath groups all these enzymes under the enzyme superfamily,. or aldolase class i; and the pria enzyme from mycobacterium tuberculosis, despite being comprehensively characterized, is still annotated as a hisa enzyme, i.e. only under the ec number., and not by both this number and., which refers to trpf or pra isomerase activity. in this paper we report that with the use of hmm profiles, and in accordance with their substrate specificities, it is possible to successfully classify sequences of unknown function as hisa or pria. we identify evolutionary intermediates between hisa and pria enzyme sub-families within a large metagenome. selected enzymes, including single amino-acid mutants constructed following evolutionary analyses, as well as a truncated enzyme, which may be useful to explore hypotheses related to the evolution of the 8-fold, were synthesized and biochemically characterized. our results demonstrate that substrate specificity within enzyme sub-families, and therefore enzyme promiscuity, can be identified with hmm profiles. methods construction of pria hidden markov model profiles the sequence alignment and modeling system program was used to predict the remote intermediates. hmm profiles were built with the w script. to score the training set of sequences with the hmm the hmmscore program was used. alignments of the hits were performed with the program align2model, and editions and redundancy analyses were performed with belvu alignment viewer. the best ten hits were aligned to the hmm as above, and the hits of the original entries were removed from the initial training set to avoid re-sampling in the following iteration. the resulting multiple sequence alignment was edited, trimmed to the aligning region, and made non redundant at 80 % identity. construction of uniprot database and profile validation a uniprot database consisting of all sequences in the size range of to residues was created. from each scoring run, the sequences with the top e-values were retrieved along with their scores and e-values for further analysis. the data of all rounds were summarized in a matrix describing the e-values of each sequence obtained per iteration. therefore, the entry mij describes the e-value. the ith is the sequence obtained when evaluated with the profile hmm from the jth iteration. construction of metagenomic database to generate an initial scoring set, sequences were retrieved from the camera database using the blast wizard tool in the webpage http://camera.calit net/. pria and subhisa sequences from actinobacteria, and hisa sequences from other eubacteria and archeae, were blasted against all metagenomic orf peptides with default parameters. pria query was p from streptomyces coelicolor, plus the enzymes annotated as hisa p9wmm from mycobacterium tuberculosis, q8g4s from bifidobacterium longum, q4jw from corynebacterium jeikenium, q0rfx from frankia alni, q5yyp from nocardia farcinica, a8lx from salinispora arenicola, a6wcu from kinecoccocus radiotolerans, a1r from arthrobacter aurescens, a7bd from actinomyces odontolyticus, and a5csk from clavibacter michiganensis. one subhisa sequence, o from corynebacterium glutamicum, was also used for blast. hisa queries included one archea: p from thermus thermofilus and six other eubacteria, two proteobacteria: p from escherichia coli and q7n8d from photorhabdus luminiscens; one cyanobacteria: b0c from acarychloris marina; and two firmicutes: o from bacillus subtilis and q2rgw from moorella thermoacetica. the top hits from each run were gathered and repeated sequences above 80 % sequence identity were removed to produce a non-redundant sequence set. phylogenetic analysis of evolutionary intermediates all msa were built using muscle within the software seaview. in order to define the best-fitting model for our data, the program prottest was used. the output of this program was used as input for mrbayes, which was run for two million generations. the sequences used, alignment and tree files were deposited in the treebase repository. construction of truncated enzyme variants a library of truncated enzyme variants was constructed by systematically removing one by one the amino acids located at the c-terminus end of the protein, covering a total region of residues. the deletion process was performed by the application of a codon-based oligonucleotide method that we developed, referred to as ‘truncagenesis’. briefly, the method relies on synthesizing a priming segment of the target oligonucleotide and removing a small fraction of the support, carrying the growing oligonucleotide, every three nucleotides during the synthesis process. the removed fractions are accumulated in a second synthesis column where they are finished with addition of an antisense codon and an appropriate restriction site for cloning purposes, when truncation is targeted at the c-terminus of a protein. due to technical limitations in the synthesis of long oligonucleotides, coverage of the amino acid targeted region was accomplished by the synthesis of five sets of truncated primers, with each set containing truncated primers. sequences of the five sets and the primers contained in set are shown in additional file 1: table s truncated variants were amplified by pcr using the forward primer ndefw and the corresponding reverse set. the pcr products, as well as the cloning vector, were double digested with the restriction enzymes ndei and hindiii, purified, and ligated. approximately colonies were obtained, which represents library coverage of 186x. plasmids from colonies were isolated and sequenced, validating that the explored region was randomly shortened, as expected. cloning, site-directed mutagenesis and functional analysis of hisa homologs sequences retrieved from camera, termed cam and cam and hisa from acidimicrobium ferrooxidans genes, were commercially synthesized. codons were optimized for e. coli heterologous expression. pria from b. longum was cloned from genomic dna generously provided by frabizio arigoni. the full sequence is shown in additional file 1: table s all genes, including cam-derived truncated library, were cloned into pask plasmid using the ndei and hindiii restriction sites. the pask plasmid is a version of the commercially available vector pask_iba3plus, from which the ndei site was deleted by site-directed mutagenesis using the oligonucleotides paskndefor and pasknderev. we have previously used pask for complementation assays involving enzymes with very low promiscuous activities. in addition, a his-tag from vector pet15b was cloned using the ecori and hindiii restriction sites. the pask final version contains an ndei site, following the n-terminal his-tag and lacks the c-terminal strep-tag. cam mutants’ a81s and a81g were constructed using a site-directed mutagenesis commercial kit and pask_cam as template. sequences of oligonucleotides used are: a81s_f and a81s_reverse a81g_forward and a81_r. the hisa enzyme from a. ferrooxidans was cloned into the pet22b vector using ndei and hindiii restriction sites, and pria from b. longum, cam1_a81g and a81s variants were cloned in the pet28a vector using the ndei and xhoi restriction sites. in vivo activity complementation assays, as well as in vitro steady-state enzyme kinetics, were performed as previously, other than pask derivatives were used and m minimal media was supplemented with anhydrous tetracycline at 20 ng/ml. construction of 3d structural models after failed attempts to obtain x-ray crystallographic structures of the proteins investigated, homology models using 2vep structure from streptomyces coelicolor as the initial template were built for cam, cam and cam1_. ten thousand decoys were produced with rosetta and clustered for each independent target. the model with the lowest energy of the widest and most populated cluster was chosen as the final structure for each sequence. protons were placed afterwards, such that an optimal hydrogen-bonding network was achieved, using the whatif package. two additional independent models were produced in silico, following the same methodology as before, in which the last residues of the cam sequence were replaced with a randomly generated sequence.  construction of hmm profiles for classification of hisa, pria and subhisa sequences our main goal was to obtain a hmm profile that could be used to classify the hisa/pria enzyme family into the three known functional groups, each representing one enzyme sub-family: hisa, pria or subhisa. moreover, in addition to classifying enzyme sub-families, the use of an hmm profile would allow detection of evolutionary intermediates between hisa and pria enzymes, where divergence is far too large to allow pinpointing key differences at the sequence level. in contrast, the differences between pria and subhisa sub-families, which show less sequence divergence despite a big functional leap, facilitate the analysis of this process by simpler sequence comparisons. the functional assignment for each enzyme sequence family, i.e. hisa, pria and subhisa, was done accordingly to the amino acid residues previously identified as specific for each sub-family, as describe in the background section and references therein. the seed sequence used for the initial hmm profile was that of the well-characterized pria enzyme from s. coelicolor, referred to as pria_scoe. hisa sequences, from bacillus subtilis and relatives, were also used for construction of hmm profiles. probably due to the larger heterogeneity of hisa enzymes in the databases, results obtained using pria gave better resolution early on during these analyses. thus, the pria hmm profile was further developed, and subsequent pria models incorporated new sequences in each iteration, with redundancy being eliminated at 80 % in each step. for instance, “q5yyp5” from nocardia farcinica and “a1t8w5” from mycobacterium vanbaaleni, which were present in the model of iteration were removed, as their similarity is above this threshold. the sequences that were included at each iteration, as well as those that were removed are shown in table  the top ten hits obtained after all iterations were aligned to the model as part of the msa that was used to build the resulting hmm. therefore, every iteration involved: construction of a hmm profile, which was used to score a training set; subtraction of the top scoring sequences from the training set; alignment of these sequences with the hmm and the previous sequences; and editing of the msa to build the model for the next iteration. this process yielded one hmm per iteration until it was stopped at the seventh iteration, as no further actinobacterial sequences were retrieved.table sequences included at each hmm profile uniprot entry1 streptomyces coelicolor thermobifida fusca kinecococcus radiotolerans acidothermus cellulolyticus nocardia farcinica mycobacterium vanbaaleni saccharopolyspora erythrea mycobacterium abscessus mycobacterium leprae arthrobacter sp. rhodococcus erythroplis salinispora arenicola nocardioides sp. corynebacterium jeikeium frankia sp. marine actinobacteria micrococcus luteus actinomyces naeslundi clavibacter michiganensis janibacter sp. leifsonia xyli bifidobacterium longum propionibacterium acnes actynomyces odontolyticus photorhabdus luminescens desulforudis audaxviator moorella thermoacetica pelotomaculum thermoprpionicum desulfotomaculum reducens heliobacterium modesticaldum geobacillus thermodentrificans bacillus subtilis 1accession number into uniprot. the pria homologs from the class actinobacteria are highlighted in bold. hisa homologues from firmicutes, proteobacteria and cyanobacteria are not shadowed. the sequences included in each hmm profile are shown as x. iter# refers to the number of iteration made to construct the hmm. rf stands for refined model obtained after iteration 6 a total of seven hmm profiles were obtained and evaluated. optimal profiles were defined as pria sequences with low variance between the assigned e-values and with significantly lower e-values for pria than for hisa. based on these criteria, profiles and were the only ones that could differentiate between the hisa and pria sub-families. for the sake of clarity, only results from profiles and are shown in fig.  however, none of the profiles could separate subhisa from the pria group. in addition, using profiles and some pria sequences coming from deep-rooted actinobacteria, as well as some hisa sequences coming from organisms belonging to the firmicutes or gram-positive with low content bacteria, could not be classified. the region that includes these outlying sequences represents the transition between these enzyme groups, labeled as transition zone in fig.  fig. performance of hmm profiles against uniprot and camera databases. hisa, pria, subhisa and ‘transition zone’ sequences are shown as red, blue, orange and green dots, respectively. uniprot codes corresponding to sequences removed to construct the hmm refined profile are indicated. the functionally characterized evolutionary intermediates from camera database are marked as cam and cam2 to improve the hmm profile, a refined profile, constructed from profile was obtained. this was done after arbitrarily removing the outlying sequences “a4j708” and “q2rgw1” of the firmicutes desulfotomaculum reducens and moorella thermoacetica, as well as “q7n8d1” and “a7bd07” from the enteric bacterium photorhabdus luminescens and the actinobacterium actinomyces odontolyticus, respectively. these outlying sequences, in particular the latter, may provide interesting candidates for future functional characterization with very unique evolutionary histories, but for the purpose of hmm construction they were removed, as previously advised. this step leads to a more homogeneous sequence sub-sample and increases the resolution of the new model. indeed, this refined pria model no longer misclassifies firmicutes sequences as pria enzymes, and instead correctly classifies them as hisa enzymes. use of the refined profile allowed us to classify uniprot hisa and pria sequences, from which selected non-redundant sequences were retrieved. these sequences were successfully grouped, with a difference of nineteen orders of magnitude, into hisa and pria. in addition, the pria sequences were slightly separated into two sub-groups. the first sub-group corresponds to dual-substrate pria enzymes from actinobacterial organisms, while the second group corresponds to subhisa enzymes from the previously identified sub-clade of the genus corynebacterium. the distribution of these enzymes was 50 % hisa,  % pria,  % subhisa. however, only  % was assigned to the tz sub-set. identification of metagenomic evolutionary intermediates between hisa and pria once our hmm refined profile was validated using uniprot, it was used to search for promiscuous evolutionary intermediates in the large and diverse metagenomic database camera. at the time when these sequence searches were done, this database had been shown to contain higher diversity, useful for sequence bioprospecting. to define the e-value threshold between pria and subhisa groups, the highest e-value of the pria group and the smallest e-value of the subhisa group obtained from uniprot were selected. the former was subtracted from the latter and the result was divided by two, giving a value that was added to the highest e-value found for pria. this is equivalent to subtracting the same value from the lowest e-value found in hisa. the result of this calculation gives the threshold e-value for classification to functional groups. thus, sequences with lower e-values were considered pria, and sequences with higher e-values were considered subhisa. the same procedure was applied to calculate the e-value threshold between subhisa or hisa and the functionally ambiguous enzymes present in the tz. use of the refined profile allowed us to classify the previously built metagenomic database, consisting of non-redundant sequences. these sequences were grouped into the three known sub-families with a distribution of  % hisa,  % pria and  % subhisa. moreover, the tz group now included up to 40 % of the sequences, increasing the likelihood of finding evolutionary intermediates. unlike the e-value distribution of the uniprot sequences, the distribution of e-values of camera sequences shows a continuum that covers the e-value distance between hisa and pria, i.e. the entire tz. as previously stated, and also shown in fig.  and additional file 1: figure s hisa and pria sequences have a glycine or a serine residue, respectively, in position however, none of the presumed subhisa sequences has the threonine in this position, which is the distinctive feature of this sub-family. this may be due to the fact that camera is based on marine samples, and species belonging to the genus corynebacterium may not be present, at least not prominently, in marine environments. this observation makes it unlikely that these sequences are in fact subhisa sequences, and may well be divergent pria homologs from uncultivated distantly related taxa. lacking of further data providing insights into the nature of these sequences, we found it unwise to further characterize these enzymes functionally.fig. phylogenetic analysis of hisa / pria evolutionary intermediates. a bayesian phylogenic tree of hisa homologs is shown. minimum branch support is. hisa, pria and intermediate sequences, as cam are grouped in different clades, labeled in red, blue and green, respectively. the key enzymes hisa_afer and pria_blon are also shown. the evolutionary intermediate, reflected at the codon level in the n-pbs motif, is shown. functionally analyzed proteins, and their first n-pbs amino acid codon usage, are marked with an asterisk the tz group includes an additional sequences, which have either glycine or serine in position making them closer to either hisa or pria, respectively. interestingly, exceptions to this observation were found in two sequences that have an alanine in this position, which we called cam and cam alanine is one of the possible transitional states from glycine to serine at the codon level. the mutational path ggc  → gcc  → ucc modifies the n-terminal pbs, from glyglygly to alaglygly and finally into serglygly. as the sequences of cam and cam suggest that they are evolutionary intermediates, these putative enzymes were investigated further by means of sequence, structural and functional analyses. at the sequence level, cam and cam were found to be almost identical, other than their c-terminal ends. as shown in the msa of additional file 1: figure s cam has amino acids, which is the average length of both hisa and pria enzymes. in contrast, cam has only amino acids, of which all except the last residues are identical to cam. thus, a phylogenetic analysis using only the complete sequence of cam and not cam was performed. this analysis also included selected sequences obtained from public databases, as retrieved at time of submission, using standard sequence similarity searches. the resulting phylogenetic tree is shown in fig.  the hisa clade encompasses hisa enzymes from proteobacteria, firmicutes and cyanaobacteria, and it serves as the root. in the transition between hisa and pria, we found the hisa enzyme from acidimicrobium ferrooxidans, an organism that belongs to the high content gram-positive bacteria, but a different order than actinobacteria. moreover, the outgroup of the pria clade is indeed cam as expected, the enzyme “q8g4s5” from bifidobacterium longum, an outgroup of the actinobacteria belonging to the order bifidobacteriales, as well as the actinobacterial deep-rooted genus actinomyces, appear as sister taxa of cam functional analysis of evolutionary intermediates we then decided to experimentally investigate the specificities, towards pra and profar, of hisa_afer and pria_blon, as well as cam and the truncated version cam. all genes were cloned into suitable plasmids for complementation, protein over-expression and purification procedures. however, despite our multiple attempts to express and purify cam this hypothetical protein could not be overexpressed to a level where purification could be attained. nevertheless, in vivo complementation of a hisa minus e. coli mutant, termed hfrg, showed that hisa_afer, pria_blon, and the putative evolutionary intermediates cam and cam all have profar isomerase or hisa activity. however, cells complemented with cam showed significantly reduced growth. moreover, based on similar in vivo complementation experiments, but using a trpf minus e. coli mutant termed fbg-wf, pra isomerase activity could only be found for pria_blon and cam but not for cam or hisa_afer.table in vivo and in vitro characterization of selected hisa/pria homologs and mutants in vitro activitya k k cat s−1 k cat/k m s−1 μm−1 k k cat s−1 k cat/k m s−1 μm−1 aeach data point comes from at least three independent determinations using freshly purified enzyme. n.d., activity not detected the aforementioned in vivo results were confirmed by in vitro steady-state enzyme kinetics, which led to kcat and km parameters for all enzymes analyzed, other than for cam which could not be purified. a decreased profar isomerase activity for hisa_afer was detected, which may be due to optimal ph  growth conditions of the organism from which this enzyme was obtained. no differences in activity, however, could be detected when the enzyme assays were done at different ph values, ranging from till. it was also noted that cam has a decreased pra isomerase activity, consistent with the proposed nature of a promiscuous evolutionary intermediate. indeed, inspection of the kinetic parameters of all these enzymes shows that cam has a very poor km value for pra, at least one order of magnitude higher than that found for pria_blon. this observation implies that cam1’s low trpf or pra isomerase activity relates to its inability to properly bind the substrate pra.fig. comparison of steady-state enzyme kinetics. a graphical comparison of the catalytic efficiencies of wild type and mutant enzymes that were characterized is shown. values for profar and pra isomerase activities, expressed as log are compared. data from hisa_afer, pria_blon, cam, cam1_a81s and cam1_a81g, as well as from pria_scoe and the ser81thr mutant of pria from s. coelicolor, labeled as pria_scoe, is included altogether, these results imply the existence of hisa enzymes able to bind pra promiscuously, or early pria enzymes with dual-substrate specificity, in the root of actinobacteria; that cam has the potential for being a hisa to pria evolutionary intermediate, as predicted by our hmm refined profile and suggested by its nascent pra isomerase activity; and despite cam having profar isomerase activity in vivo, as suggested by its truncated nature and the difficulties encountered during its characterization, this protein may be an artifact derived from sequence misannotation, and not a naturally occurring 6-barrel. mutagenesis and structural analysis of cam and cam variants serine located at the n-terminal pbs motif of pria, is part of the active site and its specific role in the binding of pra, and not profar, has been suggested by enzyme kinetics and x-ray crystallographic structural analysis. cam as previously mentioned, has an alanine in this position. site-directed mutagenesis of this residue into glycine and serine was done. as can be observed in table  and fig.  the mutants cam1_a81g and cam1_a81s retain both activities. cam1_a81g presents similar profar and pra isomerase activities to those found for cam whereas cam1_a81s activities were only reduced 8-fold, i.e. not even an order of magnitude. this observation shows that this position in cam can be mutated without major consequences upon catalytic efficiency, which contrasts with the dramatic effect that mutation of ser in pria, even into a highly similar amino acid such as threonine, has upon its isomerase activities. crystallographic structural characterization was not possible; therefore, we built homology-based structural predictions for cam and cam as expected, cam is predicted to fold as a complete barrel with two unstructured tails of ~ amino acids at the n and c terminal regions. in contrast, cam is predicted to lack the last two units of the barrel, leading to a putative 6-barrel. the model of cam shows two unstructured tails, that of the n-terminal region identical to cam and a amino acid c-terminal tail different from cam. given the different enzymatic activities found in cam and cam and starting from the sequence of cam we designed a library of truncated variants. in this library, cam was truncated by one amino acid at a time, until reaching amino acids, with the aim of isolating 6-barrels with activity.fig. structural analysis of cam and cam a. structural homology model of cam corresponding to a full 8-barrel, and cam the subunit missing in cam is shown in gray. n and c-terminal tails are shown in magenta, loops, alpha helix and beta sheet are marked in green, red and yellow, respectively. b. structural superposition of cam, cam2_ cam2_r and cam2_r where n and c-terminal tails are shown in blue, red, green and orange, respectively. c. sequence alignment of the 12-amino acid c-terminal variable region of cam cam2_, cam2_r and cam2_r the pairwise rmsd between cam and cam cam and cam2_r cam and cam2_r and cam and cam2_ was  Å,  Å,  Å,  Å, respectively a library of clones was obtained by pcr assembly and cloning into pask vector, as described in the materials and methods. the cam truncated library was screened for profar isomerase activity by complementation using the e. coli hisa− strain hfrg colonies that could grow in minimal medium were pooled together and used for screening of pra isomerase activity in the e. coli trpf− strain fbg-wf, as before. approximately positive clones were isolated and their sequences obtained and analyzed. the results show that only cam variants with lengths between and amino acids could be recovered. according to our structural model, the last amino acids correspond to a tail with no secondary structure. thus, all cam truncated variants that are active, most likely fold as a complete barrels. we sequenced cam variants immediately after screening for the profar isomerase activity. not even variants with the same length to cam but with a c-terminal region identical to cam or even slighter longer, had profar isomerase activity. moreover, deletion of the non-conserved c-terminal end of cam rendering a protein with amino acids, abolished the incipient profar isomerase activity found in cam. all together, these analyses support the idea that cam is an artifact, probably coming from misannotation, and not a naturally occurring protein. nevertheless, cam could provide an interesting protein model to investigate the evolution of the fold from smaller subunits. to further investigate cam we modeled the structure of cam2_ as can be observed in fig. 5b, the structure of the c-terminal regions of both cam and cams_ have high similarity. although these observations are based on computational models, it was interesting to find differences between these two enzymes at the n-terminal region, which corresponds to the first amino acids of this tail. cam can condense a small alpha helix at this end, while cam and cam2_ are predicted to be unable to do so. to test if the 12-residue sequence found at the c-terminal of cam is responsible for the structural variation at its n-terminal region, two additional independent structural models were produced. the last residues of the sequence of cam were replaced with a randomly generated sequence. the first variant had an amino acid distribution that mimicked that of cam while the second variant had an amino acid distribution that mimicked that observed in nature. these models were termed cam2_r and cam2_r respectively, and no alpha helix condensation could be found at their n-terminal tails. discussion our report represents an effort to classify enzyme sub-families with different substrate specificities, aiming to predict enzyme promiscuity in non-annotated databases. for this purpose, we used hmm profiles and isomerases involved in l-histidine and l-tryptophan biosynthesis, two ancestral pathways, for proof-of-concept. the use of hmm profiles allowed us to distinguish between the closely related enzyme sub-families hisa and pria, but not the evolutionarily more closely related subhisa. after validation of our pria hmm profile we classified a subset of hisa homologs contained in the camera metagenomic database. interestingly, discrete groups were no longer obtained; on the contrary, a continuum of e-values was found. this seems to correlate with a higher degree of diversity contained within non-cultivated microorganisms, as previously widely noted. by means of msa, we detected two evolutionarily intermediate sequences, cam and cam with a unique pbs sequence motif. the mutational path from glycine to serine requires a two-nucleotide substitution, with alanine being one of the plausible transition amino acids. although such amino acid transition has been previously observed in the laboratory after directed evolution experiments, to the extent of our knowledge, it is rare to find homologous enzymes that enable the study of a natural mutational path. thus, the existence of such a scenario, represented here by evolution of pria from hisa via cam, support the idea that cam and cam can be assumed to be enzymes representing possible evolutionary trajectories. along these lines, it would be interesting to take into account residues unrelated to the n-pbs but that directly mediate substrate specificity, such as the loop-embedded arg and trp involved in the metamorphosis of the active site of pria. dedicated phylogenetic analysis of hisa, cam and pria enzymes, placed cam as the outgroup of pria sub-family, between the genera acidimicrobium and bifidobacterium. in accordance with the bacterial tree of life acidimicrobium species, a group of early divergent gram positive with high content bacteria, is the outgroup of bifidobacterium. moreover, the genus bifidobacterium is the outgroup of the order actinomycetales, which includes species where pria has been found as a response to the lack of a trpf gene. in agreement with this, hisa_afer was shown to have specificity towards profar, while pria_blon can accept both pra and profar substrates. these observations suggest that evolution of pria from hisa could have occurred in the divergence events that gave rise to the gram positive with high content bacterial group, dated to billion years ago. mutation of residue could be tolerated by cam which contrasts with the effect of mutating this residue in pria from s. coelicolor. the reason why cam could withstand mutations in this key catalytic residue is not understood. it implies, however, an increased robustness, a characteristic that could be intrinsic to evolutionary intermediates, in this case between the hisa and pria sub-families. given recent reports, it is tempting to speculate that this observation has to do with the role of epistasis during early enzyme evolution, which involves promiscuous enzyme states. this could be exploited as an attribute for detection and annotation of enzyme promiscuity from sequence data using hmm profiles. given the nature of neutral mutations, the development of systematic approaches to detect enzyme promiscuity at the sequence level has been hindered. based on molecular signatures obtained by means of graph kernel support vector machines, it has been suggested previously that it is possible to make predictions with an accuracy of 85 %. more recently, random forests have been used to predict enzyme substrate-specific residues, which are intimately related to enzyme promiscuity. however, although these methods claim to predict function at the fourth level of the enzyme classification number, a situation that remains to be experimentally validated, incorrect annotation in current databases, incorrect training of the method, domain fusions, and bias in the knowledge of functional residues, among other factors, are likely to hamper their efficiency. our report expands on this repertoire of approaches, with the advantage of being computationally simpler, and possibly also more sensitive as it is based in hmms. two plausible scenarios have been hypothesized for the evolution of the 8-fold. first, based on the symmetry seen in hisa and its remote homolog hisf, whose substrates/products are also symmetrical, duplication and subsequent fusion of half barrels has been put forward. this hypothesis has been supported experimentally by reconstitution of complete barrels from two halves. second, an ‘asymmetrical’ pathway, which considers two, three or more duplication and fusion events, proposes that the barrel could have originated from and subunits. indeed, it has been shown experimentally that barrels are viable, and thus they may also provide intermediaries in the course of evolution and appearance of the. isolation of a functional hisa-like barrel in this study suggests that the latter ‘non-symmetrical’ hypothesis may also be viable for symmetrical proteins. unfortunately, despite several efforts, we could not obtain crystallographic structural data for either cam or cam this observation, together with our inability to isolate active proteins after screening of a cam truncated library, suggests that cam is not a naturally occurring protein, but rather a fictitious protein arising from misannotation of metagenomic dna. nevertheless, our homology-based structural modeling supports the idea that cam maintains some activity because its amino acid sequence at the c-terminal end allows a stable secondary structure. this predicted structural feature may stabilize the fold, and thus its active site, analogous to cases in which remote functional mutations have been found. in the absence of the c-terminal pbs, needed to bind the symmetrical profar substrate, an interesting possibility along these lines would be that cam undergoes dimerization and domain swapping reconstitute its active site.  as expected, reconstruction of the evolution of pria from hisa, using hmm profiles and phylogenetics, suggest functional shifts involving intragenic epistatic mutations. these mutations seem to be compatible with a stability-activity tradeoff that allows a broader exploration of sequence space by nature. indeed, as the evolutionary intermediates identified here behave as promiscuous enzymes, a key evolutionary raw material, intragenic epistasis is confirmed as a mechanism driving functional shifts. the use of hmm provides a convenient approach for gaining insights into these evolutionary processes, as was shown here for the evolution of pria from hisa, but not for subhisa from pria, as the latter are evolutionarily closely related and thus highly similar. isolation of a protein with hisa activity provides further evidence for the hypothesis that extant may have evolved not only from symmetrical subunits. additional file additional file 1: figure s sequence alignment of functionally characterized proteins. figure s e. coli trpf − and hisa − complementation assays with selected proteins. figure s catalytic activity ph dependence of hisa_afer. table s oligonucleotides used to construct cam truncated variants. table s protein sequence of synthetic genes. lianet noda-garcía and ana l. juárez-vázquez contributed equally to this work. competing interests the authors declare that they have no competing interests. authors’ contributions lng and aljv performed experiments, analyzed results, and helped draft the manuscript. mcaa and fbg designed and performed the hmm profile analyses. lng, aljv, and fbg designed and performed the sequence and phylogenetic analyses. mct, lng and fbg designed the structural modeling analyses, which were specifically performed by mct. eavc and gmm helped with functional analyses of enzymes. pg, lng and fbg designed and analyzed the library of truncated proteins, which was specifically constructed by pg. fbg conceived the study, participated in its design and coordination, and helped to draft the manuscript. all authors read and approved the final manuscript. 
these sequences were used to perform phylogenetic reconstructions and to identify functionally conserved amino acids. none of these mutations, presumably upon functionally essential amino acids, significantly abolished its enzyme activities. a truncated version of this enzyme predicted to adopt a 6-fold, and thus entirely lacking a c-terminus phosphate-binding site, was identified and shown to have hisa activity. we constructed a hmm profile that correctly classifies sequences of unknown function into the hisa and pria enzyme sub-families. these results are in agreement with a link between promiscuous enzymes and intragenic epistasis. using this hmm profile, we mined a large metagenome to identify plausible evolutionary intermediate sequences between hisa and pria. as expected, reconstruction of the evolution of pria from hisa with hmm profiles suggest that functional shifts involve mutations in evolutionarily intermediate enzymes of otherwise functionally essential residues or motifs. hmm provides a convenient approach for gaining insights into these evolutionary processes. we searched for functional evolutionary intermediates between the hisa and pria enzyme families in order to understand the functional divergence between these families. biochemical characterization of one selected enzyme with a mutation within the functionally essential n-terminus phosphate-binding site, namely, an alanine instead of a glycine in hisa or a serine in pria, showed that this evolutionary intermediate has dual-substrate specificity. current sequence-based approaches to identify enzyme functional shifts, such as enzyme promiscuity, have proven to be highly dependent on a priori functional knowledge, hampering our ability to reconstruct evolutionary history behind these mechanisms.
 current sequence-based approaches to identify enzyme functional shifts, such as enzyme promiscuity, have proven to be highly dependent on a priori functional knowledge, hampering our ability to reconstruct evolutionary history behind these mechanisms. hidden markov model profiles, broadly used to classify enzyme families, can be useful to distinguish between closely related enzyme families with different specificities. the 8-isomerase hisa/pria enzyme family, involved in l-histidine biosynthesis in most bacteria and plants, but also in l-tryptophan biosynthesis in most actinobacteria, has been used as model system to explore evolutionary hypotheses and therefore has a considerable amount of evolutionary, functional and structural knowledge available. we searched for functional evolutionary intermediates between the hisa and pria enzyme families in order to understand the functional divergence between these families.  we constructed a hmm profile that correctly classifies sequences of unknown function into the hisa and pria enzyme sub-families. using this hmm profile, we mined a large metagenome to identify plausible evolutionary intermediate sequences between hisa and pria. these sequences were used to perform phylogenetic reconstructions and to identify functionally conserved amino acids. biochemical characterization of one selected enzyme with a mutation within the functionally essential n-terminus phosphate-binding site, namely, an alanine instead of a glycine in hisa or a serine in pria, showed that this evolutionary intermediate has dual-substrate specificity. moreover, site-directed mutagenesis of this alanine residue, either backwards into a glycine or forward into a serine, revealed the robustness of this enzyme. none of these mutations, presumably upon functionally essential amino acids, significantly abolished its enzyme activities. a truncated version of this enzyme predicted to adopt a 6-fold, and thus entirely lacking a c-terminus phosphate-binding site, was identified and shown to have hisa activity.  as expected, reconstruction of the evolution of pria from hisa with hmm profiles suggest that functional shifts involve mutations in evolutionarily intermediate enzymes of otherwise functionally essential residues or motifs. these results are in agreement with a link between promiscuous enzymes and intragenic epistasis. hmm provides a convenient approach for gaining insights into these evolutionary processes. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. issue-copyright-statement© the author 2015 