Sequence Characterization of the MC1R Gene in Yak (Poephagus grunniens) Breeds with Different Coat Colors

Melanocortin 1 receptor (MC1R) gene plays a key role in determining coat color in several species, including the cattle. However, up to now there is no report regarding the MC1R gene and the potential association of its mutations with coat colors in yak (Poephagus grunniens). In this study, we sequenced the encoding region of the MC1R gene in three yak breeds with completely white (Tianzhu breed) or black coat color (Jiulong and Maiwa breeds). The predicted coding region of the yak MC1R gene resulted of 954 bp, the same to that of the wild-type cattle sequence, with >99% identity. None of the mutation events reported in cattle was found. Comparing the yak obtained sequences, five nucleotide substitutions were detected, which defined three haplotypes (EY1, EY2, and EY3). Of the five mutations, two, characterizing the EY1 haplotype, were nonsynonymous substitutions (c.340C>A and c.871G>A) causing amino acid changes located in the first extracellular loop (p.Q114K) and in the seventh transmembrane region (p.A291T). In silico prediction might indicate a functional effect of the latter substitution. However, all three haplotypes were present in the three yak breeds with relatively consistent frequency distribution, despite of their distinguished coat colors, which suggested that there was no across-breed association between haplotypes or genotypes and black/white phenotypes, at least in the investigated breeds. Other genes may be involved in affecting coat color in the analyzed yaks.


Introduction
The coat color of domestic animals is highly variable, and has been widely used as a unique phenotype in the morphological selection that resulted in the constitution of the modern breeds. In a variety of animals, pigmentation is mainly determined by the relative amount of two basic melanins, eumelanin (black/brown) and phaeomelanin (yellow/red), which are genetically controlled by the Extension (E) and Agouti (A) loci, respectively [1,2]. The wild alleles at the two loci would result in different coat colors, whereas dominant or recessive mutations can alter the eumelanin/phaeomelanin ratio and subsequently cause coat color variations [3]. The E locus encodes the melanocortin receptor 1 (MC1R), which is a seven-transmembrane G-proteincoupled receptor expressed in melanocytes [4]. Mutations affecting the MC1R gene fvunction determine different coat color phenotypes in a large number of vertebrates, including livestock species such as cattle [5][6][7], sheep [8], goat [9], horse [10], pig [11,12], rabbit [13], and chicken [14]. Some MC1R mutations induce eumelanin production while others prompt phaeomelanin synthesis [4].
The cattle MC1R gene is located on chromosome 18 [15] and has four principal alleles (E + , E D , E 1 , and e), which are defined by three mutation events affecting the coding region [5,7]. The wild type allele (E + ) produces a variety of coat colors. The dominant allele (E D ), differing from the E + allele by the g.296T>C substitution, is associated with black color. The recessive allele (e) is characterized by a frameshift mutation at position 310 and encodes a truncated protein with only three putative transmembrane domains. Animals homozygous for this allele show red/white coat color phenotypes across a large number of cattle breeds [5,6,16,17]. The E 1 allele is characterized by a 12 bp duplicated fragment at position 669 that was detected in a few breeds which have fawn to brown or grey coat colors 2 Journal of Biomedicine and Biotechnology  [7,16,18]. Other missense and frameshift mutations have been identified in the cattle MC1R gene, but their functional effect has not been completely investigated yet [16,19].
The yak (Poephagus grunniens) belongs to the family of Bovidae and lives in an extremely harsh habitat (>3000 m above the sea level). More than 90% of worldwide domestic yaks are distributed in the Qinghai-Tibet Plateau region in China. There are five registered yak breeds (Jiulong Yak, Qinghai Yak, Tianzhu Yak, Maiwa Yak, and Tibet Yak) in China as well as other local populations [20]. The Tianzhu yak is wellknown for the uniform white coat color even if a few subpopulations in this breed show small dark spots. The black coat is predominant in the other four yak breeds. More than 85% individuals of Jiulong Yak have uniform black coat, and the remaining individuals have a black coat with white spots/belt. The genetic basis of coat color in yaks is poorly studied compared to the cattle [1,21]. Chen et al. [22] investigated coat color inheritance in yaks by means of classical genetics and suggested that the black and yellow phenotypes were mainly determined by two loci (B and D), but the results of this report were not considered in the context of the comparative approach proposed by Searle [1]. In addition, no study has analyzed pigment or melanocyte distribution in yak, yet. Therefore it is not completely clear which could be the possible candidate genes involved in the black and yellow/white phenotype determination in yak. Here, in order to evaluate if mutations in the MC1R gene could potentially be involved in determining coat color variations in Poephagus grunniens, we characterized by sequencing this gene in 44 yaks belonging to three breeds with different coat colors: Tianzhu yaks (pure white animals), Jiulong yaks, and Maiwa yaks (black animals) ( Figure 1).

Sampling.
Ear tissue samples were collected from 44 yak individuals belonging to three pure indigenous breeds: 14 Tianzhu yaks, 16 Maiwa yaks, and 14 Jiulong yaks. The Chinese provinces and the number of herds from which the animals were sampled as well as their coat color appearances are reported in Table 1. Figure 1 shows yaks of white or black coat color types. Biological material was stored at −70 • C.

Amplification of the MC1R Gene and
Sequencing. Genomic DNA was extracted by a standard phenol/chloroform method. The whole coding region of MC1R gene was amplified using a primer pair (forward: 5 -GGACCC-TGAGAGCAAGCAC-3 ; reverse 5 -CTCACCTTCAGGGAT-GGTCTA-3 ) designed using the cattle sequence (GenBank accession number NM 174108). The amplified fragment length was 1085 bp. Amplification of MC1R gene was performed in a 50 µL reaction mixture containing 100 ng of DNA, 10 mM Tris-HCl (pH 8.3), 2.5 mM MgCl 2 , 50 mM KCl, 10 pmol of each primer, and 1 U of Taq polymerase (Takara, Dalian, China) using the following procedures: after a denaturation cycle at 94 • C for 5 minutes, amplification was performed for 35 cycles at 94 • C for 45 second, 61.5 • C for 45 second, and 72 • C for 1 minute, and ended with a final incubation at 72 • C for 10 minutes for full extension.
All PCR products were purified on spin columns and were directly sequenced by using the amplification primers and the Big Dye Terminator v. 3.1 Cycle Sequencing Kit (Applied Biosystems, Calif, USA) on an ABI Prism 3100 DNA sequencer (Applied Biosystems) according to the manufacturer's manual. Sequencing electropherograms were visualized using the SeqMan program (DNAStar Inc., Madison, Wis, USA).
When two or more heterozygous positions were recognized in the same sequencing read, in order to clarify the haplotypes, the PCR products obtained from one individual for each heterozygous pattern were cloned after separation on 1.5% agarose gel followed by purification with the Gel DNA Extraction Kit (Tiangen, Beijing, China). The target fragments were ligated to the PMD20-T vector (Takara, Dalian, China) and were transformed into E. coli strain JM109 competent cells. Seven or eight positive clones were further selected, and plasmids were isolated for insert sequencing as described above.

Sequence Data Analysis.
After editing of the sequences, the noncoding regions were truncated, and the predicted coding region was aligned with the corresponding region of the cattle MC1R wild-type allele (E + ) (GenBank accession numbers Y13957 and Y19103; [23,24]). The relationship among the yak alleles, defined according to their haplotypes, was evaluated using program Network 4.1 (http://www.fluxus-engineering.com/sharenet.htm), along with the E + , E D , and e cattle MC1R alleles. The putative transmembrane structure of the MC1R protein, conducted by the TMpred program [25], was drawn to show the position of the amino acid substitutions.  In silico functional analysis of missense mutations was obtained using PANTHER [26]. PANTHER estimates the likelihood of a particular amino acid substitution to cause a functional impact on the protein. It calculates the substitution position-specific evolutionary conservation (subPSEC) score based on an alignment of evolutionarily related proteins [26][27][28]. The probability that a given variant will cause a deleterious effect on protein function is estimated by P deleterious , such that a subPSEC score of −3 corresponds to a P deleterious of 0.5 [29]. The subPSEC score is the negative logarithm of the probability ratio of the wild-type and mutant amino acids at a particular position. PANTHER subPSEC scores are continuous values from 0 (neutral) to about −10 (most likely to be deleterious).

Results and Discussion
Very limited information is available on coat color genetics in yak despite the presence of coat color variability among and, in some cases, within breeds. In addition, to our knowledge, before our study, no information about any coat color gene was available in yak. Therefore, in order to add data in this field and try to link coat color phenotype variations in this species with polymorphisms at the DNA level, we first characterized by sequencing the MC1R gene, whose mutations in cattle and several other species are responsible for eumelanism and pheomelanism. The predicted coding region of the yak MC1R gene resulted of 954 bp, the same to that of the wild type cattle sequence [5][6][7], with identity >99%. Comparing the 44 yak sequences, a total of 5 single nucleotide polymorphisms (SNPs) were identified. These SNPs defined three haplotypes that were indicated as E Y 1 , E Y 2 , and E Y 3 . Their sequences were submitted to GenBank database with the accession numbers FJ624478, FJ624479, and FJ624480, respectively. Considering as reference the wild type cattle sequence, the E Y 1 haplotype was characterized by the c.340C>A and c.871G>A substitutions, whereas the E Y 3 haplotype included the c.375C>T and c.663T>C mutations (Figure 2). The E Y 2 haplotype differed from the E Y 3 haplotype by c.642C>G substitution (Figure 2). None of the mutation events reported in cattle was found in the sequenced yaks ( Figure 2 and data not shown). It is worth to note that in river buffalo (Bubalus bubalis) Klungland et al. [24] identified the MC1Re allele with quite high frequency (0.20), suggesting a recurrent origin of the deletion causing this allele or a mutation event that occurred before the radiation of the Bovinae subfamily.
The relationship of the three major cattle MC1R alleles with the yak haplotypes is shown in Figure 3. The E Y 1 yak haplotype has four-and five-mutation distance with the E Y 2 and E Y 3 haplotypes, respectively, passing through the differentiation node constituted by the wild-type cattle sequence. Therefore, it could be expected that another haplotype, having the same sequence of the E + cattle allele, 4 Journal of Biomedicine and Biotechnology Figure 4: Hypothetical structure of the yak MC1R protein with the two amino acid changes identified in this study. Amino acid changes are scored relative to the E Y 1 haplotype. Transmembrane domains were predicted using the TMpred program [25], with reference to the reported cattle MC1R structure. might be present (or was present) in some other yak populations. This network distribution of the yak haplotypes resembles the mtDNA data that suggested that the Chinese yak was domesticated from two distinct matrilineal sources or from a heterogeneous pool containing two divergent lineages, with occasional gene introgression from cattle [30].
Of the five identified yak SNPs, three (c.375C>T, c.642C>G, and c.663T>C) were predicted as synonymous mutations while the two polymorphisms characterizing the E Y 1 haplotype (c.340C>A and c.871G>A) were deduced as nonsynonymous substitutions causing a p.Q114K and a p.A291T amino acid change, respectively (Figures 2 and 4). These amino acid changes resided in the first extracellular loop region (p.Q114K) and in the seventh transmembrane domain (p.A291T) (Figure 4). In these regions, several amino acid substitutions, a few of which with effects on pigmentation, have been detected in other species, but none of the positions that are polymorphic in yak has been associated with any coat color variation. However, it is interesting to note that in the same transmembrane domain harboring the yak p.A291T amino acid change, a substitution at position 298 (p.Y298C) in black bear is associated with a white coat color variant [31]. Moreover, at position 294, a p.D294H substitution is associated with red hair and fair skin in human [32]. In silico functional analysis of the p.Q114K substitution did not support a putative functional role of this mutation (subSPEC = −2.60994; P deleterious = 0.4037) while a functional effects was more likely for the p.A291T amino acid change (subSPEC = −4.64874; P deleterious = 0.83872). Therefore, it could be important to confirm the putative predicted functional effect of the latter mutation with specific experiments evaluating the pharmacological activity of this variant.
The distribution of the MC1R haplotypes among the three investigated breeds is shown in Table 1. In Tianzhu yaks and Maiwa yaks, E Y 1 haplotype was predominant (0.68 and 0.75, respectively), followed by haplotype E Y 3 in Tianzhu animals (0.21) or by haplotype E Y 2 in Maiwa breed (0.22), whereas in Jiulong there was a relatively equal distribution of the three haplotypes. All six possible genotypes (derived by the combination of the three haplotypes) were identified in Jiulong with relatively even frequency distribution, while five and three genotypes were found in Tianzhu and Maiwa yaks, respectively. The E Y 2 /E Y 3 genotype, which was absent from Tianzhu yaks and Maiwa yaks, had the highest frequency in Jiulong yaks (0.43). When we tentatively grouped the yak samples according to coat color, white yaks had a higher frequency of the haplotype E Y 1 (0.68 versus 0.53) and a lower frequency of haplotype E Y 2 (0.11 versus 0.28) than black yaks. From these data it seems evident that there is no association between across-breed coat color phenotype and any haplotype or haplotype combination (genotype). However, the presence of possible regulatory mutation(s) in noncoding regions, affecting MC1R gene expression and, ultimately, yak coat colors (black versus white), was not directly excluded in this study, as we investigated only the coding region and partial 5 -and 3 -untranslated regions. Nonetheless, to account for possible regulatory mutation(s) we should hypothesize its/their occurrence only in all three haplotypes identified in one coat color type but not in the same three haplotypes observed in the other coat color type, making this possibility improbable. Therefore, based on the available data, the Extension locus in yak might not play an active role in black versus white coat color variation across breeds, at least considering the three investigated breeds.
Following these results, it remains to be clarified the nature of the white coat color in Tianzhu yaks as at present there is no study that could clarify if this phenotype is due to pheomelanism or absence of melanocytes or other mechanisms. Detailed hair and skin phenotypical evaluations and Journal of Biomedicine and Biotechnology 5 histological analyses as well as cross-breeding experiments might help in evaluating this issue.
Coat color can be considered a complex trait for which multiple genes, pathways, or networks are involved. Therefore, many genetic factors, other than the MC1R variants identified in the present study, might be more directly linked to coat color phenotypes in yak. More studies are essential to draw a complete picture of coat color determination in this species.
Other information could be extracted from the reported results. A larger level of variability seems present in the Jiulong yaks compared to the other sampled breeds. Moreover, the Maiwa yak and Jiulong yak, all with black coat color phenotype, had a significantly different distribution (P < .05, chi square test) of the E Y 1 haplotype. One possible explanation for this pattern is that the Jiulong yaks might have a different population structure and/or origin even if mtDNA data did not support any conclusion in this direction [30]. Further studies might be needed to confirm the differences of nuclear variability and, eventually, clarify the reasons of these results.