Deciphering genetic adaptations of Old World camels through comparative genomic analyses across all camelid species

Abstract

Old World camels exhibit many unique traits for adaptations to desert environments, distinct from New World camels within the same family. Here, we conducted a comparative genomic analysis of three Old World camel species relative to four New World camel species, as well as cattle, pigs, mice, and humans to search for genes showing positive selection and specific variations in Old World camels. We identified genes under positive selection, genes with specific small indels, and inactivated genes associated with the unique lipid metabolism and skin characteristics of Old World camels. Especially, through experimental validation, we confirmed the inactivation of TAS2R16, a gene encoding a prominent bitter taste receptor, which could enable Old World camels to consume certain toxic plants with a bitter taste to other mammals. These findings provide insights into the genetic basis underlying the distinctive energy metabolism, water-salt homeostasis, and dietary adaptations of Old World camels.

Keywords: Evolutionary biology; Evolutionary ecology; Zoology.

Graphical abstract

Figure 1

Construction of orthologous gene families of 11 species (A) Left: Phylogenetic tree of 11 species constructed by TimeTree, N0–N9 denoting node numbers. Right: The number of ortholog groups for each copy number pattern in 11 species. 0, 1 and x (>1) represent the copy number in an ortholog group. (B) The Venn diagram illustrates the overlap of ortholog groups among the representative species. (C) Schematic representation of PGLYRP4 and its surrounding genes in the genomes of three Old World camel species, alpacas, and humans. Framed in red is PGLYRP4, which is not found in Old World camels.

Figure 2

Genes subject to positive selection in Old World camels (A) 37 PSGs with a BEB probability of positively selected sites exceeding 0.85. In the lower panel, each dot represents a PSG identified with the branch-site model test, and the dotted line represents an FDR cutoff of 0.2. The size of each dot is the number of positively selected sites. PSGs involved in enriched GO terms are colored and highlighted. In the upper panel, each bar represents the p-value of the corresponding PSG from the drop-out test. (B) Schematic diagram of the pathway of HACD3 involved in fat metabolism. (C) Protein sequence alignment of HACD3 across multiple species, showing the evolutionary conservation and protein domains of the positively selected sites (red box) in Old World camels. (D) MutPred2 prediction of the substitution D270N for HACD3. The significance of structural changes is evaluated with empirical p-values. (E) MK test for HACD3. Divergence is estimated by comparing the ancestral sequence of all camelid species to the sequence of the Bactrian camel. Polymorphism data are derived from the genomes of 124 Bactrian camels. The p-value is calculated using the eMKT method.

Figure 3

Gene sets with specific indels in Old World camels (A) 38 genes exhibiting Old World camel-specific indels longer than one amino acid. A gene may contain multiple such indels. (B) Enrichment analysis of the indel gene set across various tissues and organs, conducted using FUMA. (C) Protein sequence alignment of HMGCR across multiple species. The evolutionary conservation of the deleted amino acids and their location within protein domains are shown in the red box. (D) MutPred-Indel prediction of the deletion for HMGCR. The significance of structural changes is evaluated with empirical p-values.

Figure 4

Gene sets with inactivating substitutions detected in Old World camels (A) Ten highly likely inactivated genes identified by both PseudoChecker and TOGA. Blue represents nonsense substitution, orange represents frameshift substitution, and the size of the dots represents the number of inactivating substitutions. (B) Left: The KEGG and GO terms significantly enriched with inactivated genes. Right: Genes involved in each term. (C) Potential inactivation events of the gene set across 423 mammalian species as annotated by TOGA. The dot color indicates the frequency of potential inactivation within each taxonomic group, while the dot size indicates the level of enrichment of potential inactivation within that group (∗ FDR <0.05, Fisher’s exact test). (D) Protein and coding sequence alignments of CES3, CPA3 and TAS2R16. The protein alignments indicate locations of the inactivating substitutions within the protein domain. The coding sequence alignments highlight these inactivating substitutions (red box) in three Old World camel species.

Figure 5

Experimental validation of TAS2R16 inactivation in Old World camels (A) The lower part shows a schematic representation of the seven-transmembrane GPCR encoded by TAS2R16. The upper part displays the ORF sequences of TAS2R16 in the human, New World camel, and Old World camel. Red highlights the frameshift deletion in the ORF1 of Old World camels. The gray shading represents the approximate location of the frameshift deletion in ORF1, and the alternative start codon of ORF2 within the GPCR structure. (B) Western blot analysis of TAS2R16 protein expression from the Bactrian camel, alpaca, and human. Two replicates were performed for each species. (C) Immunofluorescence staining of cells expressing TAS2R16 from the Bactrian camel, alpaca, and human. The nucleus is shown in blue, and the target protein in green. Images from left to right: nucleus, cell membrane, and merged view. The control group used an empty vector. Scale bar: 20 μm. (D) Dose-response curves for cells expressing TAS2R16 from the Bactrian camel, alpaca, and human to salicin. The salicin concentrations were 0.01, 0.1, 0.3, 0.6, 1, 3, and 10 mM, with the maximum response intensity observed at 10 mM. The error bars represent standard errors of estimation. Three replicates were performed for each group and concentration. Differences between the Bactrian camel and alpaca are statistically significant at all concentrations (∗∗p < 0.01, ∗∗∗p < 0.001, two-sided t-test).