Adaptive evolution and functional constraint at TLR4 during the secondary aquatic adaptation and diversification of cetaceans.

BACKGROUND: Cetaceans (whales, dolphins and porpoises) are a group of adapted marine mammals with an enigmatic history of transition from terrestrial to full aquatic habitat and rapid radiation in waters around the world. Throughout this evolution, the pathogen stress-response proteins must have faced challenges from the dramatic change of environmental pathogens in the completely different ecological niches cetaceans occupied. For this reason, cetaceans could be one of the most ideal candidate taxa for studying evolutionary process and associated driving mechanism of vertebrate innate immune systems such as Toll-like receptors (TLRs), which are located at the direct interface between the host and the microbial environment, act at the first line in recognizing specific conserved components of microorganisms, and translate them rapidly into a defense reaction.
RESULTS: We used TLR4 as an example to test whether this traditionally regarded pattern recognition receptor molecule was driven by positive selection across cetacean evolutionary history. Overall, the lineage-specific selection test showed that the dN/dS (ω) values along most (30 out of 33) examined cetartiodactylan lineages were less than 1, suggesting a common effect of functional constraint. However, some specific codons made radical changes, fell adjacent to the residues interacting with lipopolysaccharides (LPS), and showed parallel evolution between independent lineages, suggesting that TLR4 was under positive selection. Especially, strong signatures of adaptive evolution on TLR4 were identified in two periods, one corresponding to the early evolutionary transition of the terrestrial ancestors of cetaceans from land to semi-aquatic (represented by the branch leading to whale + hippo) and from semi-aquatic to full aquatic (represented by the ancestral branch leading to cetaceans) habitat, and the other to the rapid diversification and radiation of oceanic dolphins.
CONCLUSIONS: This is the first study thus far to characterize the TLR gene in cetaceans. Our data present evidences that cetacean TLR4 has undergone adaptive evolution against the background of purifying selection in response to the secondary aquatic adaptation and rapid diversification in the sea. It is suggested that microbial pathogens in different environments are important factors that promote adaptive changes at cetacean TLR4 and new functions of some amino acid sites specialized for recognizing pathogens in dramatically contrasted environments to enhance the fitness for the adaptation and survival of cetaceans.

Background

Microbial pathogens (bacteria, fungi, protozoa, and viruses) affect plants and animals of the world dramatically, including their survival, growth, development, and reproduction. In response to pathogen invasion, multicellular organisms have evolved several distinct immune-recognition systems. Unlike the adaptive immune system only found in vertebrates, the innate immune system is a universal and evolutionarily ancient mechanism existing in all multicellular organisms [1]. The innate immune system nonspecifically recognizes and kills pathogens at the first time and at the first line. The targets of innate immune recognition are called pathogen-associated molecular patterns (PAMPs), produced only by microbes and shared by a class of microorganisms. PAMPs are highly conserved because such molecular patterns are essential to the integrity, function, or replication of microbes [2]. Accordingly, PAMPs are recognized by a variety of host receptors called pattern recognition receptors (PRRs).

Toll-like receptors (TLRs) are among the best characterized PRRs that lie directly at the host-pathogen interface. Although TLRs have been regarded for a long time as a classic example of strong evolutionary conservation and intense functional constraint [3,4], a recent comparison of several Drosophila genomes showed for the first time the fast evolution between closely related species [5]. Although this contradicts the traditional view regarding innate immunity, this finding is congruent with theoretical prediction that over evolutionary time TLRs may be engaged in co-evolutionary arms races with their microbial ligands. Some recent discoveries and characterization surveys of TLRs variation in vertebrates [5-7] provide further corroboration for this prediction. To date, however, very few studies have been conducted on the evolution of TLRs in a limited number of vertebrate species, including primates [3,8-10], ungulates [11], birds [12,13], and bony fishes [14]. Furthermore, the results from different studies are incongruent with or contradict each other. For example, although Ferrer-Admetlla et al. [6] regarded balancing selection as the best explanation for sequence variation at human TLRs, Mukherjee et al. [3] did not detect any effect of natural selection on TLRs of the Indian population and thus supported the traditional viewpoint that purifying selection is the major driving force for the evolution of TLRs. In some inter-specific studies, Ortiz et al. [15] detected positive selection at the TLRs of five primate species only, whereas Nakajima et al. [8] found the action of positive selection on TLR4 when they examined a more extensive phylogenetic sampling. Recently, Wlasiuk et al. [9] and Wlasiuk and Nachman [10] detected positive selection on most TLR loci of primates, but intra-specific polymorphisms were found to be influenced mainly by population demography rather than by adaptive evolution. In other words, they found that primate TLRs are characterized by a mode of episodic evolution. Positive selection and evolutionary constraint have also been detected in birds [13] and bony fishes [14], suggesting the role of adaptive evolution in response to changes of environmental pathogens. Considering the limited number of taxa and loci examined in these studies, a clear picture of the evolution of the TLR gene family has not been painted so far, and more data are necessary to resolve this problem.

Cetaceans, including whales, dolphins, and porpoises, are a group of secondarily adapted marine mammals with a history of transition from terrestrial (land) to full aquatic habitats and subsequent adaptive radiation in waters around the world. Although the exact origin and evolutionary history of extant cetaceans remains unclear, a widely accepted view is that the direct terrestrial ancestors of cetaceans (a group of mammals called artiodactyls [16,17]) returned to the sea around 50 MYA[18-21]. The ancient cetaceans evolved gradually to conquer nearly all oceans and some rivers of the world [22-24], and finally diversified into a group of fully aquatic mammals including nearly 85 extant species that can be subdivided into two suborders (Odontoceti and Mysticeti) [25-27]. During the transition from land to sea and the radiation and diversification into various aquatic environments, cetaceans must have been confronted with formidable challenges from ever-changing environmental pathogens. For this reason, cetaceans could be one of the most ideal candidate taxa for studying the evolutionary process and the associated driving mechanisms of vertebrate innate immune systems such as TLRs.

Here, TLR4 was used as an example to reveal the evolutionary history of pattern recognition molecules across cetaceans and their closest terrestrial relatives. TLR4 is expressed on the cell membrane and is mainly responsible for the recognition of lipopolysaccharides (LPS) from Gram-negative bacteria [28] and even components of yeast, Trypanosoma, and viruses [29]. This molecule interacts with LPS indirectly aided with myeloid differentiation factor 2 (MD-2) [30] through the formation of a duplex heterodimer (TLR4-MD-2-LPS)2 that is essential to activate a signaling pathway mediating the defense against Gram-negative bacteria. It has been reported that some substitutions in the changed amino acid residues of TLR4 can alter the interaction among TLR4, MD-2, and LPS, and modify the TLR4/MD-2 immunological responses [10,13]. In this study, the open reading frames (ORF) of TLR4 from representative cetaceans and some closely related artiodactylans were sequenced to elucidate whether this innate immune gene has been the target of positive selection in cetacean evolutionary history. The aims of this study were 1) to find evidence of positive selection at TLR4 in cetacean origin and evolution, and 2) to evaluate whether the evolutionary rate of TLR4 varied in different cetacean lineages, and if so, what factors could account for this evolutionary pattern. It was interesting to find compelling evidence of positive selection acting on TLR4 throughout cetacean evolution, from their origin till the present, and it was speculated that the species-specific effects and/or the complex interaction of multiple factors (abiotic and biotic) might have played a major role in driving the heterogeneity in the evolutionary rate of cetacean TLR4.

Results

In this study, the full sequences containing 2250 bp of TLR4 open reading frame (ORF) from 17 representative cetaceans and three even-toed ungulates were obtained, 12 of which were newly determined and have been deposited in GenBank with accession nos. JN642608-JN642619 (Additional file 1: Table S1). The Bayesian analyses and Neighbor-Joining (NJ) method yielded a similar topology (Figure 1), which is basically consistent with a widely accepted hypothesis of whale phylogeny [17,31-33]. This phylogeny was then used as the working topology in the subsequent analyses. To our knowledge, this is the first study thus far to characterize a TLR locus in cetaceans and to provide some novel insights into the evolution of the innate immune system in the cetacean clade.

Figure 1

Positive selection at TLR4 across the cetacean phylogeny. Branches a to p correspond to those in supplementary Table S2. The ω value calculated by the free-ratio model is labeled along each branch. In some cases, zero synonymous substitutions lead to a ω value of infinity (n.a.). The estimated numbers of nonsynonymous and synonymous changes are shown in parentheses. The branches in red show strong evidence of undergoing positive selection. Amino acid changes were estimated by parsimony method, and every substitution of these sites is marked in blue. Six clades in which amino acid substitution occurred are filled with six different colors. The parallel amino acid changes are listed on the right of the corresponding terminal branches, while b, c, h, and l in parentheses stand for the internal branches on which parallel changes occurred. Amino acid positions (numbers) and parallel changes at each position were listed in the right part of the figure1. A = even-toed ungulates, B = river dolphins, C = oceanic dolphins, D = porpoises and white whales, E = sperm whales, F = baleen whales.

Positive selection at cetacean TLR4

The site model incorporated in Phylogenetic Analysis by Maximum Likelihood (PAML) was used to reveal whether cetacean TLR4 was subjected to positive selection. We compared nested models and found that a model including sites with ω > 1 fitted the data significantly better than did a neutral model. Model M8 detected 25 (3.3%) sites under selection with the average ω value of 3.55 in cetacean (Table 1). The specific codons identified by the Bayes empirical Bayes (BEB) approach with a posterior probability of 90% constituted an even smaller fraction (11 codons, 1.5%). With the use of Datamonkey, 17 and 13 codons were detected by fixed effects likelihood (FEL) and random effects likelihood (REL), respectively, whereas no site was detected by single likelihood ancestor counting (SLAC). When all these analyses from PAML and Datamonkey were combined, nine codons (150, 179, 183, 207, 228, 247, 272, 280, and 324) were picked out as robust sites under positive selection by at least two Maximum Likelihood (ML) methods, five (179, 207, 228, 272, 280) of which were predicted by three ML methods. In general, the more radical the amino acid substitutions are, the more likely they will affect function during evolution [34]. Most of the nine codons identified under selection made relatively conservative changes, while sites 272 and 280 were involved in radical changes in their physicochemical properties (size, polarity, and electric charge). In particular, codon 280 showed the strongest evidence of selection not only because it was detected by three ML methods, but also because it showed radical changes in three independent lineages (Table 2).

Table 1

Tests for positive selection at cetacean TLR4 using branch model and site models

Model	Models Compared	-2lnΔL	df	p value	Proportion of Sites under Selection	ω (dN/dS) of Sites under Selection
Site model	M1 versus M2	16.10	2	< 0.0001	0.033	3.45
	M7 versus M8	16.48	2	< 0.0001	0.033	3.55
	M8a versus M8	16.10	1	< 0.0001
Branch model	M0 versus full	44.13	31	0.0595

Table 2

Positive selection at amino acid sites of cetacean TLR4

AA Positionsa	PAML Site Model (M8) p > 0.9	PAML Branch-Site Modelc	FELd p < 0.2	RELd BF > 50	AA Changes	Parallel Changes	Property Changese	Protein Domain	Functional Informationf	Cladeg
25					Ser-Arg	Yes	SM, P, NEU-P.POS	LRR3		B, C
28			0.07		Leu-Trp		NP, NEU-P, NEU	LRR3		A, D
					Leu-Pro		NP, NEU-SM, NP, NEU
45			0.16		Thr-Ile	Yes	SM, P, NEU-NP, NEU	LRR4		A, F
104			0.17		Leu-Val		NP, NEU-NP, NEU	LRR6
					Leu-Ser		NP, NEU-SM, P, NEU
128		0.875			Glu-Pro		P, NEG-SM, NP, NEU	LRR7		D
133		0.723			Asn-Lys		SM, P, NEU-P, POS	LRR7		G
139		0.708			Gly-Glu		SM, NP, NEU-P, NEG	LRR8	Adjacent to site involved in interaction with MD2	G
149		0.565			Ser/Leu -Thr		SM, P, NEU/NP, NEU -SM, P, NEU	LRR8		G
150	0.995			228.23	His-Arg	Yes	P, POS-P, POS	LRR8		A, B, C
					His-Asp		P, POS-SM, P, NEG
177				61.94	Asn-ThrAsn/Thr-IleAsn-LysIle-Asn		SM, P, NEU-SM, P, NEUSM, P, NEU/SM, P, NEU-NP, NEUSM, P, NEU-P, POSNP, NEU-SM, P, NEU	LRR9	Adjacent to site involved in ligand binding and interaction with MD2	A, C, G
179	0.992		0.07	647.96	Lys-GluGlu-GlnGlu/Lys-Gln	Yes	P, POS-P, NEGP, NEG-P, NEUP, NEG/P, POS-P, NEU	LRR9		A, C, F
183			0.12	51.06	Arg-SerArg-Thr		P, POS-SM, P, NEUP, POS-SM, P, NEU	LRR9		C, D
204					Glu-His	Yes	P, NEG-P, POS	LRR10		A, D
207	1.000		0.08	1563.58	Gly/Lys-ArgArg-LysArg-ThrLys-Arg		SM, NP, NEU/P, POS-P, POSP, POS-P, POSP,POS-SM, P, NEUP, POS-P, POS	LRR10		A, G, C, E
212					Leu-Val	Yes	P, POS-NP, NEU	LRR10		A, B
221			0.1		Val-Met	Yes	NP, NEU-NP, NEU	LRR11		C, D, F
228	0.994		0.15	544.88	Asp/Ser/Cys-AsnAsp-Asn		SM, P, NEG/SM, P, NEU/SM, NP, NEU-SM, P, NEUSM, P, NEG-SM, P, NEU	LRR11		A, G
230	0.978				Gly/Glu/Asp-ArgAsp-His		SM, NP, NEU/P, NEG/SM, P, NEG-P.POSSM, P, NEG-P, POS	LRR11		A, E
239				50.32	Cys-Tyr	Yes	SM, NP, NEU-P, NEU	LRR12		B, D, G
247			0.14	86.14	Ile-ThrThr-Ile		NP, NEU-SM, P, NEUSM, P, NEU-NP, NEU	LRR12	Adjacent to site involved in interaction with ligand binding	C, G
250b	0.936				Asp/Ala-LysAsp/Lys/Ala-AsnAsn-Lys		SM, P, NEG/SM, NP, NEU-P, POSSM, P, NEG/P, POS/SM, NP, NEU-SM, P, NEUSM, P, NEU-P, POS	LRR12	Ligand binding	A, E, G
265					Phe-Leu	Yes	NP, NEU-NP, NEU	LRR13		B, E
272	0.997		0.13	188.28	Gly/Asp-HisGly-HisHis-Gly		SM, NP, NEU/SM, P, NEG-P, POSSM, NP, NEU-P, POSP, POS-SM, NP, NEU	LRR13	Adjacent to site involved in interaction with ligand binding	A, C
280	0.952		0.18	191.07	Glu-AlaGln/Glu-Ala	Yes	P, NEG-SM, NP, NEUP, NEU/P, NEG-SM, NP, NEU	LRR13		A, B, E
302		0.624			His-Arg		P, POS-P.POS	LRR14		D
304				55.05	Asp-AsnAsn-Pro		SM, P, NEG-SM, P, NEUSM, P, NEU-SM, NP, NEU	LRR14		G
324	0.996			301.87	Asn-SerAsn-LysGly-Asn		SM, P, NEU-SM, P, NEUSM, P, NEU-P, POSSM, NP, NEU-SM, P, NEU	LRR15	Adjacent to site involved in interaction with ligand binding (hydrogen bond)	C, E, G
342				53.56	Asn-SerAsn/Ser-Thr		SM, P, NEU-SM, P, NEUSM, P, NEU/SM, P, NEU-SM, P, NEU	LRR16	Adjacent to site involved in interaction with ligand binding (hydrogen bond)	A
351				0.17	Ile/Ala-Val		NP, NEU/SM, NP, NEU-NP, NEU	LRR16	Adjacent to site involved in interaction with ligand binding (hydrophobic interaction)	G
368		0.576			Ile-Thr		NP, NEU-SM, P, NEU	LRR17	Adjacent to site involved in interaction with ligand binding (hydrophobic interaction)	G
404				0.08	Leu-Met		NP, NEU-NP, NEU	LRR18		C
408					Ile-Thr	Yes	NP, NEU-SM, P, NEU	LRR19		A, G
409				0.19	Leu/Ile/Phe-Val		NP, NEU/NP, NEU/NP, NEU-NP, NEU	LRR19		A
482				0.16	Ser/Trp-PhePhe/Ser/Trp-Leu		SM, P, NEU/P, NEU-NP, NEUNP, NEU/SM, P, NEU/P, NEU-NP, NEU	LRR22		A
542	0.903				Met-Thr	Yes	NP, NEU-SM, P, NEU	LRRCT		A, D
551	0.938				Ile-ValVal-Ile	Yes	NP, NEU-NP, NEUNP, NEU-NP, NEU	Transmembrane		B, F
559				0.16	Val-Ala		NP, NEU-SM, NP, NEU	Transmembrane		G
690		0.564			Arg-Gln		P, POS-P, NEU	TIR		D
740		0.790			Glu-Asp		P, NEG-SM, P, NEG	TIR		G
742		0.697			Asn-Arg		SM, P, NEU-P, POS	TIR		G
743				0.18	Gln-Glu		P, NEU-P, NEG	TIR		A, F

a Codons identified by more than one ML method were in bold and underlined.

b Site 250 in italic was mapped onto the 3D structure of TLR4, since it directly participates in binding of LPS to TLR4.

c Codons were identified by branch-site model in PAML. Details were in Materials and Methods and Additional file 2: Table S2.

d Codons were estimated in DATAMONKEY.

e SM, small; NP, nonpolar; P, polar; NEU, neutral; POS, positively charged; NEG, negatively charged.

f Codons were in the functional regions predicted by the three-dimensional structure in Shishido et al. 2010. LRR = Leucine-rich repeat, CT = C-terminal, TIR = cytoplasmic Toll/IL-1 receptor

g Amino acid substitutions occurred in the following clades: A = even-toed ungulates, B = river dolphins, C = oceanic dolpins, D = porpoises and white whales, E = sperm whales, F = baleen whales, G = more than one equally parsimonious reconstruction

The amino acid changes reconstructed by parsimony were distributed along 42% of examined cetartiodactylan branches or 46% of examined cetacean branches. Thirteen codons (25, 45, 150, 179, 204, 212, 221, 239, 265, 280, 408, 542, and 551) showed parallel amino acid changes (Table 2), which could be regarded as candidates under selection. These codons were scattered across the entire whale phylogeny (Figure 1), rather than accumulated in just some specific lineages.

The LRT tests based on the branch model suggested that the free-ratio model fitted the data better than did the one-ratio model (Table 1), indicating that dN/dS ratios were indeed different among lineages. The ω values along three branches were found to be greater than 1 with nearly significant statistical support (p = 0.0595): branch a leading to the last common ancestor of cetaceans and hippos (ω = 4.59), branch b leading to oceanic dolphins (ω = 1.33), and branch c leading to the last common ancestor of Phocoenidae (porpoises) + Monodontidae (white whales) (ω = 1.34) (Figure 1). For all the cetacean lineages examined, ω values ranged from 0.0001 to 1.34, with an average of 0.61 (Figure 1).

When we used the branch-site model to predict positive selection acting on each branch (Additional file 2: Table S2), two lineages were detected under positive selection because likelihood ratio test (LRT) tests suggested that model A fitted the data better than did model M1a along branches a (whale + hippo) (LRT of test 2 = 5.40, df = 1, p = 0.02) and d (beluga whale) (LRT of test 2 = 8.20, df = 1, p = 0.004) (Figure 1). Six and three codons were respectively detected under positive selection along these two branches (Additional file 2: Table S2). The BEB values of the positively selected sites along these two branches were not high (0.564