Samples and DNA extraction

Samples of bloods, muscles, feces, or hair were obtained from males of twelve species in the family Canidae: Alopex lagopus (arctic fox), Fennecus zerda (fennec fox), Vulpes vulpes (red fox), Canis aureus (golden jackal), Canis mesomelas (black-backed jackal), Canis latrans (coyote), Canis lupus (gray wolf), Canis familiaris (domestic dog), Cuon alpinus (dhole), Speothos venaticus (bush dog), Chrysocyon brachyurus (maned wolf), and Nyctereutes procyonoides (raccoon dog) (Table 1).

Whole blood samples were stored at -20^°C; separated white blood cells were stored in 99% ethanol at 4^°C until use. Muscle and fecal samples were preserved in 99% ethanol at 4^°C. Hair samples with roots were preserved at 4^°C. Total DNA was extracted from whole blood (100 µl each) and pellets of white blood cells and muscle samples (about 3 × 3 × 3 mm) using the DNeasy Blood and Tissue Kit (Qiagen), from fecal samples using the QIAamp DNA Stool Mini Kit (Qiagen), and from hair samples (20–30 hairs) with the QIAamp DNA Micro Kit (Qiagen). The DNA extracts in 100–200 µl of TE buffer were preserved at 4^°C until use.

PCR amplification of the ZFX and ZFY genes

The final introns of the ZFX and ZFY genes were amplified using polymerase chain reaction (PCR) primers (U-ZF-2F/U-ZF-2R) reported by Nakagome et al. (2008) (Fig. 1). PCR amplifications were performed in 20 µl of a reaction mixture containing 1 µl of DNA extract, 0.1 µl of rTaq polymerase (5 units/µl, Takara), 2 µl of 10× buffer (Takara), 1.6 µl of dNTPs (2.5 mM, Takara), and 0.2 µl of each (25 pmol/µl) of two primers. For fecal samples, 0.4 µl of bovin serum albumin (20 mg/ml, Boehringer) was added to the above reaction mixture to eliminate some effects of PCR inhibitors.

Touchdown PCR cycling conditions were as follows: 94^°C for 10 min; 10 cycles of 94^°C for 1 min, 54^°C for 30 sec, decreasing 1^°C at every cycle, and 72^°C for 2 min; 25–35 cycles of 94^°C for 1 min, 44^°C for 30 sec, and 72^°C for 2 min; and 72^°C for 10 min.

Cloning of PCR products and nucleotide sequencing

The PCR products were cloned by using the TA Cloning Kit (Invitrogen), following the manufacturer's instructions. Positive clones containing insert DNAs of the expected molecular sizes were incubated in 2 ml of LB broth (Invitrogen) with shaking at 37^°C overnight. Plasmid DNA was then isolated with the QIAprep Miniprep Kit (Qiagen) according to manufacturer's instructions, dissolved in 50 µl of TE buffer, and stored at 4^°C.

Fig. 1.

Positions and directions of PCR and sequencing primers for the ZFX final intron (A) and the ZFY final intron (B) of Canidae.

Ten to fourteen clones were sequenced using the BigDye Terminator Cycle Sequencing Kit (Applied Biosystems), and applied to an automated sequencer (ABI3730, Applied Biosystems). Based on an alignment of obtained sequences, PCR primers specific to the ZFX final intron (C-ZFX-F: 5′-GTTCCCTTAAGGCAGGCATA-3′ and C-ZFX-R: 5′-AAAGCTTAAATCCACCTATGGAAA-3′) and the ZFY final intron (C-ZFY-F: 5′-CAAGTTAGCATAAATTTGGTTTG-3′ and C-ZFY-R: 5′-TGTCTCTGCCTCTCTGTGTCTC-3′) of Canidae were newly designed in the present study (Fig. 1).

Direct sequencing

PCR amplifications using newly designed primers (Fig. 1) for Ca. latrans and Ch. brachyurus were performed in 50 µl of a reaction mixture containing 10 µl of the DNA extract, 0.1 µl of rTaq polymerase (5 units/µl, Takara), 5 µl of 10× buffer (Takara), 4 µl of dNTPs (2.5 mM, Takara), 0.5 µl of each primer (25 pmol/µl), and 1 µl of bovin serum albumin (20 mg/ml, Boehringer). The touchdown PCR conditions for fragments X1 and X2 (Fig. 1A) were as follows: 94^°C for 10 min; 35 cycles of 94^°C for 2 min, 56^°C for 30 sec, and 72^°C for 2 min; and 72^°C for 10 min. The touchdown PCR conditions for fragment Y1 (Fig. 1B) was follows: 94^°C for 10 min; 10 cycles of 94^°C for 1 min, 65^°C for 30 sec, decreasing 0.5^°C every cycle, and 72^°C for 2 min; 25–35 cycles of 94^°C for 1 min, 60^°C for 30 sec and 72^°C for 2 min; and 72^°C for 10 min. The touchdown PCR conditions for fragment Y2 (Fig. 1B) were as follows: 94^°C for 10 min; 10 cycles of 94^°C for 1 min, 60^°C for 30 sec, decreasing 1^°C every cycle, and 72^°C for 2 min; 25–35 cycles of 94^°C for 1 min, 50^°C for 30 sec and 72^°C for 2 min; and 72^°C for 10 min. The PCR products were purified with the QIAquick PCR Purification Kit (Qiagen). Sequencing was carried out using newly designed primers, as mentioned above.

PCR amplifications using U-ZF-2F/U-ZF-2R primers (Nakagome et al., 2008) for A. lagopus were performed in 50 µl of a reaction mixture containing 2 µl of the DNA extract, 0.1 µl of rTaq polymerase (5 units/µl, Takara), 5 µl of 10× buffer (Takara), 4 µl of dNTPs (2.5 mM, Takara), and 0.5 µl of each primer. The touchdown PCR cycling conditions were as follows: 94^°C for 10 min; 10 cycles of 94^°C for 1 min, 54^°C for 30 sec, decreasing 1^°C at every cycle, and 72^°C for 2 min; 25–35 cycles of 94^°C for 1 min, 44^°C for 30 sec, and 72^°C for 2 min; and 72^°C for 10 min. The PCR products were purified with the QIAquick PCR Purification Kit (Qiagen). Sequencing was carried out using newly designed primers, as mentioned above.

Data analysis

The ZFY (accession no. AB261807) and ZFX (AB261815) final intron sequences of Ursus arctos (brown bear) as an outgroup were cited from Nakagome et al. (2008). In addition, the ZFY (AB491592) and ZFX (AB491601) final intron sequences of Mustela itatsi (Japanese weasel) were cited from Yamada and Masuda (2010).

The sequence alignment and construction of neighbor-joining (NJ) trees using Kimura's (1980) two-parameter distance model were performed with MEGA 4 (Tamura et al., 2007). Phylogenetic trees based on maximum parsimony (MP) and maximum likelihood (ML) methods were constructed with PAUP* version 4.0b10 (Swofford, 2002). Insertion or deletion (indel) sites were eliminated from calculation. Specific search conditions for the MP analysis included starting trees obtained by stepwise addition of sequences with 1000 replicates, general heuristic search and the treebisection-reconnection (TBR) branch swapping algorithm. All other settings were set by default. For the ML analysis, the program Modeltest 3.06 (Posada and Crandall, 1998) was used to select the most appropriate model of molecular evolution through a hierarchical likelihood ratio test. The model selected for sequence data for the ZFY genes was Hasegawa-Kishino-Yano (HKY) model (Hasegawa et al., 1985) following a gamma distribution shape parameter (HKY-G model), and that for the ZFX gene data was Kimura unequal base frequencies (K81uf) model (Kimura, 1981). Bootstrap values (Felsenstein, 1985) were calculated from 1000 replicates for NJ with MEGA 4, and 100 for ML and 1000 for MP with PAUP* version 4.0b10. In addition, a haplotype network was constructed by using software TCS (Clement et al., 2000), where a gap was counted as one nucleotide substitution.

Pairwise differences were calculated for sequences of ZFY and ZFX final introns, using Kimura's (1980) two parameter distance model in MEGA4. The ZFY/ZFX ratios were calculated for each ZFY matrix element with the corresponding element from ZFX, using Kimura's (1980) two parameter distance matrices. Mean ZFY/ZFX values were computed from all pairwise estimates (n = 66).

Phylogenetic relationships among the ZFY final intron sequences of Canidae

The ZFY final intron sequences (924–1146 base-pairs, bp) were obtained from all twelve canid species examined (Appendix). In addition, a sequence alignment showed that SINEs and other insertional sequences were found in the ZFY final intron, as mentioned below.

Of three equally parsimonious trees yielded by MP analysis, one tree, of which topology was the most similar to those generated by ML and NJ methods, was selected for indicating bootstrap values: this MP tree was almost identical to the ML tree, where F. zerda, A. logopus and V. vulpes were grouped into one cluster. The ML tree is shown in Fig. 2 as a representative of three methods with bootstrap values for ML, MP and NJ above each branch.

Three fox species, V. vulpes, A. lagopus, and F. zerda, formed a monophyletic group, supported by 64% (ML), 63% (MP) and 63% (NJ) bootstrap values (Fig. 2). They were included in “the red fox-like canids clade” by previous molecular phylogenetic studies (Wayne et al., 1997; Bardeleben et al., 2005). Five species of Canis, Ca. mesomelas, Ca. aureus, Ca. latrans, Ca. lupus, and Ca. familiaris, and Cu. alpinus were clustered with 73/63/64 bootstrap values (Fig. 2). This monophyletic group was congruent with “the wolf-like canids clade” defined in previous molecular phylogenetic studies (Wayne et al., 1997; Bardeleben et al., 2005). Within the wolf-like canids clade, the sequence of Ca. familiaris (domestic dog) was identical with that of Ca. lupus (gray wolf). Chrysocyon brachyurus and S. venaticus, both of which were included in “the South American canids clade” (Wayne et al., 1997; Bardeleben et al., 2005), were clustered with “the wolf-like canids clade”, supported by 92/97/99 bootstrap values (Fig. 2).

Fig. 2.

A maximum-likelihood (ML) tree as a representative phylogenetic tree of the Canidae ZFY sequences (924–1146 bp). The scale indicates genetic distances. Numbers near the internal branches are bootstrap values for ML (100 replications), MP (1000 replications) and NJ (1000 replications) in order. Because when bootstrap values are less than 50% for all three methods, the branching is bifurcated. An arrow indicates event of an insertion of SINE-I before species radiation of Canidae. An asterisk shows an independent insertion event of SINE-II into the ZFY final intron of Ca. aureus.

Fig. 3.

A maximum-likelihood (ML) tree as a representative phylogenetic tree of the Canidae ZFX sequences (834–839 bp). The scale indicates genetic distances. Numbers near the internal branches are bootstrap values for ML (100 replications), MP (1000 replications) and NJ (1000 replications) in order, where less than 50% were shown by dashes.

The position of N. procyonoides was different between the NJ tree and the other trees: Nyctereutes procyonoides was split from the other species in the NJ tree, whereas it was a sister taxon to the wolf-like and South American canids clades in the ML and MP trees.

In addition to two types of SINEs (mentioned below), five types of shorter insertional fragments were found in the ZFY final intron of Canidae (Appendix): 7-bp frangment (nucleotide sites, ns 91–97) in N. procyonoides and the three fox species; 9-bp fragment (ns 200–208) in Ca. aureus; 54-bp fragment (ns 518–572) in Cu. alpinus; 8-bp fragment (ns 644–651) in S. venaticus; and 4-bp fragment (ns 1182–1185) in N. procyonoides. Although homology searches using the DNA databases were made for the above 54-bp fragment, no resultant homologues were found. Features of the other small fragments were not clarified by the homology search.

Phylogenetic relationships among the ZFX final intron sequences of Canidae

The ZFX final intron sequences (834–839 bp) were obtained from all the twelve Canidae species (Table 1). Of the two equally parsimonious trees that resulted from MP analysis, one tree, of which topology was more similar to those of ML and NJ trees, was selected for calculating boot-strap values. Because ML, MP and NJ trees indicated the similar topologies to each other, the ML tree is shown in Fig. 3 as a representative of the three methods with bootstrap values for ML, MP and NJ, indicated in this order above each branch. The three foxes, A. lagopus, F. zerda, and V. vulpes (red fox-like canids clede), were clustered with 83/68/73 bootstrap values (Fig. 3), in agreement with the ZFY tree (Fig. 2). Four species of Canis, Ca. aureus, Ca. familiaris, Ca. latrans, and Ca. lupus, formed a monophyletic group with 61/58/55 bootstrap values (Fig. 3), in congruence with the ZFY tree (Fig. 2). Among those, the sequence of Ca. familiaris was identical with that of Ca. lupus. Canis mesomelas and Cu. alpinus clustered with Ch. brachyurus and S. venaticus with 60/60/-bootstrap values (Fig. 3). The difference between the NL tree and the other trees (ML and MP) was the position of Cu. alpinus; the ML and MP trees showed that Cu. alpinus clustered with the group consisting of Ca. mesomelas, Ch. rachyurus, and S. venaticus, whereas the NJ tree indicated that Cu. alpinus was a sister taxon of the latter group. Nyctereutes procyonoides was remote from the other species (Fig. 3).

Table 2.

Pairwise differences (Kimura's two parameter distances) of the ZFY (upper) and ZFX (below) final intron sequences among the twelve species of Canidae.

Fig. 4.

(A) Insertional positions of SINE sequences in the ZFY final intron of Canidae. Solid lines show homologous sequences among the species, and broken lines indicate the positional relationships of sequences among the species. (B) An alignment of SINE-I sequences identified from the homologous sites of the twelve Canidae species. The SINE sequence of the domestic dog (Canis familiaris, AF221611) (Pecon-Slattery et al., 2000a) is shown as a reference. (C) The SINE-II sequence (210 bp) found in Ca. aureus. RNA polymerase III boxes A and B and terminator, and (CT)n region are enclosed by solids lines. Direct repeat motifs are enclosed by broken lines. The CAN-SINE specific insertion sequence is underlined.

In contrast with the ZFY sequences, no insertional fragments were found in the ZFX final intron sequences.

Estimation of ZFY/ZFX mutation ratio

In order to evaluate the evolutionary difference between ZFY and ZFX genes, matrices composed of Kimura's (1980) two parameter genetic distances among all pairs of taxa were made separately for the two genes (Table 2). The ZFY/ZFX values were then calculated at each pair of species. The mean pairwise ZFY/ZFX value was 2.00, showing a more rapid mutation rate in the ZFY gene than the ZFX gene.

Identification of SINE sequences within the ZFY final introns of Canidae

Sequence homology searches revealed that 165–195-bp fragments (ns 810–1004) (Appendix) within the ZFY final introns of all the twelve canid species were most homologous (from 96% [117/121 bp] to 100% [176/176 bp]) to a SINE sequence (accession no. AF221611) of Ca. familiaris, which is reported to possess structural features, such as the 5′ region containing the split block promoter for RNA polymerase III common to tRNA genes, followed by a unique sequence and a CT repeat, and terminated by a 3′ poly A tail (Pecon-Slattery et al., 2000a). Therefore, we refer to the 166–195-bp fragments identified in the present study as SINE-I (Fig. 4A, B). The length difference among SINE-I sequences resulted from 30-bp insertions in the three Canis species (Ca. familiaris, Ca. latrans, and Ca. lupus) (Fig. 4B).

In addition, another insertional fragment (ns 243–453) was found in only the ZFY final intron of Ca. aureus (Fig. 4A, C; Appendix). Because this sequence also possessed the general SINE features as well as the Pol III terminator (Fig. 4C), we name this SINE-II. By contrast, no SINEs were found in the ZFX final introns of any canid species examined in the present study.

Fig. 5.

A parsimonious network of SINE-I sequences. An indel was counted as one nucleotide substitution. Open circles show presumed sequence types. One bar between circles indicates one nucleotide substitution.

Of the determined SINE-I sequences, 41 nucleotide sites were polymorphic, including 30-bp indels (Fig. 4B). Figure 5 shows a parsimonious network of SINE-I sequences. Three species of Canis (Ca. familiaris, Ca. latrans, and Ca. lupus) having the 30-bp insertional sequence were separated from the other nine species. The SINE-I sequences of Ca. mesomelas, Ch. brachyurus, and S. venaticus were identical to each other. Two foxes (A. lagopus and V. vulpes) shared identical SINE-I sequences and were closely related to the other fox, F. zerda, forming the red fox-like canids clade (Fig. 5). The phylogenetic position of N. procyonoides was remote from the other species (Fig. 5), in agreement with the ZFY (Fig. 2) and ZFX (Fig. 3) trees.

Phylogenetic relationships among Canidae species based on sex chromosomal DNA sequences

The present study investigated the phylogenetic relationships among canid species inferred from the sex chromosomal DNA sequences. The molecular phylogenetic trees obtained from Y chromosomal DNA (ZFY sequences) and X chromosomal DNA (ZFX sequences) commonly showed that the canid species examined were divided into two main clades: the red fox-like canids and wolf-like canids clades, consistent with the previous molecular phylogenetic studies (Wayne et al., 1997; Bardeleben et al., 2005). Both Y chromosomal and X chromosomal DNA data always supported the A. lagopus/F. zerda/V. vulpes cluster (red fox like-canids clade), in agreement with the previous studies based on mtDNA (Wayne et al., 1997) and autosomal DNA (Bardeleben et al., 2005; Wayne and Ostrander, 2007). In addition, the phylogenies of both sex chromosomal DNA data showed that four Canis species (Ca. aureus, Ca. familiaris, Ca. latrans, and Ca. lupus) clustered together, supporting the wolf-like canids clade previously reported on the basis of mtDNA (Wayne et al., 1997) and autosomal DNA (Bardeleben et al., 2005; Wayne and Ostrander, 2007).

On the other hand, the phylogenetic positions of Cu. alpinus, Ca. mesomelas, Ch. brachyurus, S. venaticus, and N. procyonoides were different between the ZFY and ZFX trees. Cuon alpinus and Ca. mesomelas clustered with other Canis species of the wolf-like canis clades in the ZFY tree (Fig. 2), whereas the former two species were grouped with Ch. brachyurus and S. venaticus in the ZFX tree (Fig. 3). Molecular phylogenies of mtDNA (Wayne et al., 1997) and autosomal DNA (Bardeleben et al., 2005; Wayne and Ostrander, 2007) commonly showed that the two species (Cu. alpinus and Ca. mesomelas) were clustered with the wolf-like canids clade. In addition, Zrzavý and Řičánková (2004) using the morphological and molecular data sets also revealed that both species were grouped into the wolf-like canids clade. This is not incongruent with the phylogeny shown in Fig. 2.

Diploid chromosome numbers of Cu. alpinus (n = 78) and Ca. mesomelas (n = 78) are identical with those of the other Canis species (n = 78) forming the wolf-like canids clade, but not with Ch. brachyurus (n = 76) and S. venaticus (n = 74) (Wayne, 1993). The present study indicates that the two species of South America cluster with the wolf-like canids clade, showing higher bootstrap values (more than 92%, see Fig. 2). Wayne and Ostrander (2007) reviewed the molecular phylogeny of Canidae using a large data set, and showed the earlier split of the fox-like canids clade from the others and the subsequent separation from the others into the South American canids clade and the wolf-like canids clade. The present study strongly supports the phylogenetic affinity between the latter two clades.

The phylogenetic position of N. procyonoides was different in the ZFY and ZFX trees. Molecular phylogenetic studies of mtDNA showed that N. procyonoides is not closely related to the fox-like, wolf-like, and South American canids clades (Wayne et al., 1997). On the other hand, autosomal DNA phylogenies showed that N. procyonoides was a sister taxon to the red fox-like canids clade (Bardeleben et al., 2005; Wayne and Ostrander, 2007). A morphological study indicated that N. procyonoides is included in the South American canids clade (Tedford et al., 1995). Both morphological and molecular data sets of Zrzavý and Řičánková (2004) showed that N. procyonoides was basal to the wolf-like and South American canids clades. Meanwhile, the 7-bp fragment (ns 91–97) (Appendix) was found to have been inserted into the ZFY final intron of N. procyonoides and the three fox species. Because such an insertion is thought to be cladistically informative, the common insertion of the 7-bp fragment indicates the close relatedness among N. procyonoides and the three fox species F. zerda, A. logopus, and V. vulpes, supporting the phylogenetic relationship reported by Wayne and Ostrander (2007).

The present study revealed that the mutation rate of ZFY sequences in Canidae was higher than that of ZFX, as the mean pairwise ZFY/ZFX was about 2.00. Such a greater mutation rate for the Y chromosomal genes than the X is similar to those reported in other carnivoran families such as Felidae (Pecon-Slattery et al., 1998), Ursidae (Nakagome et al., 2008) and Mustelidae (Yamada and Masuda, 2010), in concordance with male-driven evolution, reported by Haldane (1947), Miyata et al. (1987), Shimmin et al. (1993), Makova and Li (2002) and Goetting-Minesky and Makova (2006).

Notably, both the ZFY and ZFX sequences of Ca. familiaris (domestic dog) were identical to those of Ca. lupus (gray wolf). This sequence identity is not incongruent with that dogs were domesticated from gray wolfs: Vila et al. (1997) reported that the domestication time is at least 14,000 years ago. Tsuda et al. (1997) also reported that domestic dogs have maintained a large degree of mtDNA polymorphisms introduced from their ancestral wolf populations.

SINE insertions in the ZFY final introns

Carnivore-specific SINE (CAN-SINE) sequences have been identified from all carnivore families (Vassetzky and Kramerov, 2002). The SINE insertions within the ZFY final intron in Carnivora were reported in Felidae (7/34 species: Pecon-Slattery et al., 2000b, 2004) and Ursidae (8/8 species: Nakagome et al., 2008). The present study reports for the first time that the SINE-I in the ZFY final introns is shared by all the twelve canid species examined. In addition, the SINE-II found in the present study possessed characteristics specific to CAN-SINEs (Vassetzky and Kramerov, 2002): the SINE-II sequence contained two direct repeats (GAGCCTGC and GCCTGT) and insertion (TGGCCCGGGGCG) between Pol III A Box and B Box (Fig. 4B, C). From this, the SINE-II identified in the present study can be classified into a member of CAN-SINEs.

Yamada and Masuda (2010) reported an insertion of one copy of CAN-SINE in the ZFY final intron, which is specific to only Meles anakuma (Japanese badger) and Mustela erminea (ermine) among nine mustelid species investigated. We then re-examined those mustelid ZFY sequences of Yamada and Masuda (2010), and found a different insertion of another copy of SINE, which is common to all the nine mustelid species, in addition to the CAN-SINE specific to Meles anakuma and Mustela erminea. Thus, the insertional pattern of SINE-I and -II in canids identified in the present study is similar to that of the mustelid ZFY sequences of Yamada and Masuda (2010).

Canis latrans, Ca. familiaris and Ca. lupus shared the 30-bp insertion between Pol III A Box and Pol III B Box of SINE-I, compared with the other canid species. This indicates that the three Canis species had a common ancestor even after original insertion of SINE-I in canids. Similarly, among CAN-SINEs specific in the ZFY final intron of Meles spp. (Eurasian badgers), 12-bp deletions were common to Meles anakuma (Japanese badger) compared with other continental Meles species (Tashima et al., 2011). In addition to SINEs in the ZFY final intron, 4–54-bp insertions were found in single species and one lineage including multiple species in Canidae. Because the ZFY gene is located in the nonrecombining region on the Y chromosome (Pecon-Slattery et al., 2000b), it seems reasonable that more insertions are distributed in the ZFY gene than in the ZFX gene. Although the origins and functions of the small insertional sequences found in the present study are unknown yet, such insertions can be cladistically informative markers for phylogenetic analyses.

Because the SINE-I sequences occur commonly in all the twelve species examined, they could have been inserted into the ZFY final introns before species radiation in Canidae. By contract, SINE-II could have been inserted into the ZFY final intron of Ca. aureus lineage-independently after the speciation in Canidae. The studies on molecular diversity of SINEs as well as other small insertional sequences found in the ZFY gene provide a potential to further understanding the evolutionary history on diversification of Canidae.