Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Characterizing lineage-specific evolution and the processes driving genomic diversification in chordates

BMC Evolutionary Biology volume 20, Article number: 24 (2020) Cite this article

Abstract

Background

Understanding the origins of genome content has long been a goal of molecular evolution and comparative genomics. By examining genome evolution through the guise of lineage-specific evolution, it is possible to make inferences about the evolutionary events that have given rise to species-specific diversification. Here we characterize the evolutionary trends found in chordate species using The Adaptive Evolution Database (TAED). TAED is a database of phylogenetically indexed gene families designed to detect episodes of directional or diversifying selection across chordates. Gene families within the database have been assessed for lineage-specific estimates of dN/dS and have been reconciled to the chordate species to identify retained duplicates. Gene families have also been mapped to the functional pathways and amino acid changes which occurred on high dN/dS lineages have been mapped to protein structures.

Results

An analysis of this exhaustive database has enabled a characterization of the processes of lineage-specific diversification in chordates. A pathway level enrichment analysis of TAED determined that pathways most commonly found to have elevated rates of evolution included those involved in metabolism, immunity, and cell signaling. An analysis of protein fold presence on proteins, after normalizing for frequency in the database, found common folds such as Rossmann folds, Jelly Roll folds, and TIM barrels were overrepresented on proteins most likely to undergo directional selection. A set of gene families which experience increased numbers of duplications within short evolutionary times are associated with pathways involved in metabolism, olfactory reception, and signaling. An analysis of protein secondary structure indicated more relaxed constraint in β-sheets and stronger constraint on alpha Helices, amidst a general preference for substitutions at exposed sites. Lastly a detailed analysis of the ornithine decarboxylase gene family, a key enzyme in the pathway for polyamine synthesis, revealed lineage-specific evolution along the lineage leading to Cetacea through rapid sequence evolution in a duplicate gene with amino acid substitutions causing active site rearrangement.

Conclusion

Episodes of lineage-specific evolution are frequent throughout chordate species. Both duplication and directional selection have played large roles in the evolution of the phylum. TAED is a powerful tool for facilitating this understanding of lineage-specific evolution.

Background

As closely related species diverge after a speciation event, their genomes begin to accumulate changes that lead to molecular and phenotypic divergence. Speciation itself is a complex process in chordates that results from the gradual cessation of gene flow. As the isolated populations become separate species, mutations of different magnitudes affect the protein coding repertoire of the two diverging genomes. These changes include synonymous changes that only affect the nucleotide sites, nonsynonymous changes that affect the amino acid sites, and gene duplication and loss events, among other types of changes. A resource comparing chordate genomes in a phylogenetic context, The Adaptive Evolution Database (TAED) has recently been re-generated [33] extending previous versions that were released [46, 66].

The latest version of TAED contains gene families constructed systematically across chordate species as described in Hermansen et al. [33]. Gene families have been filtered for alignment quality and to prevent synonymous site saturation, with the oldest nodes in each rooted gene tree reflecting a speciation event of maximum age being the root of the chordate divergence. All pairwise alignments within each multiple sequence alignment had no more than 10% gaps and were at least 80% identical in non-gapped positions. This then created a trade-off between gene family ages (many had root nodes younger than the last common ancestor of chordates) and alignment quality, although homologous gene family relationships can still be identified through TAED. Gene families have been reconciled to the NCBI taxonomy [67] as a reference species tree and events of positive directional and diversifying selection detected using nonsynonymous to synonymous nucleotide substitution rate ratios in the branches model averaged across sites [83]. Gene families have also been used to identify duplication events using the SoftParsMap parsimony-based gene tree-species tree reconciliation software [9].

In addition to previous iterations of TAED, other studies have also sought to characterize the lineage-specific evolution of chordate genomes. This includes the generation of the Selectome Database [51] from Ensembl [2] data. Selectome extends gene family data automatically generated through the Ensembl pipeline which contains sequences from 68 different genomes. Gene families in Selectome are passed through stringent quality control steps following which tests of selection using branch-site models are implemented against tree topologies from Ensembl. While both Ensembl and Selectome examine evolution in a lineage-specific context, the method by which selection is detected varies, with Ensembl using pairwise analyses to calculate the normalized rate of nonsynonymous to synonymous substitutions (dN/dS) and Selectome using branch-site models of selection based on phylogenetic trees. Pairwise estimates of dN/dS do not account for phylogenetic information which limits the ability to understand evolution in a lineage-specific context, and prohibits detection of directional or diversifying selection on internal lineages. Branch-site models and branch models differ in their sensitivity (power) and selectivity (detection of false positives) [5, 25]. dS saturation is a potential problem for these approaches, with accuracy declining at dS ~ 3 [6].

Gene duplication is another important process to consider when assessing lineage-specific processes of evolution. As genes duplicate, they may undergo different evolutionary pressures and be either neofunctionalized, subfunctionalized, or pseudogenized [42]. In the classical model [55], duplicate gene copies can acquire mutations that lose (pseudogenize), change or gain (neofunctionalize) function mutations when the other copy retains the original function. Neofunctionalization, which can also occur to a gene subsequent to initial subfunctionalization, emerges as the dominant driver of evolution in duplicated genes in this model [35, 65]. As such it is one driver of lineage-specific differences in genome content. Subfunctionalization, the subdividing of functions from an ancestral state, can also lead to lineage-specific functional divergence of genes, without the gain of new functions in the genome as a whole. Without gene duplication as a source of genetic content unconstrained by negative selection, evolution tends to act in a conservative fashion [55].

TAED also presents a picture of lineage-specific evolution using pathway and structural information in addition to selection on individual protein encoding genes and gene duplication. Pathway level analyses of proteins may lead to understanding how proteins evolve in the context of a cell or organism, since proteins typically interact together in a pathway or network to achieve biological functions (phenotypes). Simulations have suggested that rate limiting steps are not evolutionarily stable over longer evolutionary periods [56, 57] and proteins currently involved in rate limiting steps may not remain so over long evolutionary periods. This suggests patterns that might be expected for gene-specific selective pressures in a pathway and how they relate to phenotypic evolution.

Two models for the evolution of pathways have been presented, the retrograde evolution model [34], proposing evolution to build a pathway backwards from the selected final product based upon affinity for related transition states at neighboring positions of a pathway and the patchwork model [38] suggesting that gene duplication retains catalytic mechanisms on widely distributed substrates that are dispersed throughout the network of pathways. A driver of mutational opportunity in both models is gene duplication. Analysis of protein function can identify which model is best associated with the evolution of a given pathway, with evidence suggesting that the patchwork model is more common [48]. TAED compiles duplication and selection data compiled for pathways in a lineage-specific manner that can be viewed in this light.

Understanding the structural context of substitutions within a protein may elucidate the role of individual amino acid changes in potential functional shifts under positive selection, differentiating them from compensatory or stabilizing substitutions within the protein. Modeling the effects of amino acid substitutions can demonstrate changes in structure, dynamics, allosteric regulation, and ligand binding that can be used to identify functional shifts ([19]; see also [16]). Such modeling is limited however as the process is difficult and computationally intensive, with identification of fitness effects based upon biophysical models inexact. Measurements and models based on experimental work can also contribute to our understanding [14].

The structural context of mutations also impacts the substitution rate via negative selection. Requirements for folding stability drive lower substitutions in the protein core, while binding requirements on the ligand interface slow mutation as compared to the protein surface [28]. These constraints extend to functional requirements to avoid certain alternate states, including both selection against alternate folding states and substrates that result in deleterious interactions [47]. As protein structure diverges less observably than protein sequence over equivalent units of evolutionary time [36], similar structural constraints can be assumed to be approximately equivalently applicable to sequences diverged over relatively short evolutionary times.

Understanding how genes evolve and the processes by which they lead to novel adaptations in species is fundamental to understanding the genotype-phenotype map. Here we present some new characterizations of lineage-specific evolution utilizing the TAED database; we examine specific hypotheses across lineages, as well as characterizing processes at the levels of gene duplication, pathway evolution, and of protein structure.

Results

The Adaptive Evolution Database (TAED) contains ~ 3.2 million sequences from 3214 different chordate species. The database contains 143,806 individual genes families which are mapped to the chordate species tree. Twenty-three thousand nine hundred seventy gene families contained one or more branches with dN/dS > 1, indicating positive or directional selection acting on these lineages. When the dN/dS rates are high after controlling for dS saturation, the lineages are candidates for having undergone functional shifts. It is expected that the larger the dN/dS value for a given branch, the stronger the putative selective forces were to cause functional changes to the ancestral protein [73]. A list of the lineages with the largest dN/dS values where dS > 0.01 was generated, as these proteins constituted potential strong candidates for having undergone positive selection (Table 1). Of the top 30 lineages with the largest dN/dS values, values were found to range from 88.78 to 26.57. The families that these proteins come from are putatively involved in multiple different biological processes, many of which do not map to a KEGG pathway. Interestingly strong selection was found to have occurred on the branch leading from Boreoeutherian mammals in 9 of the top 30 instances of high dN/dS. This lineage constitutes species before the split of Laurasiatheria and Euarchontoglires, following the divergence of mammals. Additionally, strong selection was seen repeatedly on the lineage leading from Laurasiatheria which is the superorder containing cetaceans, carnivores, chiropterans, and ruminants. Functional shifts in these proteins may be responsible for some of the physiological and habitat differences between these groups and shared ancestors with carnivores and primates. Strong selection was seen to occur on the lineage leading from Neognathae which comprises most avian species. Pathways under selection along this lineage may indicate some of the functional differences between flightless birds which comprise the sister order Palaeognathae and other avians. KEGG pathway mappings for the top 30 lineages with high dN/dS showed that selection may have acted on several different pathway types including metabolic pathway interactions, receptor signaling pathways, and immune response pathways. Selection can act directly on many different levels within an organism. It can occur at the DNA level, the protein level, the pathway level, and the phenotypic level. Understanding pathway evolution may ultimately be a better way to assess selection than current codon based methods [32].

Table 1 TAED gene family lineages with the largest dN/dS values where dS > 0.001
Full size table

Enrichment analysis

To gain a better understanding of pathways within TAED that are more common targets of directional selection, a test to determine which pathways were over or under represented for instances of putative positive selection was undertaken. Table 2 shows the list of the top 25 enriched KEGG pathways within TAED for directional selection. From the top 25 pathways that are over-represented in the database, 8 of the pathways are involved in metabolic reactions (the pathway labeled “Metabolic pathways” contains proteins from all metabolic pathways, and therefore is not a unique pathway). Metabolism, or the process of constructing useful cellular molecules, is essential for life. Given the vast array of different physiological and environmental conditions that exist within chordate species, it is plausible that developing different metabolic strategies is a primary way for organisms to cope with their surroundings. As such, seeing that these pathways are often targets for directional selection is not surprising. Furthermore, it is evident from the list that pathways involved in immune response and cellular health have also been directly impacted by selection. Over-represented pathways involved in immune response included: Herpes simplex infection, Influenza A, Toxoplasmosis, and Th17 cell differentiation. It has been documented in the literature that selection against pathogens is a constant arms race that requires novel adaptations to overcome the constant pressures of pathogenic infection [15, 44, 78]; that these pathways should be over-represented for putative positive selection is not surprising. Additionally, pathways which alleviate physiological stress also appear to be over-represented for directional selection as seen in the pathways: fluid shear stress and atherosclerosis, non-alcoholic fatty liver disease, and chemical carcinogenesis. Cellular components were also found to be under selective pressure to evolve as seen in the pathways, protein processing in endoplasmic reticulum, RNA transport, lysosome, and peroxisome. Lastly, many lineages were found to have evolved under directional selection relating to olfactory transduction. Olfactory genes are the most duplicated genes within the human genome and are known to be largely expanded in other chordate species [54]. Olfactory sense is a primary means of communication, predation, and foraging for many species and thus is unsurprising that many lineages relating to this pathway have instances of dN/dS > 1.

Table 2 Pathways present in lineages under positive selection
Full size table

Of the pathways found within TAED to be under-represented for functional shifts, surprisingly phototransduction was found to be included within the top 25 (Table 3). The ability to visually see pigments is important in both sexual selection and predation. In birds [12, 84], fish ([72, 74, 79];) and cetaceans [24] instances of positive selection have been discovered relating to selection on opsin and rhodopsin genes. Therefore, it is surprising that selection on this KEGG pathway would be under-represented within TAED. However, KEGG pathways for zeatin biosynthesis, penicillin and cephalosporin biosynthesis, bacterial secretion systems, and MAPK signaling pathway – plant, should be underrepresented in the database as these pathways are primarily involved in either plant or microbial systems and do not constitute meaningful pathways in chordates although orthologous proteins to some of the components to these pathways do exist in chordates, but may have different functions. RNA polymerase is a highly-conserved protein found throughout all domains of life, and therefore is unsurprising that the pathway for RNA polymerase would be under-represented for functional shifts within chordate species.

Table 3 Pathways absent in lineages under positive selection
Full size table

Another interesting question which was generated from structural elements contained in TAED was if some functional protein domains are more likely to experience elevated rates of evolution compared to others. To determine if this is true a systematic search was performed to determine what functional domain topologies are enriched within lineages in TAED that have signals for functional change (Table 4). Functional domains were annotated from the CATH database which assigns each domain a CATH classification. Annotations for this analysis looked at the topology level as it contains a wide array of functional domain annotations. The most over-represented domain/fold within TAED was the Rossmann fold which constituted approximately a quarter of all lineages in TAED with dN/dS > 1 that could map to a domain (the analysis did normalize for abundance in the database). The Rossmann fold is a common fold comprised of a b-a-b-a-b (b – beta sheet, a – alpha helix) subunit motif and is commonly found within nucleotide-binding proteins [63]. Proteins that include this fold type include kinases, guanine nucleotide binding proteins (G proteins), proteins that bind cyclic adenosine monophosphate (cAMP), and NAD(P)-binding proteins [31]. These proteins are abundant within a cell and therefore proteins in which these domains reside are likely candidates for directional selection. However due to the nature and importance of nucleotide binding, it is unlikely that the Rossmann fold is under selection, but other domains within the same protein are as this domain is likely under strong negative constraint unless there are selective pressures on binding affinity or specificity. More structural analyses of the lineages under selection that contain the Rossmann fold would be warranted to examine this in more detail. The second most over represented domain topology was the Jelly Rolls fold which a subset of the beta-barrels superfamily. This fold type is composed of 8 beta-sheets which fold into a roll shape [1]. These folds are commonly found in viral capsid proteins [64]. It is possible that since these folds are commonly found in viral proteins that they evolve quickly and are prone to high mutation rates. This would suggest that protein families which contain this domain would be over-represented. The third most over-represented domain topology was TIM barrel folds. These are very common folds found with proteins that share alpha-beta structures. The TIM barrel folds are known to be highly promiscuous in sequence with many different sequences able to generate the TIM barrel fold. Therefore, there is biophysical flexibility for amino acids within these domains to be substituted while still maintaining the same domain structure [82]. These folds are in some cases known over longer evolutionary periods as folds that are structurally adaptable and evolve under relaxed selective constraint [17, 27, 45], consistent with their observation here in divergence among closely related species.

Table 4 Domains present in lineages under positive selection
Full size table

From the list of the top under-represented domain topologies (Table 5), two of the most under-represented domains were derived from the SMAD3 (mothers against decapentaplegic homolog 3) protein (smad3 chain A and Smad anchor for receptor activation chain B). The SMAD3 protein is involved in the signal trafficking of TGF-β which plays an important role in cell growth and death. This protein structure is known to contain two different domains, a DNA-binding domain and a protein-protein interacting domain. These two domains have been shown to be conserved across many species and play an essential role in the function of SMAD proteins [52, 53]. Accordingly, it is expected that these domains would be very limited in the rate at which they evolve and that they would evolve mostly under strong negative selection. Another interesting protein domain that was under-represented within the database was the fold for cAMP-dependent protein kinase. The primary enzyme which contains this domain is protein kinase A (PKA) which is involved in many different cellular pathways and plays a role in cell growth and differentiation, signaling, and migration [21]. As a central hub protein within a protein interaction network, it would be expected that this would be highly negatively constrained [58] and therefore domains that are essential to this protein are also under strong negative selection.

Table 5 Domains absent in lineages under positive selection
Full size table

Duplication analysis

One important element of lineage-specific evolution is the expansion and contraction of genes within the genome. As genes duplicate they may undergo different evolutionary pressures and be either neofunctionalized, subfunctionalized, or pseudogenize [42]. Following the completion of the TAED database, it was interesting to determine if some gene families are more likely to undergo gene duplication events than others and what pathways these genes reside in. Are some pathways more flexible to gene duplication and dosage balance constraints [76] than others? A systematic examination of TAED gene family duplications was performed by scaling the number of duplication events detected within a family by the amount of time over which the family evolved. Three different proxies for time were used in the analysis, the maximum phylogenetic tree length measured in substitutions per site (Additional file 1: Figure S1), the median tree length measured in substitutions per site (Additional file 1: Figure S2), and the relative age of each family found by mapping the root of each gene tree to the chordate species tree (Fig. 1). Each analysis determined that there is a positive correlation between the number of duplications within the family and the amount of time over which the family evolved. Outliers from the regression line identified families that were highly duplicated over a shortened timespan. These families are also those with a high rate of duplication compared to other gene families. Table 6 shows the Cook’s distance calculations for the analysis using family node age as a proxy for time and the corresponding gene families that were calculated to be furthest from the regression line. Cook’s distances for the maximum tree length and median tree length are found in Additional file 1: Tables S1 and S2, respectively. From the families with the largest Cook’s distance the number of times a highly duplicable family mapped to a give KEGG pathways was counted (Table 7). Pathway counts for the maximum tree length and median tree lengths were also calculated (Additional file 1: Tables S3 and S4).

Fig. 1
figure 1

Duplication analysis regression plot using family node ages as a proxy for time – The x-axis is measured in MYA based on the root node for each TAED gene family. The best Pearson’s r coefficient was found when neither axes were log transformed. The upper left half (shaded orange) of the scatterplot was used to determine TAED gene families that were statistically different from the regression line using Cook’s distance

Full size image
Table 6 TAED gene families with many duplications based on family node age from summed branch lengths
Full size table
Table 7 TAED KEGG pathways based on duplication analysis using family node age from summed branch lengths
Full size table

The data shows metabolic pathways and olfactory receptors are consistently the top pathways where duplications occur. Olfactory receptors are known to be the largest expanded gene family [26], aligning our study with the currently known data.

Additionally, the top 25 most highly duplicable gene families included serine/threonine-protein phosphatase 2A 56 kDa regulatory subunit epsilon isoform, abl interactor 1 - partial, aldolase B, guanine nucleotide-binding protein G(i) subunit alpha-1 - partial, and myosin regulatory light polypeptide 9. A further examination of the structural components and pathway components of these families may explain why they are more tolerable to duplication events and the mechanisms that are causing large gene family expansions. Interestingly, many of the most duplicated gene families mapped to KEGG pathways involved in immunity (HTLV-I infection; Herpes simplex infection; Epstein-Barr virus infection; Influenza A) and cancer (Pathways in cancer; Proteoglycans in cancer; Transcriptional misregulation in cancer; Viral carcinogenesis), possibly suggesting that duplication plays a strong role in this arms race.

Protein structure based analysis

The combination of gene families and information from the Protein Databank allows examination of how selection acts on a protein structural level. Gene families with associated protein structures were collated and aligned to the PDB alongside maximum likelihood ancestral sequences calculated by PAML.

The resulting profile is significantly different than the profile of non-substituted sites in the background on those lineages (Table 8). For both positively and negatively selected lineages, fewer substituted sites are buried relative to all sites on the protein; this is true both looking at all sites, and sites of any specific secondary structure, except for β-Sheet (p = 0.0361) and β-Bridge (p = 0.0081) sites on positively selected lineages, which was not significant after a multiple testing correction. The result in β-Bridge sites may simply be a matter of lower power due to the relatively small number of residues compared to most other secondary structures. β-Sheet sites are the most commonly substituted buried site on positive lineages (14.2744% vs 13.1684% for all helices), though α-Helix sites, as well as helices in general, are more common amongst all sites (15.9368 and 17.6017% vs 14.5822% for β-Sheet).

Table 8 Sitewise substitution rates in TAED lineages sorted by selective pressure and structural features
Full size table

Negatively selected lineages consistently have an increase in the prevalence of exposed residues across all secondary structures, but this is not universal for positively selected lineages. α-Helix sites are the most frequent in the dataset and show no change in prevalence of exposed sites compared to non-substituted sites under positive selection. 310 Helix sites show an overall increase in substitution rates in negatively selected lineages, unlike other helixes but consistent with bends, turns and coil sites. This is likely linked to their lower stability and higher proportion of exposed vs buried sites.

In terms of secondary structure when both exposed and buried regions are considered together, substitutions are more likely to occur across less structured regions (Turns, Bends, and Coil areas) that are more likely to be exposed than buried on both positively and negatively selected lineages, but also β-Sheet sites on positively selected lineages and 310 Helix sites on negatively selected lineages. The changes in prevalence for each secondary structure is strongly related to the buried/exposed ratio of their own residues (particularly in negatively selected sites), so solvent exposure, while a significant factor, is not the only one. This corresponds with observations seen in other studies ([18] and studies cited therein).

The lack of significant change in β-Sheet buried sites on positively selected lineages, suggests that positive selection is freer to act on it than comparable α-Helix sites, which have a considerable drop in frequency amongst substituted (13.1684%) rather than all (17.6017%) sites. The β-Sheet site changes also point at differences between positive and negative selection. Unlike in positively selected lineages, in negatively selected lineages, a smaller proportion of substituted sites are buried β-Sheet sites compared to all sites. This suggests the difference on positively selected lineages is not simply due to lower fragility in β-Sheet structure, but an active role for β-Sheet internal structure in driving evolution of new functionality. It should also be considered that, in general, positively selected lineages have fewer α-Helix (30.1108% vs 32.7617%) and more β-Sheet (21.7820% vs 19.8385%) sites compared to negatively selected lineages. Since, as discussed earlier, certain gene families and pathways are under more frequent positive selection than others, the lower selective constraint on β-Sheet sites has a long term impact on protein structure.

β-Bridge sites did not show a reduction in prevalence for substitutions on positively selected lineages. As these sites are used to hydrogen bond, particularly between β-sheets, the most likely source for these substitutions is to allow for protein restructuring. Purely compensatory driven changes are a less likely explanation, as negatively selected lineages where they are more likely than positively selected ones show a reduction in β-Bridge prevalence amongst substituted sites.

It should be noted that the same PDB structure is assumed to be applicable to all sequences in a gene family. As sequence pairs with divergence > 20% were split into separate families and as the median pairwise comparison among family members was 85% identity, the slow divergence of structural RMSD makes this a reasonable approximation [36]. Over longer evolutionary times [68, 69] and especially after lateral transfer events [60], repeated regions are known to lead to structural divergence.

Gene family analysis of ornithine decarboxylase

Lastly TAED can be a valuable resource in understanding the lineage-specific evolution of individual gene families. To examine this, one gene family was selected based on criteria that it contained KEGG pathway mappings and structural information. The gene family that was analyzed encoded a putative ornithine decarboxylase. Ornithine decarboxylase is responsible for the decarboxylation of L-ornithine to putrescine. L-ornithine is a key component to the urea cycle and the decarboxylation of L-ornithine signals the irreversible reaction of forming putrescine which is the first step in polyamine synthesis [59]. Polyamines are polycations able to bind negatively charged molecules such as DNA and RNA. Three primary polyamines are important regulators of the MAPK pathway which plays a role in cell proliferation: putrescine, spermidine, and spermine. Spermidine is produced from putrescine which can further impact apoptosis [50]. As these molecules play an important role in cell growth and cellular death, the committed step in the synthesis of polyamines would be hypothesized to be evolving under strong negative constraint.

An analysis of the TAED gene family showed six lineages with dN/dS > 1. These rates varied from a dN/dS rate of 2.0096 to 1.5451 (Table 9). Directional selection was found to have occurred on the lineage leading to Afrotherian mammals which are primarily localized to the continent of Africa and include: moles, elephants, manatees, and aardvarks. Other lineages with elevated rates of evolution were found for both Macaca mulatta (Rhesus macaque) and Dasypus novemcinctus (Nine-banded armadillo). Lastly, three different lineages involved cetacean species which may reflect the evolutionary pressures of moving from a terrestrial to an aquatic lifestyle. It was found that these instances of positive selection occurred following a duplication event, suggesting that the ornithine decarboxylase duplicate gene may have been under relaxed selective constraint following the duplication and not under the same strong constraints imposed by the polyamine synthesis pathway (Fig. 2). Although, since this protein was maintained and not lost over the 34 MYA of divergence between Orcinus orca (Killer whale) and Balaenoptera acutorostrata scammoni (Minke whale), it is likely that it has retained some functionality within these organisms.

Table 9 Lineages with dN/dS > 1 in Ornithine decarboxylase family
Full size table
Fig. 2
figure 2

Gene tree for cetacean lineages of ornithine decarboxylase – Presented here is the gene tree taken from the TAED Tree Viewer for the TAED gene family 557. Lineages not associated with Cetaceans are collapsed. Internal nodes labeled with a while box are duplication events found within the tree. Nodes with solid grey dots represent speciation events. Nodes labeled in black indicate a leaf node. Lineages labeled in red have a dN/dS > 1 and the numbers along each branch are the associated dN/dS value for the given branch. Image was generated from the TAED Tree Viewer

Full size image

To better understand the molecular mechanisms associated with the increased rate of evolution detected within the evolution of ornithine decarboxylase in cetaceans an examination of the ancestral changes mapped to the extant version of human ornithine decarboxylase was performed. For the changes on the branch Cetacea, it was seen that a nonsynonymous substitution occurred at site 238 with an asparagine substituting to an aspartic acid (N238D). This substitution is situated one residue from site 237 which is a known pyridoxal phosphate binding site [22] (Fig. 3. The decarboxylation of L-ornithine to putrescine is known to be a pyridoxal 5′-phosphate dependent reaction [37] and therefore changes to this site in the protein may impact the rate or ability to catalyze L-ornithine. The N238D substitution caused a substitution for an uncharged amino acid to be replaced by a negatively charged amino acid which could potentially impact the pyridoxal phosphate binding site (Fig. 3).

Fig. 3
figure 3

Pyridoxal phosphate binding site for ornithine decarboxylase along the lineage of Cetacea – A protein homology model of the ancestral protein leading to Cetacea was created. Template for the model was from human ornithine decarboxylase (PDB:2OO0; chain A). Ancestral changes occurring on the lineage for Cetacea have been mapped to the model, sites colored in red indicate nonsynonymous changes in the ancestral protein, sites colored in dark grey are synonymous site changes. The site indicated in green is the pyridoxal phosphate binding site 238. The site adjacent to the binding site is the substitution N238D found on the ancestral lineage. Image was generated from Swiss-PdbViewer

Full size image

The active site of ornithine decarboxylase in humans is at residue 357 (Cystine - 357) [3]. While no substitutions were found at the active site, four different nonsynonymous substitutions were localized on the beta-sheets surrounding the active site. The substitutions P368Q, R375C, I376M, and R379H were all proximally close to the active site and may have been involved in remodeling of the active site for the cetacean duplicate of ornithine decarboxylase (Fig. 4). These mutations have impacted the ability of the protein in several ways, by either helping to stabilize the active site, change the specificity of the binding pocket, change the rate of the reaction, or cause the active site to become inert. Further experimental validation would be necessary to understand how the N238D substitution and the putative remodeling of the active site may impact the function of the protein. However, evidence from TAED does suggest that cetacean ornithine decarboxylase has undergone functional shifts in several different sites which may impact the efficacy of the decarboxylation of L-ornithine to putrescine. Why this enzyme would be under selection within Cetaceans is also an unanswered question, but understanding the lineage-specific evolution of ornithine decarboxylase may help to decipher the mechanistic reasons for how cetaceans were able to readapt to life in the water.

Fig. 4
figure 4

Active site remodeling for ornithine decarboxylase along the lineage of Cetacea – A protein homology model of the ancestral protein leading to Cetacea was created. Template for the model was from human ornithine decarboxylase (PDB:2OO0; chain A). Ancestral changes occurring on the lineage for Cetacea have been mapped to the model, sites colored in red indicate nonsynonymous changes in the ancestral protein, sites colored in dark grey are synonymous site changes. The site indicated in gold is the active site cysteine-357. Remodeling of the active site can be seen in the changes P368Q, R375C, I376M, and R379H which are positioned around the loop containing the active site

Full size image

Discussion

Understanding the mechanistic reasons that species diverge is of central importance to the field of molecular evolution. Gaining insight into how individual proteins evolve in context of the pathways in which they occur may help elucidate the underlying molecular mechanisms of speciation. Placing evolutionary events in the context of a species tree allows the interpretation of understanding how selective forces have varied across species. Here we have presented findings from The Adaptive Evolution Database (TAED) that have attempted to characterize the lineage-specific evolution of chordates. We know that selection can act on multiple levels within an organism, from the level of individual nucleotides to phenotypic traits in a population. We therefore have examined the effects of directional selection at the domain level, gene level, and pathway level to better understand the dynamics of lineage-specific evolution. Examination of high level trends within TAED have confirmed that some pathways including those that are related to metabolism, immunity, and cell signaling have been repeated targets for functional change and may play important roles in species divergence. Additionally, we have shown that some protein families have undergone many duplication events which have impacted the evolutionary constraints of the duplicate pairs. These duplicated genes may evolve to new functions within the genome and develop new links within pathways. Tools developed on TAED can be utilized to find gene families that have undergone instances of adaptive evolution and help to propose hypotheses for how these genes have evolved.

Not all parts of a protein are under the same selective constraints and residues located on the outside or surface of a protein may be more likely to evolve, and evolve at a different rate, than a residue which comprises the hydrophobic core of the protein. Our comparison of the solvent accessible surface area (SASA) and dN/dS showed that this holds for both positively selected and negatively selected lineages. It distinguishes differences between the action of the two kinds of selection beyond this by showing that while solvent accessibility is more exclusively the primary driver of changes in the nature of substituted sites on negatively selected lineages, positively selected lineages show relaxed selective constraint on β-Sheet and strengthen constraints on α-Helix sites.

Additionally, the relationship between the energetics of different substitutions and how they interplay with dN/dS could be explored by comparing dN/dS to the change in the change of free energy (ΔΔG) of a protein when different substitutions are introduced. Studies of this nature have examined how the thermodynamics of a protein influence the rate of dN/dS and how compensatory substitutions affect protein stability [61, 70]. Current evolutionary tests do not consider epistatic relationships within proteins, treating each site as acting independently from a statistical perspective.

Further, it is known that when Ne is large, selection is more efficient and the chance of an allele being lost from the population is small. However, when Ne is small the effects of genetic drift are greater and selection is less efficient [49, 75]. As such selection has limited ability to eliminate deleterious variants in chordates or fix advantageous changes, as chordate species have low effective populations sizes. Weber, et al. [80] found an unexpected negative correlation between Ne and dN/dS in bird populations, but found expected signals when considering the magnitude of biophysical effects of changes [80, 81].

TAED as a tool and resource in detecting episodes of lineage-specific evolution may also be useful in helping to understand the differences between directional selection and intra- and inter- molecular forces. Not all amino acid substitutions are the direct result of directional selection acting on a protein to functionally evolve. When physical changes within a molecule do occur, corresponding compensatory changes can occur which alleviate the deleterious effects of a mutation. These compensatory changes ensure that the newly substituted amino acid becomes the preferred amino acid for the residue in which it is located [61, 70]. Using traditional approaches of dN/dS it is difficult to differentiate between directional selection and compensatory changes as both aggregate across the branch. However, by examining changes in a lineage-specific context and determining when each substitution occurred along the lineage, it may be possible to begin to differentiate between these two processes.

The secondary structure analysis raises questions about the nature of the selective pressures on a protein-structure level, and points to the need for further investigation of β-sheet, α-helix, and 310 Helix structures and their role in protein evolution in particular.

Conclusions

TAED is a useful tool for understanding lineage-specific evolution and provides a source of data to develop further hypothesis-based inquiries into the mechanisms that drive diversification. In addition to providing an example of lineage-specific evolution in cetaceans, this work examined gene family evolution through the lenses of protein structure, co-evolution in pathways, as well as characterizing the duplication process within families. At the structural level, the study utilized the database to understand the differential patterns of amino acid substitution, including filtering by secondary structure, in comparing proteins under negative and positive selection. Overall, this work provides a further empirical window into the lineage-specific processes of evolution.

Methods

Database construction

The TAED database was constructed following the pipeline outlined in Hermansen et al. [33]. The pipeline includes generation of gene families from single-linkage clustering of BLAST results from chordate genes found on GenBank. A point accepted mutation (PAM) distance threshold of 120 was used for gene family construction. Gene families were refined for quality using an iterative method controlling for pairwise percent identity (> 80%) and the fraction of pairwise aligned gaps (< 10%). Gene families where then aligned using MAFFT [41] and phylogenetic trees were constructed using PhyML [30]. Gene tree – species tree reconciliation against the NCBI chordate taxonomy was implemented to determine putative duplication events and gene tree roots using SoftParsMap. Gene families were defined phylogenetically by the species tree except in cases where alignment quality prohibited this, as described here and in Hermansen et al. [33] (see [4] for a recent discussion of gene family construction methodology). Putative rates of evolution were then calculated using the branches model from PAML and dN/dS rates was computed. BLAST was then performed on TAED gene families against the KEGG database [40] to determine KEGG pathway relatedness and against PDB [10] to determine protein structure for each gene in TAED. All branches, including specifically those found to have a dN/dS > 1 (putatively evolving under positive selection) were mapped to the corresponding chordate species tree to determine along what lineage the elevated rates of evolution occurred and which proteins evolved rapidly on the same species tree lineage. Roots of all genes families were additionally mapped to the chordates species tree. To determine the approximate family root age for each gene family, information from TimeTree [43] was collected and root ages determined in MYA (millions of years ago). Domain classification information was gathered from the CATH database [71]. Putative functional annotations were assigned to each gene family based on NCBI nomenclature and KEGG pathway annotations when available.

Enrichment analysis

Over/under-represented KEGG pathway and domain analyses were performed with a BLAST search against the KEGG database of TAED gene families. KO numbers were assigned to each individual protein in TAED that contained a BLAST hit with an e-value <1e− 10. This threshold was set so that all putative hits would be the result of orthologous descent instead of chance. The KO number from the top BLAST result was assigned to each TAED gene. KO numbers were then used to assess each putative biological pathway in which the protein is known to play a role. Over/under-representation of these pathways was then calculated using Fisher’s Exact test [23] and significance was estimated using an α-level of 0.05. The resulting p-values were corrected for multiple testing by performing a false discovery rate (FDR) analysis [8] with an FDR threshold of 0.05 and using a Bonferroni correction [13]. The FDR calculation was computed using the R statistical programming package [62]. A similar method was used to determine the over/under-representation of CATH domain topologies. The topology level classification was used as it represented a broad enough group that multiple topologies were found throughout TAED.

Duplication analysis

For each gene family in TAED, the root node of the family was mapped to its associated lineage on the chordate species tree. Nodes were then given approximate dates in MYA based on estimates from the TimeTree database [43]. The number of duplication events that occurred in each gene family was used as inferred by SoftParsMap [9] through reconciliation with the NCBI taxonomy for chordates. A linear regression was performed on the resulting comparison between the family root node ages and the number of duplication found within each gene family. The Pearson’s r coefficient was calculated for the resulting linear regression with a Pearson’s r = 0.59. Log scaled transformations of the data did not yield a strong regression coefficient.

Since families were sought that showed a high propensity for duplicability in a short amount of time, families that fell below the regression line were filtered out (Fig. 1). We also filtered out all families whose length was below the 5th percentile, since evolutionary forces may not have had time to act on families with so few substitutions. Outliers in the resulting set of families were detected using Cook’s distance [20], which measures the change in regression coefficients due to the removal of a data point, and is often used as a proxy for the influence of that point. Gene families were then sorted according to this distance (Table 6). Finally, the top quartile of families was measured using this distance and the number of times they occur in each KEGG pathway was counted (Table 7).

Additionally, to test how different proxies of time impacted the duplication analysis, two additional proxies for time were generated: the maximum tree length, and the median tree length. The maximum tree length estimated in substitutions per site was calculated for all gene tree topologies by taking the maximal tree length from root to leaf node for every TAED gene family as estimated by PhyML. The median tree length was calculated in a similar manner by taking the median of all distances between the root and leaf of the phylogenetic tree for each gene family. Additional file 1: Figures S1 and S2 illustrate the differences in the duplication distribution of the families based on the change of the time component to the analysis. Each axis was of the analysis was given the transformation y = log (1 + x) and the Pearson’s r coefficient was calculated. The resulting best coefficients for both the maximum tree length and the median tree were found when both axes were log-transformed. Cook’s distance was calculated for each proxy of time and the families with pathways from the families with the largest Cook’s distance to the regression line were tabulated.

Protein structure based analysis

Protein information was determined from stored PDB information associated with each gene family. To show that sites in different locations and belonging to different structures evolve at different rates, DSSP [39] values were used to ascertain the relative solvent accessibility (RSA) and secondary structure of individual sites within the protein was obtained. While newer and less approximate, but more computationally intensive methods than DSSP are available, a pilot analysis suggested that DSSP and more computationally intensive methods gave similar results for the purposes of this study. Membrane proteins and multimers were removed from the dataset based on identifying information in the PDB data. Sites were binned based on RSA using maximum surface areas from Tien et al. [77]; sites with a ratio greater than 0.20 were marked as exposed and buried otherwise, and then further categorized according to secondary structure. PAML analysis was used to determine the maximum likelihood ancestral sequence for each gene associated with a protein and the results controlled for lineages with dN/dS > 1 and lineages with a dN/dS <  0.5. dN/dS values of 0 or between 0.5 and 1 were ignored, as were any sites that did not align with the PDB sequence or were not one of the most common 20 amino acids. To determine the significance of the calculated values, two-tailed non-parametric bootstrapping was performed. For each lineage, simulated datasets of size matching the total substituted residue count were generated, using the distribution of all sites on the respective lineages as a baseline.

Gene family analysis of ornithine decarboxylase

To demonstrate the application of lineage-specific analyses of evolution on specific gene families using TAED data, a gene family was selected for analysis based on the criteria that the gene family contained 3 or more lineages with dN/dS > 1 and it contained lineages that mapped to KEGG pathways and to a PDB structure. Using these criteria, the TAED gene family 554 (ornithine decarboxylase) was selected for further examination of lineage-specific evolution. dN/dS estimates of each lineage were taken from the TAED database. A homology model was generated using Swiss-Model [11], with the automated build method. The top template used in the homology model was PDB entry 2OO0 chain A. Ancestral amino acids were mapped to the model. Active site and binding site information was taken from the PDB website for the same entry. Uniprot [7] data for ornithine decarboxylase was also used to make inferences into important catalytic sites within the molecule. Images of the homology model were generated using Swiss-PdbViewer [29].

Availability of data and materials

The TAED Database is accessible at https://liberles.cst.temple.edu/TAED/index.html. Code available upon request.

References

  1. Aik W, McDonough MA, Thalhammer A, Chowdhury R, Schofield CJ. Role of the jelly-roll fold in substrate binding by 2-oxoglutarate oxygenases. Curr Opin Struct Biol. 2012;22(6):691–700. https://0-doi-org.brum.beds.ac.uk/10.1016/J.SBI.2012.10.001.

    Article  CAS  PubMed  Google Scholar 

  2. Aken BL, Achuthan P, Akanni W, Amode MR, Bernsdorff F, Bhai J, et al. Ensembl 2017. Nucleic Acids Res. 2017;45(D1):D635–42. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkw1104.

    Article  CAS  PubMed  Google Scholar 

  3. Almrud JJ, Oliveira MA, Kern AD, Grishin NV, Phillips MA, Hackert ML. Crystal structure of human ornithine decarboxylase at 2.1 å resolution: structural insights to antizyme binding. J Mol Biol. 2000;295(1):7–16. https://0-doi-org.brum.beds.ac.uk/10.1006/JMBI.1999.3331.

    Article  CAS  PubMed  Google Scholar 

  4. Altenhoff AM, Boeckmann B, Capella-Gutierrez S, Dalquen DA, DeLuca T, Forslund K, Huerta-Cepas J, Linard B, Pereira C, Pryszcz LP, Schreiber F, da Silva AS, Szklarczyk D, Train CM, Bork P, Lecompte O, von Mering C, Xenarios I, Sjölander K, Jensen LJ, Martin MJ, Muffato M, Quest for Orthologs consortium, Gabaldón T, Lewis SE, Thomas PD, Sonnhammer E, Dessimoz C. Standardized benchmarking in the quest for orthologs. Nat Methods. 2016;13(5):425–30. https://0-doi-org.brum.beds.ac.uk/10.1038/nmeth.3830.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Anisimova M, Yang Z. Multiple hypothesis testing to detect lineages under positive selection that affects only a few sites. Mol Biol Evol. 2007;24(5):1219–28. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msm042.

    Article  CAS  PubMed  Google Scholar 

  6. Anisimova M, Liberles DA. Detecting and understanding natural selection. In: Codon evolution: Oxford University Press; 2012. p. 73–96. https://0-doi-org.brum.beds.ac.uk/10.1093/acprof:osobl/9780199601165.003.0006.

    Chapter  Google Scholar 

  7. Apweiler R, Bairoch A, Wu CH, Barker WC, Boeckmann B, Ferro S, et al. UniProt: the universal protein knowledgebase. Nucleic Acids Res. 2004;32(90001):115D–19. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkh131.

    Article  CAS  Google Scholar 

  8. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J Roy Sta Soc. 1995;57:289–300.

    Google Scholar 

  9. Berglund-Sonnhammer A-C, Steffansson P, Betts MJ, Liberles DA. Optimal gene trees from sequences and species trees using a soft interpretation of parsimony. J Mol Evol. 2006;63(2):240–50. https://0-doi-org.brum.beds.ac.uk/10.1007/s00239-005-0096-1.

    Article  CAS  PubMed  Google Scholar 

  10. Berman HM, Westbrook J, Feng Z, Gilliland G, Bhat TN, Weissig H, et al. The protein data bank. Nucleic Acids Res. 2000;28(1):235–42. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/28.1.235.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Biasini M, Bienert S, Waterhouse A, Arnold K, Studer G, Schmidt T, et al. SWISS-MODEL: modelling protein tertiary and quaternary structure using evolutionary information. Nucleic Acids Res. 2014;42(Web Server issue):W252–8. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gku340.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Bloch NI, Morrow JM, Chang BSW, Price TD. SWS2 visual pigment evolution as a test of historically contingent patterns of plumage color evolution in warblers. Evolution. 2015;69(2):341–56. https://0-doi-org.brum.beds.ac.uk/10.1111/evo.12572.

    Article  PubMed  PubMed Central  Google Scholar 

  13. Bonferroni C. Teoria statistica delle classi e calcolo delle probabilita. Pubblicazioni Del R Istituto Superiore Di Scienze Economiche e Commericiali Di Firenze. 1936;8:3–62 Retrieved from https://ci.nii.ac.jp/naid/20001561442/.

    Google Scholar 

  14. Boucher JI, Bolon DNA, Tawfik DS. Quantifying and understanding the fitness effects of protein mutations: laboratory versus nature. Protein Sci. 2016;25(7):1219–26. https://0-doi-org.brum.beds.ac.uk/10.1002/pro.2928.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Brockhurst MA, Chapman T, King KC, Mank JE, Paterson S, Hurst GDD. Running with the red queen: the role of biotic conflicts in evolution. Proc R Soc B Biol Sci. 2014;281(1797):20141382. https://0-doi-org.brum.beds.ac.uk/10.1098/rspb.2014.1382.

    Article  Google Scholar 

  16. Brown SD, Babbitt PC. New insights about enzyme evolution from large scale studies of sequence and structure relationships. J Biol Chem. 2014;289(44):30221–8. https://0-doi-org.brum.beds.ac.uk/10.1074/jbc.R114.569350.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Cheng S, Brooks CL 3rd. Viral capsid proteins are segregated in structural fold space. PLoS Comput Biol. 2013;9(2):e1002905. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pcbi.1002905.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Chi PB, Kim D, Lai JK, Bykova N, Weber CC, Kubelka J, Liberles DA. A new parameter-rich structure-aware mechanistic model for amino acid substitution during evolution. Proteins. 2018;86(2):218–28. https://0-doi-org.brum.beds.ac.uk/10.1002/prot.25429.

    Article  CAS  PubMed  Google Scholar 

  19. Chi PB, Liberles DA. Selection on protein structure, interaction, and sequence. Protein Sci. 2016;25(7):1168–78. https://0-doi-org.brum.beds.ac.uk/10.1002/pro.2886.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Cook RD. Detection of influential observation in linear regression. Technometrics. 1977;19(1):15–8. https://0-doi-org.brum.beds.ac.uk/10.1080/00401706.1977.10489493.

    Article  Google Scholar 

  21. Deming PB, Campbell SL, Stone JB, Rivard RL, Mercier AL, Howe AK. Anchoring of protein kinase A by ERM (ezrin-radixin-moesin) proteins is required for proper netrin signaling through DCC (deleted in colorectal cancer). J Biol Chem. 2015;290(9):5783-96. https://0-doi-org.brum.beds.ac.uk/10.1074/jbc.M114.628644.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Dufe VT, Ingner D, Heby O, Khomutov AR, Persson L, Al-Karadaghi S. A structural insight into the inhibition of human and Leishmania donovani ornithine decarboxylases by 1-amino-oxy-3-aminopropane. Biochem J. 2007;405(2):261–8. https://0-doi-org.brum.beds.ac.uk/10.1042/BJ20070188.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Fisher RA. On the interpretation of χ 2 from contingency tables, and the calculation of P. J R Stat Soc. 1922;85(1):87. https://0-doi-org.brum.beds.ac.uk/10.2307/2340521.

    Article  Google Scholar 

  24. Gatesy J, Geisler JH, Chang J, Buell C, Berta A, Meredith RW, et al. A phylogenetic blueprint for a modern whale. Mol Phylogenet Evol. 2013;66(2):479–506. https://0-doi-org.brum.beds.ac.uk/10.1016/J.YMPEV.2012.10.012.

    Article  PubMed  Google Scholar 

  25. Gharib WH, Robinson-Rechavi M. The branch-site test of positive selection is surprisingly robust but lacks power under synonymous substitution saturation and variation in GC. Mol Biol Evol. 2013;30(7):1675–86. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/mst062.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Glusman G, Bahar A, Sharon D, Pilpel Y, White J, Lancet D. The olfactory receptor gene superfamily: data mining, classification, and nomenclature. Mamm Genome. 2000;11(11):1016–23 Retrieved from http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/pubmed/11063259.

    Article  CAS  PubMed  Google Scholar 

  27. Goldman AD, Beatty JT, Landweber LF. The TIM barrel architecture facilitated the early evolution of protein-mediated metabolism. J Mol Evol. 2016;82(1):17–26. https://0-doi-org.brum.beds.ac.uk/10.1007/s00239-015-9722-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Grahnen JA, Nandakumar P, Kubelka J, Liberles DA. Biophysical and structural considerations for protein sequence evolution. BMC Evol Biol. 2011;11(1):361. https://0-doi-org.brum.beds.ac.uk/10.1186/1471-2148-11-361.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Guex N, Peitsch MC, Schwede T. Automated comparative protein structure modeling with SWISS-MODEL and Swiss-PdbViewer: a historical perspective. Electrophoresis. 2009;30(S1):S162–73. https://0-doi-org.brum.beds.ac.uk/10.1002/elps.200900140.

    Article  PubMed  Google Scholar 

  30. Guindon S, Dufayard J-F, Lefort V, Anisimova M, Hordijk W, Gascuel O. New algorithms and methods to estimate maximum-likelihood phylogenies: assessing the performance of PhyML 3.0. Syst Biol. 2010;59(3):307–21. https://0-doi-org.brum.beds.ac.uk/10.1093/sysbio/syq010.

    Article  CAS  PubMed  Google Scholar 

  31. Hanukoglu I. Proteopedia: Rossmann fold: a beta-alpha-beta fold at dinucleotide binding sites. Biochem Mol Biol Educ. 2015;43(3):206–9. https://0-doi-org.brum.beds.ac.uk/10.1002/bmb.20849.

    Article  CAS  PubMed  Google Scholar 

  32. Hermansen RA, Mannakee BK, Knecht W, Liberles DA, Gutenkunst RN. Characterizing selective pressures on the pathway for de novo biosynthesis of pyrimidines in yeast. BMC Evol Biol. 2015;15(1):232. https://0-doi-org.brum.beds.ac.uk/10.1186/s12862-015-0515-x.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Hermansen RA, Oswald BP, Knight S, Shank SD, Northover D, Korunes KL, et al. The adaptive evolution database (TAED): a new release of a database of Phylogenetically indexed gene families from chordates. J Mol Evol. 2017;85(1–2):46–56. https://0-doi-org.brum.beds.ac.uk/10.1007/s00239-017-9806-8.

    Article  CAS  PubMed  Google Scholar 

  34. Horowitz NH. On the evolution of biochemical syntheses. Proc Natl Acad Sci U S A. 1945;31(6):153–7. https://doi.org/10.1073/PNAS.31.6.153.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Hughes T, Liberles DA. The pattern of evolution of smaller-scale gene duplicates in mammalian genomes is more consistent with neo- than subfunctionalisation. J Mol Evol. 2007;65(5):574–88. https://0-doi-org.brum.beds.ac.uk/10.1007/s00239-007-9041-9.

    Article  CAS  PubMed  Google Scholar 

  36. Illergård K, Ardell DH, Elofsson A. Structure is three to ten times more conserved than sequence--a study of structural response in protein cores. Proteins. 2009;77(3):499–508. https://0-doi-org.brum.beds.ac.uk/10.1002/prot.22458.

    Article  CAS  PubMed  Google Scholar 

  37. Jackson LK, Baldwin J, Akella R, Goldsmith EJ, Phillips MA. Multiple active site conformations revealed by distant site mutation in ornithine decarboxylase. Biochemistry. 2004;43(41):12990–9. https://0-doi-org.brum.beds.ac.uk/10.1021/bi048933l.

    Article  CAS  PubMed  Google Scholar 

  38. Jensen RA. Enzyme recruitment in evolution of new function. Annu Rev Microbiol. 1976;30(1):409–25. https://0-doi-org.brum.beds.ac.uk/10.1146/annurev.mi.30.100176.002205.

    Article  CAS  PubMed  Google Scholar 

  39. Kabsch W, Sander C. Dictionary of protein secondary structure: pattern recognition of hydrogen-bonded and geometrical features. Biopolymers. 1983;22(12):2577–637. https://0-doi-org.brum.beds.ac.uk/10.1002/bip.360221211.

    Article  CAS  PubMed  Google Scholar 

  40. Kanehisa M, Sato Y, Kawashima M, Furumichi M, Tanabe M. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Res. 2016;44(D1):D457–62. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkv1070.

    Article  CAS  PubMed  Google Scholar 

  41. Katoh K, Standley DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013;30(4):772–80. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/mst010.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Konrad A, Teufel AI, Grahnen JA, Liberles DA. Toward a general model for the evolutionary dynamics of gene duplicates. Genome Biol Evol. 2011;3:1197–209. https://0-doi-org.brum.beds.ac.uk/10.1093/gbe/evr093.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Kumar S, Stecher G, Suleski M, Hedges SB. TimeTree: a resource for timelines, Timetrees, and divergence times. Mol Biol Evol. 2017;34(7):1812–9. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msx116.

    Article  CAS  PubMed  Google Scholar 

  44. Kurtz J, Schulenburg H, Reusch TBH. Host–parasite coevolution—rapid reciprocal adaptation and its genetic basis. Zoology. 2016;119(4):241–3. https://0-doi-org.brum.beds.ac.uk/10.1016/j.zool.2016.06.011.

    Article  PubMed  Google Scholar 

  45. Laurino P, Tóth-Petróczy Á, Meana-Pañeda R, Lin W, Truhlar DG, Tawfik DS. An ancient fingerprint indicates the common ancestry of Rossmann-fold enzymes utilizing different ribose-based cofactors. PLoS Biol. 2016;14(3):e1002396. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pbio.1002396.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Liberles DA, Schreiber DR, Govindarajan S, Chamberlin SG, Benner SA. The adaptive evolution database (TAED). Genome Biol. 2001;2(8):RESEARCH0028 Retrieved from http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/pubmed/11532212.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Liberles DA, Tisdell MDM, Grahnen JA. Binding constraints on the evolution of enzymes and signalling proteins: the important role of negative pleiotropy. Proc Biol Sci. 2011;278(1714):1930–5. https://0-doi-org.brum.beds.ac.uk/10.1098/rspb.2010.2637.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Light S, Kraulis P. Network analysis of metabolic enzyme evolution in Escherichia coli. BMC Bioinformatics. 2004;5:15. https://0-doi-org.brum.beds.ac.uk/10.1186/1471-2105-5-15.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Lynch M, Bobay L-M, Catania F, Gout J-F, Rho M. The repatterning of eukaryotic genomes by random genetic drift. Annu Rev Genomics Hum Genet. 2011;12:347–66. https://0-doi-org.brum.beds.ac.uk/10.1146/annurev-genom-082410-101412.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Minois N, Carmona-Gutierrez D, Madeo F. Polyamines in aging and disease. Aging. 2011;3(8):716–32. https://0-doi-org.brum.beds.ac.uk/10.18632/aging.100361.

    Article  PubMed  PubMed Central  Google Scholar 

  51. Moretti S, Laurenczy B, Gharib WH, Castella B, Kuzniar A, Schabauer H, et al. Selectome update: quality control and computational improvements to a database of positive selection. Nucleic Acids Res. 2014;42(Database issue):D917–21. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkt1065.

    Article  CAS  PubMed  Google Scholar 

  52. Newfeld SJ, Wisotzkey RG, Kumar S. Molecular evolution of a developmental pathway: phylogenetic analyses of transforming growth factor-beta family ligands, receptors and Smad signal transducers. Genetics. 1999;152(2):783–95 Retrieved from http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/pubmed/10353918.

    CAS  PubMed  PubMed Central  Google Scholar 

  53. Newfeld SJ, Wisotzkey RG. Molecular evolution of Smad proteins. In: Smad signal transduction. Dordrecht: Springer Netherlands; 2006. p. 15–35. https://0-doi-org.brum.beds.ac.uk/10.1007/1-4020-4709-6_1.

    Chapter  Google Scholar 

  54. Niimura Y, Nei M. Extensive gains and losses of olfactory receptor genes in mammalian evolution. PLoS One. 2007;2(8):s. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pone.0000708.

    Article  CAS  Google Scholar 

  55. Ohno S. Evolution by Gene Duplication. Berlin: Springer-Verlag; 1970.

  56. Orlenko A, Chi PB, Liberles DA. Characterizing the roles of changing population size and selection on the evolution of flux control in metabolic pathways. BMC Evol Biol. 2017;17(1):117. https://0-doi-org.brum.beds.ac.uk/10.1186/s12862-017-0962-7.

    Article  PubMed  PubMed Central  Google Scholar 

  57. Orlenko A, Teufel AI, Chi PB, Liberles DA. Selection on metabolic pathway function in the presence of mutation-selection-drift balance leads to rate-limiting steps that are not evolutionarily stable. Biol Direct. 2016;11(1):31. https://0-doi-org.brum.beds.ac.uk/10.1186/s13062-016-0133-6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Pang K, Cheng C, Xuan Z, Sheng H, Ma X. Understanding protein evolutionary rate by integrating gene co-expression with protein interactions. BMC Syst Biol. 2010;4:179. https://0-doi-org.brum.beds.ac.uk/10.1186/1752-0509-4-179.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Pegg AE. Regulation of ornithine decarboxylase. J Biol Chem. 2006;281(21):14529–32. https://0-doi-org.brum.beds.ac.uk/10.1074/jbc.R500031200.

    Article  CAS  PubMed  Google Scholar 

  60. Persi E, Wolf YI, Koonin EV. Positive and strongly relaxed purifying selection drive the evolution of repeats in proteins. Nat Commun. 2016;7:13570. https://0-doi-org.brum.beds.ac.uk/10.1038/ncomms13570.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Pollock DD, Thiltgen G, Goldstein RA. Amino acid coevolution induces an evolutionary stokes shift. Proc Natl Acad Sci U S A. 2012;109(21):E1352–9. https://0-doi-org.brum.beds.ac.uk/10.1073/pnas.1120084109.

    Article  PubMed  PubMed Central  Google Scholar 

  62. R Core Team. R: A Language and Environment for Statistical Computing: R Foundation for Statistical Computing; 2016.

  63. Rao ST, Rossmann MG. Comparison of super-secondary structures in proteins. J Mol Biol. 1973;76(2):241–56 Retrieved from http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/pubmed/4737475.

    Article  CAS  PubMed  Google Scholar 

  64. Raoult D, Forterre P. Redefining viruses: lessons from Mimivirus. Nat Rev Microbiol. 2008;6(4):315–9. https://0-doi-org.brum.beds.ac.uk/10.1038/nrmicro1858.

    Article  CAS  PubMed  Google Scholar 

  65. Rastogi S, Liberles DA. Subfunctionalization of duplicated genes as a transition state to neofunctionalization. BMC Evol Biol. 2005;5:28.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Roth C, Betts MJ, Steffansson P, Saelensminde G, Liberles DA. The Adaptive Evolution Database (TAED): a phylogeny based tool for comparative genomics. Nucleic Acids Res. 2005;33(Database issue):D495–7. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gki090.

    Article  CAS  PubMed  Google Scholar 

  67. Sayers EW, Cavanaugh M, Clark K, Ostell J, Pruitt KD, Karsch-Mizrachi I. GenBank. Nucleic Acids Res. 2020;48(D1):D84–D86. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkz956.

  68. Schaeffer RD, Kinch LN, Liao Y, Grishin NV. Classification of proteins with shared motifs and internal repeats in the ECOD database. Protein Sci. 2016;25(7):1188–203. https://0-doi-org.brum.beds.ac.uk/10.1002/pro.2893.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Schüler A, Bornberg-Bauer E. Evolution of protein domain repeats in Metazoa. Mol Biol Evol. 2016;33(12):3170–82.

    Article  PubMed  PubMed Central  Google Scholar 

  70. Shah P, McCandlish DM, Plotkin JB. Contingency and entrenchment in protein evolution under purifying selection. Proc Natl Acad Sci U S A. 2015;112(25):E3226–35. https://0-doi-org.brum.beds.ac.uk/10.1073/pnas.1412933112.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Sillitoe I, Lewis TE, Cuff A, Das S, Ashford P, Dawson NL, et al. CATH: comprehensive structural and functional annotations for genome sequences. Nucleic Acids Res. 2015;43(Database issue):D376–81. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gku947.

    Article  CAS  PubMed  Google Scholar 

  72. Spady TC, Seehausen O, Loew ER, Jordan RC, Kocher TD, Carleton KL. Adaptive molecular evolution in the opsin genes of rapidly speciating cichlid species. Mol Biol Evol. 2005;22(6):1412–22.

    Article  CAS  PubMed  Google Scholar 

  73. Spielman SJ, Wilke CO. The relationship between dN/dS and scaled selection coefficients. Mol Biol Evol. 2015;32(4):1097–108. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msv003.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Sugawara T, Imai H, Nikaido M, Imamoto Y, Okada N. Vertebrate rhodopsin adaptation to dim light via rapid meta-II intermediate formation. Mol Biol Evol. 2010;27(3):506–19. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msp252.

    Article  CAS  PubMed  Google Scholar 

  75. Sung W, Ackerman MS, Miller SF, Doak TG, Lynch M. Drift-barrier hypothesis and mutation-rate evolution. Proc Natl Acad Sci U S A. 2012;109(45):18488–92. https://0-doi-org.brum.beds.ac.uk/10.1073/pnas.1216223109.

    Article  PubMed  PubMed Central  Google Scholar 

  76. Teufel AI, Liu L, Liberles DA. Models for gene duplication when dosage balance works as a transition state to subsequent neo- or sub-functionalization. BMC Evol Biol. 2016;16:45. https://0-doi-org.brum.beds.ac.uk/10.1186/s12862-016-0616-1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Tien MZ, Meyer AG, Sydykova DK, Spielman SJ, Wilke CO. Maximum allowed solvent accessibilites of residues in proteins. PLoS One. 2013;8(11):e80635. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pone.0080635.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Van Valen L. Molecular evolution as predicted by natural selection. J Mol Evol. 1974;3(2):89–101 Retrieved from http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/pubmed/4407466.

    Article  PubMed  Google Scholar 

  79. Weadick CJ, Loew ER, Rodd FH, Chang BSW. Visual pigment molecular evolution in the Trinidadian pike cichlid (Crenicichla frenata): a less colorful world for neotropical cichlids? Mol Biol Evol. 2012;29(10):3045–60. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/mss115.

    Article  CAS  PubMed  Google Scholar 

  80. Weber CC, Nabholz B, Romiguier J, Ellegren H. K r /K c but not d N /d S correlates positively with body mass in birds, raising implications for inferring lineage-specific selection. Genome Biol. 2014;15(12):542. https://0-doi-org.brum.beds.ac.uk/10.1186/s13059-014-0542-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Weber CC, Whelan S. Physicochemical amino acid properties better describe substitution rates in large populations. Mol Biol Evol. 2019;2019. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msz003 in press.

    Article  CAS  PubMed  Google Scholar 

  82. Wierenga R. The TIM-barrel fold: a versatile framework for efficient enzymes. FEBS Lett. 2001;492(3):193–8. https://0-doi-org.brum.beds.ac.uk/10.1016/S0014-5793(01)02236-0.

    Article  CAS  PubMed  Google Scholar 

  83. Yang Z. Likelihood ratio tests for detecting positive selection and application to primate lysozyme evolution. Mol Biol Evol. 1998;15(5):568–73. https://0-doi-org.brum.beds.ac.uk/10.1093/oxfordjournals.molbev.a025957.

    Article  CAS  PubMed  Google Scholar 

  84. Zhang G, Li C, Li Q, Li B, Larkin DM, Lee C, et al. Comparative genomics reveals insights into avian genome evolution and adaptation. Science. 2014;346(6215):1311–20. https://0-doi-org.brum.beds.ac.uk/10.1126/science.1251385.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgments

Russell A. Hermansen was originally the first author on this work and contributed at an author level. He has left academic science and became non-responsive in our attempts to contact him about this manuscript. After discussion with the Editorial Advisory Group at BMC on how to handle this, he has been moved from the author list to the acknowledgments. We also thank Peter Chi, Tristan Stark, Mark Walker, Kelly Clark, and Jordan Zehr for helpful discussions.

Funding

The authors would like to thank National Science Foundation Grants DBI-0743374 and DBI-1355846 for support as well as the University of Wyoming INBRE Award P20 RR016474. DBI-0743374 and P20 RR016474 supported RAH in this analysis and in the creation of the original database. SDS received support from DBI-1355846.

Author information

Authors and Affiliations

  1. Department of Biology and Center for Computational Genetics and Genomics, Temple University, Philadelphia, PA, 19122, USA

    David E. Northover, Stephen D. Shank & David A. Liberles

  2. Department of Molecular Biology, University of Wyoming, Laramie, WY, 82071, USA

    David A. Liberles

Authors
  1. David E. Northover
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Stephen D. Shank
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. David A. Liberles
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

DAL conceived the work in this study together with Russell A. Hermansen (RAH). RAH, DEN, and SDS performed the analyses in this work. RAH wrote the initial version of the manuscript, with additional writing by DEN and DAL. All authors read and approved the final manuscript.

Corresponding author

Correspondence to David A. Liberles.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Additional file 1: Figure S1.

The logarithm of the maximum phylogenetic tree length of TAED gene families, regressed against the logarithm of the number of duplications in a given family. Maximum length consists of the maximum cumulative branch length from the corresponding phylogenetic trees, considering all paths from root to tips. A log-scale was chosen since it resulted in a higher correlation coefficient. Figure S2. The logarithm of the median phylogenetic tree length of TAED gene families, regressed against the logarithm of the number of duplications in a given family. Meidan length consists of the median cumulative branch length from the corresponding phylogenetic trees, considering all paths from root to tips. A log-scale was chosen since it resulted in a higher correlation coefficient. Table S1. TAED gene families with many duplications based on maximum tree length. Table S2. TAED gene families with many duplications based on median tree length. Table S3. KEGG Pathways with many duplications based on maximum tree length. Table S4. KEGG Pathways with many duplications based on median tree length.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Northover, D.E., Shank, S.D. & Liberles, D.A. Characterizing lineage-specific evolution and the processes driving genomic diversification in chordates. BMC Evol Biol 20, 24 (2020). https://0-doi-org.brum.beds.ac.uk/10.1186/s12862-020-1585-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://0-doi-org.brum.beds.ac.uk/10.1186/s12862-020-1585-y

Keywords

BMC Ecology and Evolution

ISSN: 2730-7182

Contact us