Buffalo species identification and delineation using genetic barcoding markers.

Enrichment of barcode databases with mitochondrial cytochrome c oxidase subunit I (COI) barcode sequences in different animal taxa has become important for identification of animal source in food samples to prevent commercial fraud. In this study, COI barcode sequence in seventy one river buffalo samples were determined, analyzed and deposited in Genbank barcode database and barcode of life database (BOLD) to contribute for construction of public reference library for COI barcode sequence in river buffalo. Moreover COI barcode sequence was used to identify the closely related buffalo groups: river buffalo, swamp buffalo, lowland anoa and African buffalo. Results indicated the success of the COI barcode in the identification of each of the tested groups. Whereas a suggested sequence of other mitochondrial segment representing two successive transfer RNA (tRNA) genes; tRNA-Threonine (MT-TT) and tRNA-Proline (MT-TP) was failed to be used as a barcode marker for differentiation between the tested buffalo groups.

Introduction

Animal species identification and delineation using genetic barcoding method has become essential in food quality control procedures and for the detection and identification of animal source in food samples to prevent commercial fraud [1], [2], [3]. The nucleotides sequence of approximately 650 base pairs (bp) of the 5′ half of the mitochondrial (mt) cytochrome c oxidase subunit I (COI) gene was initially proposed as a standard barcode for rapid animal species identification [4]. The choice of this gene was impelled by its low intraspecific variability within same species and high interspecific variability between different species [4], [5], [6], [7].

Domesticated water buffalo species (Bubalus bubalis) belongs to genus Bubalus, subfamily Bovinea and family Bovidea. It is a very common farm animal species which is particularly wide spread in tropical and subtropical countries with hot and humid climate [8]. Water buffalo species includes two subspecies known as river type (Bubalus bubalis bubalis) and swamp type (Bubalus bubalis carabanesis) [9], [10], [11]. The efficacy of COI barcode was previously examined, in silico for the identification of water buffalo species in few studies [12], [13]. Results of these studies indicated the presence of high intraspecific divergences in COI barcode sequence within this buffalo species. These results based on the available limited number of sequence data in the Genbank for river and swamp buffalo at that time. Oraby and colleagues [13] in their report also suggested the use of sequences of two successive mitochondrial transfer RNA genes; tRNA-Threonine (MT-TT) and tRNA-Proline (MT-TP), as a complementary marker for COI barcode marker.

In view of the above findings, this work was planned to study the COI barcode sequence in larger population of river buffalo. A total of seventy one river buffalo samples were collected from the Egyptian river buffalo and the Egyptian-Italian river buffalo hybrids. The nucleotide sequence of 690 bp form the 5′ half of COI gene in these samples were determined, analyzed and further submitted to both databases; the Genbank barcode database (https://www.ncbi.nlm.nih.gov/WebSub/?tool=barcode) and the DNA barcode of life database (BOLD) (https://www.boldsystems.org). Sequences of the haplotypes that obtained from the analysis of the seventy one river buffalo COI sequences in the present study were used along with other COI homologous sequences, retrieved from the database, for river buffalo, swamp buffalo and two other closely related buffalo species from the subfamily Bovinea for verification of the efficacy of COI gene barcoding in the identification and delineation of these four buffalo taxa groups. The two related buffalo species are lowland anoa (Bubalus depressicornis) from genus Bubalus and African buffalo (Syncerus caffer) from genus Syncerus. The usefulness of the sequence of the two successive mitochondrial transfer RNA genes; MT-TT and MT-TP (MT-TT&MT-TP) suggested by Oraby and colleagues [13] as a complementary marker for COI barcode marker has also been investigated and its ability was compared to the efficacy of COI gene barcode in identification and delineation of the above mentioned buffalo taxa groups.

Materials and methods

Blood samples collection and DNA extraction

A total of sixty Egyptian river buffalo blood samples representing the northern and southern subpopulations were collected from different locations in Egypt. Northern subpopulation samples were twenty nine samples from Kafr El-Sheikh and twenty six samples from Ismailia while Southern subpopulation samples were five samples from Sohag. Eleven blood samples were also collected from Egyptian-Italian river buffalo hybrids from Ismailia. Genomic DNA was extracted from whole blood by salting out method [14].

Primers design

Two sets of primers were designed using NCBI/Primer-Blast program and the river buffalo whole mitochondrial DNA sequence (accession no. AF547270.1). First set of primers (COI primers set) was designed for amplification of mtDNA segment (709 bp) from the 5′ half of the mitochondrial cytochrome c oxidase subunit I (COI) gene. The second set of primers (CYTB-TT-TP primers set) was designed for amplification of mtDNA segment (228 bp) which represents partial sequence of cytochrome b (CYTB) gene from nucleotide one to nucleotide 90 (90 bp), and the two successive tRNA genes; MT-TT and MT-TP (MT-TT&MT-TP), from nucleotide 95 to nucleotide 228 (134 bp). Table 1 indicates primers sequences, annealing temperatures and sizes of PCR amplicons.

Table 1

Sequences of primers, annealing temperatures and amplicons size.

Designed primers	Forward and reverse primers sequences	Annealing temperature	PCR amplicon size
COI primers set	Forward 5′-TCTCAACCAACCATAAAGATATCGG-3′Reverse 5′-TATACTTCAGGGTGTCCGAAGAATCA-3	58 °C	709 bp
CYTB-TT-TP primers set	Forward 5′-GGACAACTAGCATCCATCAC-3′Reverse 5′-CAGGGAGTAGTTTAAATAGAA-3′	59 °C	228 bp

Polymerase chain reaction conditions and profile

The polymerase chain reaction (PCR) was performed for each DNA sample in 25 µl reaction mixture consisting of 0.2 mM dNTPs, 10 mMTris, 50 mMKCl, 1.5 mM MgCl2, 0.01% galatin (w/v), 1U Taq polymerase, 1 µM forward primer, 1 µM reverse primer and 100 ng DNA. Reaction mixture was overlaid with sterile mineral oil and was run in an MJ research PTC-100 Thermocycler. PCR amplifications were performed by means of 30–35 cycles each consisted of 1 min for denaturing at 94 °C, 2 min for annealing at the suitable temperature for each set of primers (Table 1), and 2 min for extension at 72 °C. PCR reaction was completed by a final extension at 72 °C for 7 min.

Following completion of cycling reaction, each PCR product was analyzed by electrophoreses separation of a 5 µl aliquot on a 2% (W/V) agarose gel containing ethidium bromide (1 µg/ml ethidium bromide). Stained agarose gels were analyzed using SYNGENE Bio Imaging Gel Documentation System, for the presence of a fluorescent band of the expected level for both investigated segments.

Purification of PCR amplicons and sequencing

PCR products were purified using the Exo SAP-IT PCR Purification Kit (Applied Biosystems) following the manufacturer’s recommended protocol. Sequencing reactions were performed using Big Dye TM terminator Cycle Sequencing Kit (Applied Biosystems). Sequences were determined in Macrogen company using ABI3730 automated DNA sequencer (Applied Biosystems).

Multiple sequence alignments, diversity analysis and construction of the phylogenetic trees

Homologous sequences for the two investigated segments in the present study (COI and MT-TT&MT-TP) were retrieved from Genbank database of four taxa groups representing three buffalo species from subfamily Bovinea. The four buffalo taxa, their relationships and accession numbers of all retrieved sequences of COI and MT-TT&MT-TP segments from Genbank are presented in Table 2. All retrieved sequences of COI and MT-TT&MT-TP were from whole mitochondrial genome sequences in Genbank for the four buffalo taxa except five sequences (HM461918.1-HM461922.1) that represent mitochondrial sequences for MT-TT&MT-TP genes in Egyptian river buffalo. Multiple alignments of the determined sequences of river buffalo COI and MT-TT&MT-TP in this study and the homologous sequences retrieved from Genbank were performed by CLUSTALW (http://www.genome.jp/tools/clustalw/) program [15].

Table 2

Relationships between the four buffalo taxa groups and Genbank accession numbers of the retrieved sequences of COI and MT-TT&MT-TP.

Subfamily	Genus	Species	Subspecies	Investigated mt genes	Accession number in Genbank
Bovinae	Bubalus	Bubalis (water buffalo)	Bubalis (river buffalo)	COIMT-TT&MT-TP	AF547270.1
				MT-TT&MT-TP	HM461918.1-HM461922.1
			Carabanesis (swamp buffalo)	COIMT-TT&MT-TP	AY702618.1NC_006295.1JN632607.1KX758296.1KX758331.1KX758352.1KX758353.1KX758363.1KX758371.1KX758389.1KX758397.1KX758399.1KX758402.1
		Depressicornis (lowland anoa)		COIMT-TT&MT-TP	EF536351.1NC_020615.1
	Syncerus	Caffer (African buffalo)		COIMT-TT&MT-TP	EF536353.1JQ235505.1-JQ235547.1NC_020617.1

Genetic diversity calculations and phylogenetic trees construction were performed using sequences of COI and MT-TT&MT-PT that were investigated and analyzed in the present work (four haplotypes for COI in the investigated Egyptian and Egyptian-Italian river buffalo samples and five sequences for MT-TT&MT-PT in Egyptian-Italian river buffalo) along with homologous sequences that retrieved from Genbank for the four buffalo taxa groups: river buffalo, swamp buffalo, lowland anoa and African buffalo (Table 2).

Mean intraspecific diversity within samples of each of the tested buffalo taxa groups, mean intraspecific diversity within all taxa and mean interspecific diversity between all taxa were calculated, for COI and MT-TT&MT-TP markers sequences, by MEGA7 software using the same algorithm, the maximum composite likelihood (MCL) [16], that was also used for construction of the phylogenetic trees.

Two phylogenetic trees were created based on COI and MT-TT&MT-TP segments sequences to provide the graphic representation of sequence divergences of the two segments within and among the four tested buffalo groups. Homologous sequences for the COI and MT-TT&MT-TP investigated segments in goat (Capra hircus) whole mitochondrial genome (accession no. GU229278.1) was used in trees construction as an out group. The trees were constructed by Molecular Evolutionary Genetics Analysis version 7.0 (MEGA7) software [17] using Unweighted Pair Group Method with Arithmetic Mean (UPGMA) [18] and maximum composite likelihood (MCL) method [16].

Results

COI sequence analysis and alignment

PCR assay was applied on seventy one DNA samples from Egyptian river buffalo (sixty samples) and Egyptian-Italian river buffalo hybrids (eleven samples) for amplification of 709 bp from the 5′ half of COI gene.

The recovered nucleotide sequence (690 bp) from the COI amplified segment (709 bp) was determined, assembled and manually corrected for all buffalo samples. Sequences were submitted to the Genbank barcode database (accession nos. from KU932052.1 to KU932122.1) and to the BOLD database (accession nos. from GBMIN46384-16 to GBMIN46454-16).

Multiple sequence alignment of COI sequences revealed the presence of five variable nucleotide sites between the investigated sequences in accordance with nucleotide positions in the whole mitochondrial DNA sequence of river buffalo (accession no. AF547270.1). Presence of these five variable nucleotide sites between sequences of the mitochondrial COI investigated segment in the different buffalo samples resulted in segregation of the buffalo samples into four haplotypes. Haplotype 1 was the prevalent haplotype that was shared by 57 samples. These samples comprise 47 samples of the Egyptian river buffalo and 10 samples of the Egyptian-Italian hybrid river buffalo. Haplotype 2 was shared by 12 samples of the Egyptian river buffalo. Haplotype 3 was only presented by one Egyptian river buffalo sample while haplotype 4 was presented by only one sample of the Egyptian-Italian hybrid river buffalo. Nucleotide variations characterizing the four haplotypes are illustrated in Table 3.

Table 3

Nucleotide variations between the four COI haplotypes in accordance with nucleotide positions in the whole mitochondrial DNA sequence of river buffalo (Accession No: AF547270.1). The dots means: identical nucleotides.

Number of segregated haplotypes	Nucleotide positions/Nucleotide variations					Number of samples sharing the haplotype
	5816	5840	5867	5952	6174
Haplotype 1	C	A	A	T	T	57 samples
Haplotype 2	·	G	·	C	C	12 samples
Haplotype 3	·	G	G	C	C	1 sample
Haplotype 4	T	·	·	·	·	1 sample

Multiple alignment of river buffalo COI haplotypes (in the present work) and all the homologous segment sequences in river buffalo, swamp buffalo, lowland anoa and African buffalo from Genbank database revealed the presence of 104 nucleotide variable sites between sequences of these four investigated taxa groups.

Analysis of the COI diversity and barcoding efficacy

The calculated mean ± standard error (SE) and percentage values of intraspecific diversity within each of the tested taxa groups, intraspecific diversity within all taxa and interspecific diversity between the taxa for the investigated COI segment are presented in Table 4.

Table 4

Genetic diversity of the investigated segments COI and MT-TT&MT-TP.

Studied taxa	COI				MT-TT&MT-TP
	Intraspecific diversity		Interspecific diversity		Intraspecific diversity		Interspecific diversity
	Mean ± SE	%	Mean ± SE	%	Mean ± SE	%	Mean ± SE	%
River buffalo	0.005 ± 0.002	0.5	0.049 ± 0.013	4.9	0.000 ± 0.000	0.0	0.002 ± 0.001	0.2
Swamp buffalo	0.004 ± 0.002	0.4			0.005 ± 0.004	0.5
Lowland anoa	0.000 ± 0.000	0.0			0.000 ± 0.000	0.0
African buffalo	0.004 ± 0.002	0.4			0.002 ± 0.002	0.2
All taxa	0.003 ± 0.001	0.3			0.002 ± 0.001	0.2

The phylogenetic tree that based on COI barcode sequences in the four buffalo taxa groups is presented in Fig. 1. Graphic presentation of different sequences in the tree indicated that the four buffalo groups were separated independently from each other and distributed in two clades. First clade included African buffalo sequences and the second clade was re-branched into two sub-clades. First sub-clade included all river buffalo samples and the second sub-clade was re-branched into two branches; one for swamp buffalo and the other one for lowland anoa.

Fig. 1

COI barcode-based phylogenetic tree constructed using sequences from four buffalo taxa groups and a sequence of goat (Capra hircus) as an out group. Hap.: haplotype.

MT-TT and MT-TP sequence analysis and alignment

Mitochondrial DNA segment (228 bp) representing partial sequence of CYTB gene and two successive tRNA genes (MT-TT&MT-TP) was amplified in five Egyptian-Italian river buffalo DNA samples. Sequence of the amplified segment (228 bp) was determined, assembled and manually corrected in all samples. Sequences were further deposited in Genbank database (accession nos.: from MG544232.1 to MG544236.1). Nucleotide sequence of the MT-TT&MT-TP genes (134 bp) was selected from the whole amplified segment (228 bp) in the five buffalo samples. Alignment of MT-TT&MT-TP sequences revealed that the nucleotide sequence of this DNA segment was 100% similar in the five investigated river buffalo samples.

Multiple alignment of MT-TT&MT-TP genes segment in Egyptian-Italian river buffalo samples (present study) and the homologues sequences for this segment in river buffalo, swamp buffalo, lowland anoa and African buffalo that retrieved from Genbank revealed the presence of 7 nucleotide variable sites between sequences of the four taxa groups.

Analysis of MT-TT and MT-TP diversity and barcoding efficacy

The calculated mean ± standard error (SE) and percentage values of intraspecific diversity within each of the tested taxa groups, intraspecific diversity within all taxa and interspecific diversity between the taxa for MT-TT&MT-TP are presented in Table 4. Results indicated that mean intraspecific diversity and mean interspecific diversity had the same value: 0.002 ± 0.001 (0.2%).

Graphic presentation of different sequences in the phylogenetic tree that based on MT-TT&MT-TP sequences in the four buffalo taxa groups (Fig. 2) indicated that only lowland anoa species was separated independently from the other tested species in a separate clade. Graphic presentation also indicated that the three buffalo groups: river buffalo, swamp buffalo and African buffalo shared the same clade which was differentiated into two sub-clades with mixed origins for the three taxa groups.

Fig. 2

MT-TT&MT-TP-based phylogenetic tree constructed using sequences from four buffalo taxa groups and a sequence of goat (Capra hircus) as an out group. Eg-It: Egyptian-Italian hybrid.

Discussion

Currently, the 5′ half of COI gene has been selected as a standard barcode marker for animal groups [1], [2], [3], [4], [19]. Since DNA barcode approach depends on a good library of COI sequences being established prior to testing the species identity [20], enrichment of Genbank database and BOLD system with COI barcode sequences in different animal taxa is important for developing the global biodiversity inventory.

Scarcity of river and swamp buffalo COI barcode sequences in the previously mentioned databases necessitates submission of more COI sequence data for the two subspecies. In the present study, seventy one sequences (690 bp each) form 5′ half of COI gene in Egyptian and Egyptian-Italian river buffalo samples were determined and deposited in Genbank barcode database and BOLD system as a contribution in construction of a public reference library for COI barcode sequence of the river buffalo subspecies. Recently, new 107 sequences for COI barcode in swamp buffalo subspecies are available in Genbank database through the project of sequencing swamp buffalo whole mitochondrial genome [21].

Analysis of river buffalo COI sequences that obtained in the present work revealed the presence of only five variable nucleotide sites between investigated sequences that resulted in segregation of the buffalo samples into four haplotypes. These results indicate that number of nucleotide variations within sequences of mitochondrial 5′ half of COI gene in Egyptian and Egyptian-Italian river buffalo is very low compared to number of nucleotide variations that was detected in other study using other mitochondrial segment. Hassan et al. [11] detected 77 variable nucleotide sites, in the mitochondrial D-loop sequence of 48 Egyptian river buffalo samples, that resulted in segregation of the buffalo samples into 28 haplotypes.

Unlike distinct species, closely related species constitute a great challenge for phylogeny construction and species identification with DNA barcoding due to their often overlapping genetic variation [22]. In the present study, the efficacy of COI barcode sequence was examined in the identification and delineation of three closely related buffalo species: water buffalo (Bubalus bubalis), lowland anoa (Bubalus depressicornis) and African buffalo (Syncerus caffer) that are distributed in two genera of the subfamily Bovinea: genus Bubalus and genus Syncerus. Our study was further expanded to verify the efficiency of the use of COI barcode to differentiate between the two water buffalo subspecies: river buffalo (Bubalus bubalis bubalis) and swamp (Bubalus bubalis carabanesis).

For successful barcoding, sequence variation within species must be adequately low within species and sufficiently high between species so that they can be differentiated from one another [23]. Our present results revealed the presence of marked difference between the mean intraspecific diversity (0.3%) and the mean interspecific diversity (4.9%) of COI barcode sequence. It has been previously reported that >2% sequence divergence is possible to identify and discriminate between animal species [24].

Results obtained from the COI barcode-based phylogenetic tree also revealed that COI sequence was successfully able to delineate the three tested buffalo species (water buffalo, lowland anoa and African buffalo). COI sequence was also efficiently able to delineate the two subspecies of water buffalo (river buffalo and swamp buffalo). Our results are contrary to those reported in two previous studies [12], [13]. The authors of the two studies signified to the presence of intraspecific divergences within water buffalo species and as a result they could not differentiate between samples from the two subspecies. Their findings were based on the available limited number of sequences at that time in Genbank database. On the other hand, COI barcode-based phylogenetic tree constructed in the present study did not reflect real taxonomic relationships between river buffalo, swamp buffalo and lowland anoa. These results are acceptable as DNA barcoding based on a short DNA sequence of the mitochondrial COI gene is a method to correctly identify species and is not an alternative tool for taxonomic classification [5]. Moreover, it has been also reported in a previous study that the phylogenetic tree which based on the whole mitochondrial genome sequence was not able to reflect the real relationships between the same three closely related buffalo taxa groups [25].

The mitochondrial DNA sequence of MT-TT and MT-TP genes has been suggested previously by Oraby et al. [13] to be used as a complementary marker for the COI barcode marker. The authors reported a 100% success of this DNA sequence in the differentiation between water buffalo, cattle, sheep and goat species. In the present study efficacy of the sequence of MT-TT and MT-TP genes was also compared with the efficacy of COI gene barcode in identification and delineation of the investigated species and subspecies. Diversity analysis of MT-TT&MT-TP genes segment revealed that there were no difference between percentages of intraspecific and interspecific divergences of the four buffalo taxa groups, both were equal to 0.2%. This is due to the high similarity between sequences of this segment in the investigated buffalo groups. These results were also reflected on the graphic presentation of different sequences in the phylogenetic tree that was based on MT-TT&MT-TP sequences. The phylogenetic tree indicated that MT-TT&MT-TP sequence was only able to delineate the lowland anoa species and failed to delineate any of the other three buffalo groups: river buffalo, swamp buffalo and African buffalo.

Conclusions

In the present study, we contributed to enrich Genbank barcode database and BOLD system for constructing public reference library for COI barcode sequence in river buffalo. We also reported success of COI gene barcode in identification and delineation of the samples from each of the three closely related buffalo species: water buffalo, lowland anoa and African buffalo. Moreover, COI barcode efficiently discriminated between samples of the two subspecies; river and swamp buffaloes. Results also indicated that mean interspecific diversity of the MT-TT&MT-TP genes segment between buffalo groups is very low and it failed to be used as a barcode marker for differentiation between these closely related buffalo taxa groups.