Association of SNP Haplotypes at the Myostatin Gene with Muscular Hypertrophy in Sheep *

The myostatin gene of seven important meat (Beltex (Australia), Beltex×Huyang (F1), Meat and Multi-Prolific Chinese Merino Fine Wool, Meat Chinese Merino Fine Wool and Dorper (South Africa)) and non-meat (Huyang and Kazak) sheep breeds was analyzed to study the genetic basis of muscular hypertrophy (double muscling) phenotype in sheep. SNPs, four in regulatory regions and several in the introns in the myostatin gene, were identified, and the former four SNPs were used for further studies. Twelve haplotypes were predicted by PHASE program, of which four main haplotypes (1, 3, 7, 9) were present in 90% of the 364 sheep in the study. Haplotypes 1-4 were mainly present in meat breeds while haplotypes 7 and 9 dominated the non-meat breeds. The association between haplotypes and average daily gain (ADG) was analyzed among 116 sheep with production data, Haplo2 (CGAA) and Haplo8 (TGAA) were identified to have significant (p<0.05) effect on ADG by the model (JMP5.1 software) taking into account the effects of breed, family background, haplotype, birth weight and sex. ADG of these haplotype groups also correlated well (r = 0.82) with hypertrophic phenotype scores. In conclusion, the mutations -956 (T→C), -41 (C→A) and 6223 (G→A) involved in Haplo2 and 8 may be associated with the double-muscling trait by influencing myostatin function and be suitable markers in selecting meat sheep. (


INTRODUCTION
The explanation for double-muscled animals has puzzled breeders since its identification in the 1880s.Only recently, developments in modern genetics have allowed the genetic elements underlying the phenomenon to be revealed.In 1997, the myostatin gene (MSTN) was first demonstrated in mice to be a potent negative regulator of the growth and development of skeletal muscle, and loss of myostatin in mice lead to a dramatic increase in skeletal muscle mass (McPherron et al., 1997).In cattle, mutations and deletion(s) identified in the coding region of the myostatin gene of two well-known double-muscled cattle breeds Belgian Blue and Piedmontese were found to be highly associated with the double muscling trait (Kambadur et al., 1997;McPherron et al., 1997).In pigs, mutations in the myostatin promoter were found to be associated with higher average daily gain (ADG) of the double-muscled pig in contrast to regular breeds (Jiang et al., 2002).In chickens, mutations of the myostatin gene were also reported to be associated with production traits (Gu et al., 2002).In goats, four haplotypes were identified in the intron 2 of the myostatin gene in 17 Chinese indigenous goat breeds (Li et al., 2006a), and myostatin genotypes were associated with differences in body weight (Li et al., 2006b).In dogs, a two-base deletion in a myostatin exon was recently found to result in premature stopping of myostatin translation and dogs with a heterozyote genotype had an enhanced muscle development and racing performance (Mosher et al., 2007).These observations indicate that the myostatin gene is an important candidate gene associated with double-muscling in animals.
For sheep, however, a convincing functional connection of myostatin variation with double-muscling has remained elusive.The whole coding region (CDS) of the myostatin gene in Texel double-muscling sheep was sequenced (Marcq et al., 1998), but no sequence differences were found in the coding region.However, a strong quantitative trait locus (QTL) for the double muscling trait was found at the OAR2 locus in the region where the myostatin gene is located (Walling et al., 2004).This evidence implies that mutations in non-coding regions of the sheep myostatin gene may also be associated with the double muscling trait.Mutation(s) in the myostatin 3'UTR have been found at the molecular or cellular level in a recent study by Clop et al. (2006).The mutations create an illegitimate miRNA binding site which might affect the double muscling trait of the Texel sheep.In addition, mutations in the promoter regions of the swine myostatin were reported to affect the stability of myostatin mRNA, and therefore affect animal phenotype (Cieślak et al., 2003).
Based on the above facts, we carried out a study to systematically examine polymorphisms in the promoter, coding region, introns and UTRs of the myostatin gene in seven important production sheep breeds, and sought to identify the associations between these SNPs and ADG among different sheep breeds.To our knowledge, so far there has been no study performed on the association between myostatin polymorphisms and the double muscling of sheep.

Animals and samples
Experimental sheep breeds were selected based on meat production characteristics (meat and non-meat).A total of 364 sheep from seven pure and crossed breeds were analyzed, including Meat and Multi-Prolific Chinese Merino Fine Wool (MPMF), Meat Chinese Merino Fine Wool (MF), Huyang, Kazak, Beltex, Dorper and B×H (a crossbreed of Beltex and Huyang).Huyang and Kazak are classified as non-meat breeds and MPMF, MF, Beltex, B×H and Dorper are typical meat (cross)breeds.In particular, the double-muscling Beltex is known for its exceptional ability in meat production.Moreover, the sheep breeds were ranked by the degree of muscular hypertrophy phenotype using a four category scale from 0 (no hypertrophic sign) to 3 (high double muscling trait).The double-muscling sheep Beltex and B×H were scored 3; the non-meat breeds Huyang and Kazak were scored 0; Dorper and MF were scored 2 and MPMF was scored 1.The score of each haplotype group was calculated using the sum of scores of haplotype group divided by the number of the sheep.
Blood and ear tissue were collected from all the breeds.Skeletal muscle and heart muscle were collected from a small number of Beltex and Huyang sheep in a manner which complied with International Guiding Principles for Animal Research.Samples for DNA extraction were stored at -20°C, and samples for RNA extraction were stored in liquid nitrogen.

Isolation of genomic DNA and total RNA
Genomic DNA was extracted according to the phenolchloroform method and dissolved in 50 μl TE buffer (pH 8.0).Total RNA was extracted from muscle samples by the Trizol method according to the manufacturer's instructions (Tiangen RNA extract kit) and each preparation was dissolved in 30 μl of RNase-free water.

Identification of SNPs
PCR-SSCP method (PCR-Single Strand Conformation Polymorphism) (Chu et al., 2007) was used to detect the polymorphisms in each individual sheep with the primers listed in Table 1.For PCR-SSCP, 2 μl of PCR product triggered by each pair of primers was mixed with 8 μl loading buffer (95% formamide, 0.5 EDTA pH 8.0 and 0.005% bromophenol blue); the mixture was denatured at 98°C for 10 min, chilled on ice for 5 min, and loaded on polyacrylamide gel of proper concentration (8%-12%) according to the size of PCR products.Electrorphoresis was performed in 1×TBE buffer at 120 V for 12 h at 20°C.The gels were visualized after Silver Staining.The genotypes were determined according to different band types on SSCP gels.SNPs identified by SSCP were further confirmed by direct sequencing of the PCR products.

Haplotype and statistical analyses
Haplotype analysis has higher precision and quality than single nucleotide polymorphism analysis when analyzing SNP data sets (Rieder et al., 1999;Stephens et al., 2001).Haplotypes were predicted for the sheep populations in this study by a PHASE (v2.1) program (Stephens et al., 2003;Scheet et al., 2006), which can be downloaded freely from http://www.stat.washington.edu/stephens/home.html.
In order to estimate the association between haplotype and degree of muscular hypertrophy, we established a general linear model to analyze fixed effects of breed, family within breed, haplotypes within family within breed, birth weight, sex and age on ADG (from birth to six month).Using the Fit Model Test, all effects in our model were analyzed with statistical software JMP5.1 (SAS Institute Inc, 2000) simultaneously.The effect in a model is tested for significance by comparing the sum of squared residuals to the sum of squared residuals of the model with that effect removed.Residual errors that are much smaller when the effect is included in the model confirm that the effect has a significant contribution.The distance from a point to the line of fit shows the actual residual.The distance from the point to the horizontal line of the mean shows what the residual error would be without the effect in the model.The statistical model was designed as follows: Y opqmn = μ+Mstn o +Breed p +Family r +BW q +Sex m +Age n +ε opqmn, where Y opqmn is ADG value of the analyzed trait; μ is mean value in a population; Mstn o represents fixed effect of different haplotype within family within breed; Breed p represents fixed effect of breeds; Family r means fixed effect of family within breed; BW q is fixed effect of birth weight; Sex m and Age n represent fixed effect of sex and age(the same age); ε opqmn represents random error effect.

Sequence analysis and SNPs identification
A genomic region (4,826 bp) containing partial coding sequence of sheep (Beltex) myostatin was obtained by long range PCR.The sequence was submitted to Genbank (DQ990914).The size of intron I and intron II was identified as 1,833 bp and 2,030 bp, respectively, by aligning the new sequence with a previously published sequence AF019622 (McPherron et al., 1997).The sheep myostatin was reported to have a putative 3'UTR region (Clop et al., 2006), but the exact sequence has not been identified.In this study, sheep myostatin 3'UTR sequence was successfully cloned by 3'RACE based on the bovine myostatin 3'UTR sequence.The complete 3'UTR was identified as 1,575 bp by aligning myostatin genomic sequence DQ530260 (Clop et al., 2006) with AF019622 and the newly cloned sequence EF439841 (submitted).
Genotyping by PCR-SSCP was performed on selected genomic regions (the promoter region, the 5'UTR, the coding region, part of the introns and the 3'UTR) of sheep breeds.Altogether 22 SNPs were detected: two in the promoter region, one in the 5'UTR, none in the coding region, 18 in the introns and one in the 3'UTR.Additionally, 29 new SNPs were identified in two introns from sequencing results (positions identified as "N" by CHROMAS software, Technelysium Pty Ltd., Australia) of the PCR products obtained using mixed genomic DNA from different breeds as the template.The detailed SNP sites (referencing the A in the translation initiation codon ATG as +1) are displayed in Table 2. Of these 51 SNPs, 41 were confirmed to be new with respect to the published sequence (DQ530260).The density of SNPs was about 1 per 100 bases, which was similar to that of bovine myostatin in beef cattle (Dunner et al., 2003).

Alleles and haplotypes distribution
Since no polymorphism was found in the myostatin coding region in our study, polymorphism analysis focused on the four SNPs located in the regulatory regions (-956 (T→C), -781 (G→A) in promoter, -41 (C→A) in 5'UTR and 6223 (G→A) in 3'UTR).The data on allele frequency are summarized in Table 3.The most remarkable allele distribution was in the bi-allelelic (6223 (G→A)), which had frequencies of 86% and 55% for the A allele in the double-muscling sheep breeds Beltex and B×H, respectively, but the same allele did not exist in the other less hypertrophic breeds.The next noteworthy distribution was that of the C allele at the -956 (T→C) locus in the meat breeds (Beltex, B×H, Dorper, MPMF, and MF) (ranging from 63% to 85%).The C allele was much rarer in the nonmeat breeds (Huyang and Kazak, 35% and 24% respectively).For the -41 (C→A) locus, the A allele was more common in the five meat breeds (ranging from 54% to 89%) than in the non-meat breeds (Huyang and Kazak, 46% and 50%, respectively).The smallest difference was observed at the bi-allelelic (-781 (G→A)) where the G frequency was 100% in all breeds except in the MPMF breed which had 89% G and 11% A. Based on the four SNPs, twelve haplotypes were expected in seven sheep breeds or crosses by PHASE (v2.1) (Figure 1).The predicted frequencies were compared with the actual haplotype frequencies in the sheep population examined (Table 4) and were found to be very close to each other.
According to the result of haplotype analysis, the haplotypes Haplo 1, Haplo 3, Haplo 7 and Haplo 9 were the four main haplotypes, accounting for almost 90% of haplotypes found in the studied sheep populations, and Haplo 1 (CGAG) accounted for more than 30% (Table 3).
In contrast, the rare Haplo 10 (0.1% predicted by PHASE) was not found in the populations.Haplo2, 4 and 8 were only detected in double-muscling sheep (Beltex and B×H), and Haplo1, 2, 3 and 4 were prevalent in meat breeds (including double muscled breed or cross).The frequencies of Haplo1, 2, 3 and 4 combined represented approximately 85%, 72%, 76% , 69% and 63%, respectively, of the Beltex, B×H, Dorper, MPMF and MF flocks, but were found in less than 35% of the non-meat breeds (Huyang and Kazak).On the contrary, the pooled percentage of Haplo 7 and Haplo 9 were much higher in the two non-meat sheep breeds (65% and 76% for Huyang and Kazak, respectively) than in the five meat breeds (range from 14% to 37%).
Associations of ADG with haplotypes and phenotypic score 116 out of 364 samples with production records were analyzed, including 28 MPMF, 23 MF, 30 Hu yang, 28 B×H and 7 Beltex sheep, the four SNPs (-956 (T→C), -781   , 2, 3, 5, 6, 7, 8 and 9) in our study (Table 5), and the ADG of different haplotypes was calculated from the birth and 6-month weight.The effect of haplotype was evaluated to have significant (p = 0.0463) effect on ADG in the fixed model considering the effects of breed, birth weight, the family background and sex on ADG (Figure 3).The results of association analyses between haplotypes and ADG show that the ADG of Haplo2 (232.61 ±SE = 4.85) g/d) and 8 (233.28±0.00)were significantly higher than that of other haplotypes.The Haplo2 (CAGAA) and 8 (TGAA) were mainly comprised of the three transitions (T→C at position -956, C→A at position -41 and G→A at position 6223).The mutant allele C (SNP-956 (T→C)) and the mutant allele A (SNP -41 (C→A), SNP (G→A)) occupied almost 97% and 100% respectively in Haplo 2 and 8 while no mutant allele A (SNP -781 (G→A)) was found in them.
The association of the double-muscling phenotype and sheep haplotypes was further analysed for phenotypic score.A strong correlation (r = 0.8172) was found between phenotypic score and daily weight gain (Table 5, Figure 2).The general trend was that the higher the ADG, the higher the phenotypic score.The double-muscling haplotypes (Haplo2, 8) scored the highest (3.0) in all eight haplotypes, and their ADGs were highest.Moreover, the haplotypes 2 and 8 belonged exclusively to the double-muscling flocks, while the other haplotypes with low phenotypic score were present mainly in the breeds with less muscular hypertrophy.

DISCUSSION
Myostatin was considered as the candidate gene for double-muscling animals because mutations in its coding region are responsible for double-muscled Belgian Blue and Piedmontese cattle (Kambadur et al., 1997).The fact that no variation in the coding region was found in doublemuscling sheep myostatin has confused many researchers (Marcq et al., 1998), but it does not exclude the functional involvement of myostatin in double-muscling sheep because polymorphisms in other regions of myostatin could be equally important.Mutations in regions upstream of initial codon ATG may control the transcript level of mRNA   by increasing or destroying the number of transcription factor binding sites; and mutations in 3'UTR may create microRNA targets, whose binding to some well-paired miRNAs might suppress or block myostatin mRNA translation and lead to differences in sheep phenotypes.
In our study, based on sequence analysis of the noncoded region of the sheep myostatin gene, four SNPs were identified in the 5' promoter region, 5'-UTR and 3'-UTR, and dozens of new or previously reported SNPs were identified in the introns.The main emphasis of our research has shifted to analyze the SNPs in the regulatory region (putative promoter and UTRs) because they are most likely to influence the meat traits by modulating the level of MSTN expression.The two SNPs (-956 (T→C), -781 (G→A)) were identified in meat sheep for the first time; the SNP (-41C→A) and the SNP (6223 G→A) have been reported by Clop et al. (2006).However, until now, no work has been carried out to study the role of the three SNPs in meat sheep except for the SNP (6223 G→A).In this study, the association of the haplotypes consisting of the former four SNPs with sheep ADG was further analyzed, and the Haplo2 (CGAA) and 8 (TGAA), containing three main base transitions (a T→C at position -956, a C→A at position -41 a G→A at position 6223), were identified to be strongly associated with muscular hypertrophy traits.The high correlation of the haplotypes with breed (in the form of phenotype score here) is present in our sheep flocks, which may imply that the three mutations (-956 (T→C), -41 (C→A) and 6223 (G→A)) make a concerted contribution to control double muscling phenotype in sheep.Compared to SNP 6223, SNP -41 and SNP -956, the mutation of SNP -781 existing only in MPMF sheep breed (common meat breed) was not found in double muscling haplotypes (Haplo2 and 8), showing insignificant association with ADG.Hence there is no evidence that it works in meat sheep.
SNP(6223 (G→A)), which was identified as a QTN (quantitative trait nucleotide) located in myostatin 3'UTR and reported in Texel sheep by Clop et al. (2006), was also found in this study in Beltex and B×H sheep flocks with a high percentage but was not found in the other five breeds.Beltex and Texel sheep are both renowned for their extraordinary ability for meat production (Banks, 1997;Busboom et al., 1999).In the study by Clop et al. (2006), a G to A transition found in the 3'UTR of myostatin was verified to create a target site for miR-1 and miR-206 which are highly expressed in skeletal muscle.The miRNA binding results in translational inhibition of myostatin mRNA and thus in the reduction of myostatin protein level, which weakens the inhibition of myostatin protein in skeletal muscle development and growth (Clop et al., 2006).This may provide the biological mechanism for doublemuscling sheep.
In pigs, mutations in the promoter region and 5'UTR were reported to affect the stability of myostatin mRNA and influence the weight of fattening pigs to some extent (Jiang et al., 2002;Cieślak et al., 2003).The two SNPs (-41 and -956) maybe act similarly in meat sheep.To discover the potential biological function of the mutations in the promoter and 5'UTR, the SNPs at position -41 and -956 were analyzed with an online transcription factor prediction software (http://www.cbrc.jp/research/db/TFSEARCH.html).
The results showed that the SNP (-41 (C→A)) could create a putative transcription binding factor HSF (heat shock factor) which may have a role in development and differentiation; the novel SNP (-956 (T→C)) might destroy a putative C\EBPb (CCAAT) enhancer binding protein site which may serve in fat deposition.These suggestions are in line with the fact that myostatin knock-out mice are characterized by less fat and more muscle (McPherron et al., 1997).These, however, are conjectures, and the roles and mechanisms of these three mutations should be studied in future work.Furthermore, a strong correlation (r = 0.8172) between ADG and phenotypic score of sheep grouped by haplotypes provided some support for a possible role of the SNPs other than the 6223 SNP.Haplotypes with the highest ADG (Haplo 2 and 8) had the highest phenotypic score (3.0), indicating a perfect (100%) association between Haplo 2 and 8 with the double-muscling phenotype trait.This is easy to understand and also quite remarkable since the two breeds with double-muscling phenotype (Beltex and BxH) were hardly occupied by Haplo2 and 8. Haplotypes with relatively high ADG also had relatively high phenotypic scores, suggesting that the SNPs in the promoter region and 5'UTR might also be taken as meaningful genetic markers of meat sheep.
Combining the analysis of the SNPs, haplotypes, phenotypic scores and ADG in our study, the conclusions were that the 6223 G→A SNP and -41 C→A, which involved Haplo2 and 8, were strongly associated with double muscling and may be responsible for the trait of muscular hypertrophy.The SNPs (-956 (T→C), -781 (G→A)) might act as meaningful markers in selecting meat sheep breeds and their effects should be confirmed in further studies.Moreover, our study contributed to understanding of the biological functions of myostatin polymorphisms at the level of quantitative genetics and to future breeding applications of the myostatin gene.

Figure 1 .
Figure 1.Haplotypes frequencies of the four SNPs in sheep myostatin gene.The unfilled blocks represent frequency (%) predicted by PHASE software, and the filled blocks represent actual frequency (%) of the haplotypes in the study populations.

Figure 2 .
Figure 2. Correlation between ADG and phenotypic score of the haplotype groups.

Figure 3 .
Figure 3. Leverage Plot of the haplotype effect in fixed model.Leverage plot showing the effect of haplotype from the JMP5.1 analysis.The effect of haplotype is significant (p = 0.0463).

Table 1 .
Part of primers for SSCP and sequencing for the sheep myostatin gene

Table 4 .
Myostatin haplotype distribution in different breeds (cross) "/" means haplo10 was not present in the studied sheep populations, and animal numbers are in brackets after the breed.Values are given in percentages.

Table 5 .
Mean and standard error of daily gain (g/d) in different haplotypes **" indicate significant (p<0.01)difference of ADG between the indicated haplotype and the unmarked haplotype(s), statistical comparison was by t test.The phenotypic scores is calculated as the sum of individual phenotypic scores divided by the total number of sheep in that group.Unless indicated, 116 sheep with production data (28 MPMF, 23 MF, 30 Huyang, 28 B×H and 7 Beltex) were pooled and then grouped according only to the specified phenotype. "