Short Communication
Short Communication
The Caucasian Toad, Bufo verrucosissimus (Pallas, 1814) in the Levant: evidence from mitochondrial DNA
expand article infoDaniel Jablonski, Riyad A. Sadek§
‡ Comenius University in Bratislava, Bratislava, Slovakia
§ American University of Beirut, Beirut, Lebanon
Open Access


We present the first molecular analysis of the Bufo bufo (Linnaeus, 1758) complex population in the Levant. This complex contains four species distributed through the Western Palearctic. Based on the analysis of the widely used mitochondrial marker 16S rRNA, populations from the Levant belong to B. verrucosissimus (Pallas, 1814) and have a close affiliation to populations from the Transcaucasian region rather than southern Turkey. Such identified relationships may imply past rapid colonisation from north to south and support the need for further research.

Key Words

common toad, Bufonidae, mtDNA, biogeography, taxonomy, Middle East

The Caucasian Toad, Bufo verrucosissimus (Pallas, 1814), is a member of the Western Palearctic B. bufo (Linnaeus, 1758) complex along with three other species: B. bufo, B. eichwaldi and B. spinosus (García-Porta et al. 2012; Recuero et al. 2012; Arntzen et al. 2013). Some authors consider B. verrucosissimus as a subspecies of B. bufo due to overlapping distribution, allele-sharing and possible hybridisation events (García-Porta et al. 2012). However, individual-based analyses of enzyme data, presented by Arntzen et al. (2013), support B. verrucosissimus as a full species until better evidence for its taxonomic position is presented. This species is larger in size than other members of the complex (70–190 mm snout-vent length), has a rugged skin texture, different colouration and marked sexual dimorphism in body size (Litvinchuk et al. 2008). The Caucasian Toad prefers wet and shaded sites in mountain coniferous, mixed and deciduous forests (Dufresnes 2019) in the Transcaucasian region (Armenia, Azerbaijan, Georgia, southern Russia and north-eastern Turkey). However, recent molecular data also confirmed this species in western and southern Turkey (García-Porta et al. 2012; Arntzen et al. 2013).

Although information about B. bufo complex in the Levant region (Syria, Lebanon) is absent in recent phylogeographic studies (García-Porta et al. 2012; Recuero et al. 2012), its occurrence had been documented there in the past (Disi and Böhme 1996; Hraoui-Bloquet et al. 2001). Potential isolated distribution of the complex within southern Turkey to northern Israel was also demonstrated through paleoclimatic reconstructions by García-Porta et al. (2012; fig. 7, 125 p.). However, recent information concerning their distribution is very incomplete due to ongoing political turmoil and the existence of fragmented and low density populations of these toads in the region.

We, therefore, document the species’ status in the Levant, specifically in Lebanon. The presence of B. bufo complex in this country was first recorded from the locality Jannat-Kartaba near Nahr Ibrahim (Abraham River) valley (Hraoui-Bloquet et al. 2001; Suppl. material 1: Fig. S1) and represents the southernmost known distribution limit for the complex in the eastern part of the range. However, for nearly two decades, no information about these toads was reported in literature as they are rarely encountered and are restricted to several river valleys. Although Kuzmin (2008) expected the presence of B. verrucosissimus in Lebanon due to its close proximity to known populations in the Anatolian part of Turkey, the taxonomic status (possible independent phylogenetic clade), genetic affiliation and origin (possible existence of the relict population) of this toad population in the Levant are unknown. Therefore, we conducted the first molecular-based study analysing the species/population affiliations of this toad.

Total genomic DNA was extracted from eight blood samples (see Forzán et al. 2012 for details on blood sampling) collected on 29 March 2018 in one population from south-central Lebanon (Moukhtara, around Awali river, 33.66N, 35.61E, 689 m a.s.l.; voucher specimens were not collected; Suppl. material 1: Fig. S2), using the Qiagen DNeasy Blood and Tissue Kit, following the manufacturer’s protocols. The majority of sequences for the B. bufo complex in GenBank (Suppl. material 1: Table S1) was available for the mitochondrial DNA (mtDNA) fragment of the commonly used 16S ribosomal RNA (16S rRNA). We used this standard phylogenetic marker in amphibian DNA barcoding and species classification (Vences et al. 2012) for DNA amplification of the Lebanese population (primers 16Sar-L and 16Sbr-H, using the laboratory protocol of Palumbi et al. 1991). PCR products were sent to Macrogene Inc. (Amsterdam, The Netherlands) for sequencing. For basic species affiliation, the new sequences were checked via BLAST ( and combined with available 16S rRNA GenBank sequences of B. verrucosissimus published in previous studies by Kutrup et al. (2006), Van Bocxlaer et al. (2009), Recuero et al. (2012), García-Porta et al. (2012), Litvinchuk et al. (2012) and Liedtke et al. (2017; Suppl. material 1: Table S1 and see the comment and the map of fig. 3, 1206 p. in Arntzen et al. 2013). The 16S rRNA fragment (~550 bp) was aligned using Clustal W algorithm (Thompson et al. 1994) as implemented in BioEdit (Hall 1999). The alignment was checked by eye and low quality ends trimmed. Ambiguously aligned regions/gaps were ignored for the analysis. We used a network approach (Posada and Crandall 2001) to infer inter-population relationships. We constructed a haplotype network for B. verrucosissimus and allopatric representatives of its related species (JN647044B. bufo, JN647239B. eichwaldi, JN647207B. spinosus) as outgroups, using the 95% limit of parsimony (TCS algorithm; Clement et al. 2000). This was implemented and visualised in the software PopArt ( The DnaSP 5.10 (Librado and Rozas 2009) was used to estimate the number of haplotypes (h), haplotype diversity (Hd), number of segregating sites (S), average number of nucleotide differences (k), nucleotide diversity (π) and Watterson’s theta per site (θW) for all available sequences of the species. The same programme was used to estimate uncorrected p-distances between haplotypes. The eight novel sequences were deposited in GenBank under accession numbers MN094870-MN094877.

Overall 32 sequences (16S rRNA) covering the species range in Turkey and Transcaucasia were combined with eight new ones from Lebanon (Fig. 1A, B; Suppl. material 1: Table S1). The network analysis supports the inclusion of the Lebanese population to B. verrucosissimus which is the first molecular-based confirmation of this species in Lebanon and in the Levant. As indicated by the network analysis (Fig. 1B), the sequences obtained from eight specimens form three different haplotypes in the studied Lebanese locality. The maximal distance amongst these haplotypes is two mutational steps (Fig. 1B). These haplotypes are only one to three mutation steps distant from the geographically most common haplotype (e.g. FJ882807) of the species occurring in Georgia and Russia (0.2 to 0.6% in p-distances) and surprisingly two to five steps (0.4 to 1.0%) from haplotypes detected in localities of southern Turkey (JQ348811 from Alanya and AY840240-1 from Mersin; Fig. 1A, B). Our analysis showed 12 detected haplotypes overall, forming a star-like pattern with a nucleotide diversity of 0.34% (Fig. 1, Suppl. material 1: Table S2).

Figure 1. 

(A) Map displaying the localities of genetically investigated (16S rRNA) populations of Bufo verrucosissimus. The expected range of the species, according to Kutrup et al. (2006), García-Porta et al. (2012), Litvinchuk et al. (2012), Recuero et al. (2012) and Arntzen et al. (2013), is highlighted in orange. The newly investigated population is represented by a red circle with a black dot in the middle. Geographically closest (Alanya and Mersin, Turkey) and furthest (Kyurdzhinovo, Russia) investigated localities, relative to Lebanon, are marked by GenBank accession numbers. The question mark represents parts where species distribution needs further investigation. (B) Haplotype network inferred from a 16S rRNA showing the intraspecific relations of B. verrucosissimus. Black dots inside the red circle represent population from Lebanon and arrows following geographically closest and furthest localities given to Lebanon as presented in Fig. 1A. Sequences of other species of the B. bufo complex were used as the outgroup. The pictured specimen originates from the locality Moukhtara, Lebanon.

The Levant region is known as a speciation centre with endemic phylogenetic lineages and cryptic taxa hidden within morphologically similar sister species that also have fully or partially isolated distributions (e.g. Van Riemsdijk et al. 2017, Dufresnes et al. 2019a, b). Our results are, therefore, surprising in that the Levant population (i) does not form any significantly divergent lineage as would be expected from their likely isolated distribution there, (ii) is genetically very close to populations from the Transcaucasian region rather than the population from southern Turkey (Fig. 1). The northernmost (Kyurdzhinovo, Russia; JQ348764 and JX218095) and southernmost genetically investigated localities (Moukhtara, Lebanon) are separated by approximately 1300 km (straight line distance). However, the genetic distance between these localities is only one or two mutations steps (average 0.4%), similar to some other haplotypes inside the Transcaucasian region (Fig. 1). In the context of geographical distance, one mutation step between the most common Transcaucasian haplotype (accession number, e.g. FJ882807) and the haplotype from Alanya, Turkey (JQ348811) is also noteworthy. The star-like pattern of haplotypes thus suggests a probable origin in Transcaucasia and species radiation from northern to western and southern regions of the current range through Anatolia and the Levant. Such results imply most likely past rapid but different colonisation scenarios. As we used here only mtDNA barcoding, that does not resolve strictly maternal inheritance, the lack of recombination or possible massive mitochondrial introgression (Ballard and Whitlock 2004), denser sampling in Anatolia and Transcaucasia with a combination of genome-wide multilocus analyses and compilation of phenotypic datasets is needed to accurately delineate the phylogeography of B. verrucosissimus and origin of populations in the Eastern Mediterranean.

Lebanon is a mountainous country with deep, wet valleys that could regionally provide suitable microrefugia for B. verrucosissimus. This species prefers wet and forested areas (Kuzmin 2008) and has a clear preference for river basins with particularly rocky riparian habitats. However, related Bufotes toads are also found here but in various different habitats. Bufo verrucosissimus is probably less adapted to the drier climate of the Levant than genus Bufotes. This could be a reason for its rarity in the region and its restriction to the western slopes of Mt. Lebanon which are wetter than other regions in Lebanon. Additional distribution research is needed to identify other localities for B. verrucosissimus in the Levant (see Tron and Rocha 2005) in order to establish local conservation priorities for these southernmost populations of the species.


We thank Daniel Gruľa for his help with the fieldwork in Lebanon, Jana Polakova for laboratory work, Chace Holzheuser for language improving and anonymous reviewers for their important suggestions. This study was supported by the Slovak Research and Development Agency under contract no. APVV-15-0147.


  • Disi AM, Böhme W (1996) Zoogeography of the amphibians and reptiles of Syria, with additional new records. Herpetozoa 9: 63–70.
  • Dufresnes C (2019) Amphibians of Europe, North Africa and the Middle East. Bloomsbury Wildlife, 224 pp.
  • Dufresnes C, Mazepa G, Jablonski D, Sadek RA, Litvinchuk SN (2019a) A river runs through it: tree frogs genomics supports the Dead Sea Rift as a rare phylogeographic break. Biological Journal of the Linnean Society 128: 130–137.
  • Dufresnes C, Strachinis I, Suriadna N, Mykytynets G, Cogălniceanu D, Székely P, Vukov T, Arntzen JW, Wielstra B, Lymberakis P, Geffen E, Gafny S, Kumlutas Y, Ilgaz C, Candan K, Mizsei E, Szabolcs M, Kolenda K, Smirnov N, Geniez P, Lukanov S, Crochet P-A, Dubey S, Perrin N, Litvinchuk S, Denoël M (2019b) Phylogeography of a cryptic speciation continuum in Eurasian spadefoot toads (Pelobates). Molecular Ecology 28: 3257–3270.
  • Forzán MJ, Vanderstichel RV, Ogbuah CT, Barta JR, Smith TG (2012) Blood Collection from the Facial (Maxillary)/Musculo-Cutaneous Vein in True Frogs (Family Ranidae). Journal of Wildlife Diseases 48: 176–180.
  • García-Porta J, Litvinchuk SN, Crochet PA, Romano A, Geniez PH, Lo-Valvo M, Lymberakis P, Carranza S (2012) Molecular phylogenetics and historical biogeography of the west-palearctic common toads (Bufo bufo species complex). Molecular Phylogenetics and Evolution 63: 113–130.
  • Hall TA (1999) BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symposium Series 41: 95–98.
  • Hraoui-Bloquet S, Sadek R, Geze R (2001) Les Amphibiens du Liban: inventaire, répartition géographique et altitudinale. Bulletin de la Société Herpétologique de France 99: 19–28.
  • Kutrup B, Yilmaz N, Canakci S, Belduz A, Doglio S (2006) Intraspecific variation of Bufo bufo, based on 16S Ribosomal RNA sequences. Amphibia-Reptilia 27: 268–273.
  • Liedtke HC, Müller H, Hafner J, Penner J, Gower DJ, Mazuch T, Rödel MO, Loader SP (2017) Terrestrial reproduction as an adaptation to steep terrain in African toads. Proceedings of the Royal Society B, Biological sciences 284: 20162598.
  • Litvinchuk SN, Borkin LJ, Skorinov DV, Rosanov JM (2008) A new species of Common Toads from the Talysh Mountains, south-eastern Caucasus: genome size, allozyme, and morphological evidences. Russian Journal of Herpetology 15: 19–43.
  • Litvinchuk SN, Mazepa GO, Kami HG, Auer M (2012) Taxonomic status and distribution of common toads in Iran. Herpetological Journal 22: 271–274.
  • Palumbi SR, Martin A, Romano S, McMillan WO, Stice L, Grabowski G (1991) The simple fool’s guide to PCR, version 2. Zoology Department, University of Hawaii, Honolulu, 45 pp.
  • Recuero E, Canestrelli D, Vörös J, Szabó K, Poyarkov NA, Arntzen JW, Crnobrnja-Isailović J, Kidov AA, Cogalniceanu D, Caputo FP, Nascetti G, Martínez-Solano I (2012) Multilocus species tree analyses resolve the radiation of the widespread Bufo bufo species group (Anura, Bufonidae). Molecular Phylogenetics and Evolution 62: 71–86.
  • Thompson JD, Higgins DG, Gibson TJ (1994) Clustal W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Research 22: 4673–4680.
  • Van Bocxlaer I, Biju SD, Loader SP, Bossuyt F (2009) Toad radiation reveals into-Indian dispersal as a source of endemism in the Western Ghats-Sri Lanka biodiversity hotspot. BMC Evolutionary Biology 9: 131 pp.
  • Van Riemsdijk I, Arntzen JW, Bogaerts S, Franzen M, Litvinchuk SN, Olgun K, Wielstra B (2017) The Near East as a cradle of biodiversity: a phylogeography of banded newts (genus Ommatotriton) reveals extensive inter- and intraspecific genetic differentiation. Molecular Phylogenetics and Evolution 114: 73–81.
  • Vences M, Nagy ZT, Sonet G, Verheyen E (2012) DNA barcoding amphibians and reptiles. In: Kress WJ, Erickson DL (Eds) DNA Barcodes: Methods and Protocols. Methods in Molecular Biology (vol. 858). Humana Press, Totowa, 79–107.

Supplementary material

Supplementary material 1 

Tables S1–S2, Figures S1–S2

Daniel Jablonski, Riyad A. Sadek

Data type: docx file

Explanation note: Table S1. Dataset of 16S sequences of Bufo verrucosissimus used in this study. Table S2. Summary of genetic polymorphism of Bufo verrucosissimus 16S rRNA sequences. Figure S1. The map of Lebanon with an overview of Bufo verrucosissimus distribution records in the country. Figure S2. Adult individuals of Bufo verrucosissimus from locality Mukhtara, Lebanon.

This dataset is made available under the Open Database License ( The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this Dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited.
Download file (1.50 MB)
login to comment