Microsatellite polymorphism and its association with body weight and selected morphometrics of farm red fox (Vulpes vulpes L.)

Polymorphism of 30 canine-derived microsatellites was studied in a group of 200 red foxes kept on 2 Polish farms. 22 out of 30 microsatellites were selected to study association between marker genotypes and body weight (BW), body length (BL), body circumference (BC), tail length (TL), ear height (EH), length of the right front limb (FRLL), length of the right rear limb (RRLL), length of the right front foot (FRFL) and length of the right rear foot (RRFL). A total of 112 alleles and 243 genotypes were found at 22 autosomal microsatellite loci. Three monomorphic loci deemed as uninformative were excluded from the study. The association between marker genotypes and the studied traits was analysed using general linear model (GLM) procedure and least squares means (LSM). Linkage disequilibrium (LD) was estimated to assess non-random association between microsatellite loci. Out of 19 microsatellites studied four markers showed no association with the studied traits, three markers had a significant effect on one trait, and another three markers had significant effect on two traits. Among ten microsatellites with significant effect on four economically important traits (BW, BL, BC, TL) four were associated with two characters: marker FH2613 with BW and BC, marker FH2097withBL and BC, marker ZUBECA6 with BW and BC, whereas marker REN75M10 was associated with BL and TL. The strongest LD (r2 ranged from 0.15 to 0.33) was estimated between nine loci with significant effect on economically important traits (BW, BL, BC, TL).


Introduction
The red fox (Vulpes vulpes) had a common ancestor with the domestic dog (Canis lupus familiaris) about 12-15 million years ago (Bardeleben et al. 2005). Both species belong to the dog family Canidae, a group that has a long history spanning the last 50 million years (Ostrander and Wayne 2005). However, in contrast to Canis lupus familiaris, the red fox kept in captivity for more than 100 years has not been fully domesticated. Evolutionary processes of the red fox kept on farms have been influenced by artificial selection. This process, controlled by humans, is focused on retaining traits regarded as desirable (Kukekova et al. 2004). As a consequence of breeding programs carried out on fox farms not only morphological, physiological and behavioural changes have taken place, but also genetic structure of the selected fox populations has been shaped according to breeders needs.
In order to understand the genetic background of the red fox selective breeding, a genome map of the fox is needed. Evolutionary proximity to the dog enables to make good use of progress in studying the canine genome to the molecular genetics of the red fox (Sidjanin et al. 2002;Zhang et al. 2002). Development of a fox map is facilitated by the known cytogenetic homologies between the dog and fox, and by the availability of high resolution canine genome maps and sequence data. Alignment of the fox meiotic map against the canine genome sequence revealed high conservation of marker order between homologous regions of the two species (Kukekova et al. 2007). However, comparative analyses performed by mapping of DNA microsatellite markers or genes in chromosomes of the red fox, arctic fox and Chinese raccoon dog revealed that the number of the FISH mapped loci is still low (Switonski et al. 2009).
It is known that a majority of the primer sequences of canine microsatellites amplify analogous sequences in the red fox (Zajac et al. 2000). Consequently, molecular tools developed from the dog genome sequencing project (Lindblad-Toh et al. 2005) are likely to be applicable to the fox-like canids. Kukekova et al. (2004) reported, that over 60 % of tested canine-derived microsatellite markers robustly amplified fox DNA were polymorphic in foxes and were thus applicable for genotyping fox pedigrees. Because of the close phylogenetic relationship between the dog and the red fox, genomic resources developed previously in the dog have proven useful in the construction of the fox meiotic linkage map (Spady and Ostrander 2007). The first meiotic linkage map of the silver fox (which is a colour variant of the red fox) has been published by Kukekova et al. (2007).
The application of the latest advances in identifying microsatellite markers in the fox genome to red fox breeding is still in its infancy. Marker assisted selection (MAS), which combines traditional genetics and molecular biology would become a valuable tool in selecting foxes for traits of interest, such as body size, fur quality, fertility, tame behaviour or disease resistance. The information from the DNA testing, combined with the observed performance records for foxes could improve the accuracy of selection and increase the possibility of identifying individuals carrying desirable traits at an earlier stage of breeding.
The genome of the dog has the most developed set of tools among the Canidae, including an assembled genome sequence (Lindblad-Toh et al. 2005), which is lacking for other members of this family. Close evolutionary proximity between the dog and the red fox and high conservation of marker order between homologous regions of the two species (Kukekova et al. 2011a) offers a valuable resource for association studies in the red fox.
The aim of this study was to estimate the association of 30 canine-derived microsatellites with body weight and selected morphometrics in the farm red fox.

Material
A total of 200 unrelated silver foxes (115 females and 85 males) kept on two farms (the Leszno fox farm-120 foxes, the Prochy fox farm-80 foxes) located in a region of the Leszno city (West Poland) were studied. To ensure unrelatedness of foxes the pedigree of each animal was carefully checked. Only unrelated animals originating from different litters (one fox from a litter produced by unrelated parents) were sampled for the study. The analysed traits were: body weight (BW), body length (BL, measured from the tip of nose to the base of tail), body circumference (BC, measured behind the front limbs), tail length (TL), ear height (EH), length of the right front limb (FRLL), length of the right rear limb (RRLL), length of the right front foot (FRFL) and length of the right rear foot (RRFL). The traits measurements were taken at the end of the farm season (December) after technological slaughter of animals.

Sampling and molecular analysis
Tissue samples for DNA extraction were taken from tongues of the studied animals. Thin slices (approximately 1-2 mm) of tongues were taken post mortem using a razor blade. The DNA was isolated in accordance to ARK Genomic selfprotocol.
Thirty microsatellite markers described and localized in the dog genome were chosen for the study (Table 1). The multiplex PCR method (based on Qiagen Multiplex PCR Kit) was used to amplify these sequences. A whole set of microsatellites was divided into four pools for multiplex analysis (Table 2). Only three markers were amplified separately. Thirty multiplex PCR cycles at annealing temperature of 55°C were run. After multiplex PCR, markers were genotyped using automated sequencer 3730xl (Applied Biosystems), and then obtained data were analysed using GeneMapper v. 4.0 (Applied Biosystems).

Statistical analyses
The association between microsatellite genotypes and the studied traits was analysed using general linear model (GLM) procedure and the following linear model: where: y ijkl is the studied phenotype μ is the overall mean F i is the ith farm effect S j is the jth sex effect G k is the kth microsatellite marker genotype effect e ijkl is the random error associated with y ijkl th phenotype, N(0, σ 2 ). The least squares means (LSM) for the studied traits in relation to microsatellite genotype were computed using 'lsmeans' procedure implemented in the R package. LSM were estimated only for genotypes with the frequency higher than five. Multiple comparisons of LSM were performed with the Bonferroni-Dunn test. The studied microsatellite loci were tested separately (one at a time) and the significant differences between microsatellite genotypes for each locus were indicated.
To assess non-random association between microsatellite loci linkage disequilibrium (LD) was estimated. The squared coefficient of correlation (r 2 ) computed using MIDAS (Gaunt et al. 2006) was used as a measure of LD between loci.  (Guyon et al. 2003). The locations were calculated with the TSP/CONCORDE program. The program calculates the intermarker distances in arbitrary TSP units, and the TSP units are translated to Mb based on the known physical size of each chromosome, as determined by cytofluorimetry

Results and discussion
Microsatellites polymorphism and heterozygosity Table 2 lists the loci producing a PCR product. Out of 30 autosomal microsatellites 22 were selected for further analysis. Four microsatellites (REN88H03, FH2613, ZUBECA6 and FH3824) were observed in all studied silver foxes, while markers FH3970 and FH3771 were shared by the lowest number of foxes (163 and 169, respectively). A total of 112 alleles and 243 genotypes were found at 22 autosomal microsatellite loci. The number of alleles found in the studied loci ranged from one (loci UOR4101, REN248F14 and REN252E18) to 12 (locus ZUBECA6) with an average of 5.09 per locus. The three monomorphic loci (all of them were found in pool 3) deemed as uninformative were excluded from further study. The highest heterozygosity (Table 2) was found at locus FH3824 (0.758) followed by loci REN75M10, REN75M10 and REN258F18 (0.753, 0.738 and 0.705, respectively). The lowest heterozygosity was found at loci REN64E19 and FH2613 (0.255 and 0.245, respectively). An average heterozygosity per locus was estimated at 0.485. Monomorphic loci as mentioned above were not considered for further study.
The microsatellite polymorphism in foxes has not often been studied. Knowing that a majority of primers used for PCR amplification of the canine microsatellites can be successfully applied in the fox genome studies (Zając et al. 2000), researchers started to search for the linkage of the caninederived microsatellites in the red fox and arctic fox genomes (Klukowska et al. 2002) or to construct a genetic map of the silver fox (Kukekova et al. 2004(Kukekova et al. , 2007Spady and Ostrander 2007). Klukowska et al. (2002) who studied 19 canine-derived microsatellites in 14 red fox families (65 offspring) and 17 arctic fox families (113 offspring) kept on a commercial farm reported high microsatellites polymorphism in the investigated animals. The number of alleles varied between two and 14 with the mean value of 6.4 alleles in the red fox, and between one and 14 with the mean value of 8.2 alleles in the arctic fox. Kukekova et al. (2004) tested 700 canine microsatellites searching for a set of markers useful to study the fox genome. In the set of selected 30 microsatellites the calculated polymorphism information content (PIC) ranged from 0.06 to 0.77. The number of alleles per locus varied from two to 12, with a mean allele number of 5.1. The results of the cited studies are comparable to those presented in Table 2.
Association of microsatellite genotypes with the studied traits Least squares means and their standard errors of the silver fox morphometrics and body weight for different genotypes of ten microsatellites with significant effect on four economically important traits (BW, BL, BC, TL) are presented in Table 3. The number of foxes with a given genotype ranged from five (marker FH2097, genotype 275/283) to 150 (marker FH2623, genotype 180/180). Out of 19 microsatellites studied four markers (AHT137, REN258F18, FH3824 and FH3287) showed no association with the studied traits, three markers had a significant effect on one trait (REN135K06 with BL, FH3713 with FRFL and FH3775 with FRLL), and another three markers had significant effect on two traits (REN88H03 with BW and FRLL, FH3241 with FRLL and FRFL, REN210I14 with FRLL and RRLL). A significant effect on the highest number of traits (five out of nine studied: BL, BC, FRLL, RRLL and FRFL) was estimated for marker FH2097, whereas markers REN75M10 and FH3970 were significantly associated with four traits (BL, TL, FRLL, RRLL and TL, FRLL, RRLL, RRFL, respectively).
To our knowledge association between the polymorphism of canine-derived microsatellites and traits under selection pressure in red fox farming has not yet been studied. In the fox breeding the most important traits are those of economic importance. These are usually litter size, pelt size and pelt quality characters. According to Peura et al. (2004) who carried out the study in the Finnish blue fox population the highest relative economic weight was estimated for pelt size (0.37), followed by litter size (0.31) and pelt quality (0.26). The relative economic weights derived for the Polish blue fox population, depending on the interest rate per year assumed, ranged from 0.46 to 0.48 for litter size, 0.15 for body size, 0.35 to 0.36 for fur quality and 0.02 to 0.03 for colour type (Wierzbicki et al. 2007). In the present study four out of nine traits analysed, namely BW, BL, BC and TL, can be considered as economically important. These traits are genetically correlated with pelt or body size, the characters which significantly affect pelt prices in the international trading system (Wierzbicki and Jagusiak 2006).
The analysis of association between the studied microsatellites and four economically important traits (Table 3) indicates that significant effect on BW have five markers (FH2613, FH3771, REN64E19, REN88H03 and ZUBECA6), while the remaining three traits are significantly affected by three markers: FH2097, REN135K06 and REN75M10 are associated with BL; FH2097, FH2613 and ZUBECA6 are associated with BC; FH3970, REN30J23 and REN75M10 are associated with TL. Out of the ten microsatellites with significant effect on the four economically important traits, four were associated with two characters: marker FH2613 was associated with BW and BC, marker FH2097 was associated with BL and BC, marker ZUBECA6 was associated with BW and BC, whereas marker REN75M10 was associated with BL and TL. Thus, in further studies on implementation of MAS in red fox breeding efforts should be focused on the set of ten microsatellite loci significantly associated with traits of breeders' interest. In recent years most of the association studies carried out in red foxes have mainly focused on searching for genetic markers related to body weight. Skorczyk et al. (2010), who conducted the association study in the red fox, showed that the melanocortin 3 receptor (MC3R) is a promising candidate gene for body weight. The authors reported that two polymorphisms in the red fox, i.e. a silent substitution c.957A>C and c.*185C>T in the 3'-flanking sequence, showed a significant association (P<0.01) with body weight.
The polymorphism of another gene encoding melanocortin receptors, melanocortin 4 receptor gene (MC4R), in relation to body weight of the red fox was studied by Nowacka-Woszuk et al. (2011). In the analyzed sequence the authors found three already known SNPs, two novelindels and three novel SNPs. The 11 bp indel and four SNPs segregated as two haplotypes. The comparative in silico search for functional sequences in the studied fragment revealed two uORFs. However, in contrast to Skorczyk et al. (2010) the authors found no significant association between polymorphic variants within a putative upstream open reading frame of the MC4R gene and body weight of farmed red foxes. The polymorphism of the MC4R gene and its potential use in animal breeding was also investigated in the Chinese raccoon dog (Nycterute sprocyonoides procyonoides) and the arctic fox (Vulpes lagopus) (Skorczyk et al. 2012).
The polymorphisms of two candidate genes, the insulin induced gene 2 (INSIG2) and the fat mass and obesity associated gene (FTO), in four species belonging to the family Canidae (the dog, red fox, arctic fox and Chinese raccoon dog) were studied by Grzes et al. (2011). Two synonymous SNPs, one in the FTO gene (−28T>C in the 5′-flanking region) and one in the INSIG2 (10175C>T in intron 2), were used for the association studies in the red fox. The evidence was observed for their association with body weight (FTO, p<0.08) and weight of raw skin (INSIG2, p<0.05). These associations indicate that both genes are potential candidates for growth or adipose tissue accumulation in the red fox. Figure 1 shows LD between the studied marker loci. It is interesting that the highest values of r 2 were estimated between nine loci with significant effect on economically important traits (BW, BL, BC and TL). The strongest LD was found between ZUBECA6 and FH3771 (r 2 =0.33), followed by ZUBECA6 and FH3970 (r 2 =0.25), FH3771 and FH3790 (r 2 =0.25), FH3824 and REN75M10 (r 2 =0.20). The FH2097 microsatellite was strongly linked with three loci: REN210I14 (r 2 =0.16), REN258F18 (r 2 =0.16) and REN64E19 (r 2 =0.15). The pattern of LD found between the studied microsatellites may indicate that many decades of artificial selection towards genetic improvement of economically important traits in the red fox, favoured those marker loci which were tightly associated with preferred phenotypes. This strong footprint of selection in regions of the genome that are important for controlling important traits, opens up opportunities to successfully utilize LD-markers in MAS, which may markedly improve genetic progress in red fox farming. A similar conclusion was drawn by Boyko et al. (2010) who generated a high density map of canine genetic variation by genotyping 915 dogs from 80 domestic dog breeds, 83 wild canids, and ten outbred African shelter dogs across 60,968 single-nucleotide polymorphisms (SNPs). The authors identify 51 regions of the dog genome associated with phenotypic variation among breeds in 57 traits, and concluded that artificial selection has played an important role in altering the genetic architecture of key traits in the studied canids. Identification of 44 chromosomal regions (associated, among others, with body size and ear morphology) in 46 breeds of the dog, bearing signatures of strong human-directed selection was also reported by Vaysse et al. (2011).
A closer look at the canine chromosomal regions where the studied microsatellite markers are located reveals a number of genes positioned in close distance to the markers (Guyon et al. 2003). In Table 4 the studied markers and possible candidate genes located nearby (±5 Mb) are presented. In the close proximity of eight markers located on eight canine chromosomes, 11 possible candidate genes are positioned. The selected genes are directly or indirectly involved in processes underlying behaviour (e.g. HTR1B), body weight (e.g. GHR), formation of hair structure and fur coat (e.g. KRT17).
If the marker loci were located in chromosomal regions near the selected genes they might also create a similar genetic architecture and plausibly play a similar role in the red fox genome. However, to elucidate this the population-based genome screens together with candidate genes studies are needed.
In a group of the nine studied traits, five (EH, FRLL, RRLL, FRFL and RRFL) can be considered less important in selective breeding of foxes. These traits do not directly and significantly influence efficiency of fur production, thus the genetic improvement of these traits is not included in the breeding goal in fur farming ). However, for evolutionary geneticists the morphological changes associated with domestication are of great interest. Domesticated animals or those being domesticated can be easily distinguished from their wild relatives by skull shape and other skeletal and morphological features (Wayne 2001). Preliminary results confirming significant differences in morfometrics between farm silver foxes and wild red foxes were reported by Zatoń-Dobrowolska et al. (2012). Development of a fox genetic map and the mapping morphological phenotypes segregating in selected red fox populations may provide insight into the mechanisms underlying the red fox evolution and domestication (Kukekova et al. 2004(Kukekova et al. , 2011b.

Conclusions
The present study confirms that canine-derived microsatellites can be successfully applied in red fox breeding. MAS using a set of selected LD-markers with significant effects on economically important traits may increase genetic progress, and thus profitability of fur farming. Further studies, including the candidate gene approach are needed to fully exploit possibilities of MAS in shaping the genetic structure of the farm red fox populations. Table 4 Description of chromosomal regions with studied markers and possible candidate genes located nearby (±5 Mb)

Chromosome
Marker and its location Gene and its location Gene function CFA02 FH2613 (49.9 Mb) FGF1 (46.5 Mb) Fibroblast growth factor 1 (acidic). The protein encoded by this gene is a member of the fibroblast growth factor (FGF) family. FGF family members possess broad mitogenic and cell survival activities, and are involved in a variety of biological processes, including embryonic development, cell growth, morphogenesis, tissue repair NR3C1 (47 Mb) Nuclear receptor subfamily 3, group C, member 1, which can function both as a transcription factor that binds to glucocorticoid response elements in the promoters of glucocorticoid responsive genes to activate their transcription, and as a regulator of other transcription factors. Can act as a coactivator for STAT5-dependent transcription upon growth hormone stimulation and could reveal an essential role of hepatic GR in the control of body growth CFA04 FH2097 ( Collagen, type I, alpha 1 encodes the pro-alpha1 chains of type I collagen whose triple helix comprises two alpha1 chains and one alpha2 chain. Type I is a fibril-forming collagen found in most connective tissues and is abundant in bone, cornea, dermis and tendon CFA12 FH3713 (52.3 Mb) HTR1B (57.1 Mb) 5-Hydroxytryptamine (serotonin) receptor 1B, G protein-coupled 1, previously known as 5-HT1D receptors, are primarily located in the basal ganglia, striatum, hippocampus and vascular smooth muscle. 5-HT1B receptors play a role in appetite control, sexual behaviour, aggression and anxiety CFA15 REN307J23 (40.6 Mb) MYF5 (41.2 Mb) Miogenic factor 5, involved in muscle differentiation (myogenic factor). Induces fibroblasts to differentiate into myoblasts CFA20 FH3771 (52 Mb) DAG1 (47.8 Mb) Dystroglycan 1 (dystrophin-associated glycoprotein 1). Dystroglycan is a laminin binding component of the dystrophin-glycoprotein complex which provides a linkage between the subsarcolemmal cytoskeleton and the extracellular matrix. Dystroglycan 1 is a candidate gene for the site of the mutation in autosomal recessive muscular dystrophies CFA34 REN64E19 (25.4 Mb) SST (28.9 Mb) Somatostatin hormone. This hormone is an important regulator of the endocrine system through its interactions with pituitary growth hormone, thyroid stimulating hormone and most hormones of the gastrointestinal tract Source: Guyon et al. (2003)