Zoosyst. Evol. 98 (2) 2022, 275-284 | DOI 10.3897/zse.98.79213 eee BERLIN On the edge of the Shivaliks: An insight into the origin and taxonomic position of Pakistani toads from the Duttaphrynus melanostictus complex (Amphibia, Bufonidae) Daniel Jablonski!, Rafaqat Masroor?, Sylvia Hofmann? 1 Department of Zoology, Comenius University in Bratislava, Ilkovicova 6, Mlynska dolina, 842 15 Bratislava, Slovakia 2 Zoological Sciences Division, Pakistan Museum of Natural History, Garden Avenue, Islamabad 44000, Pakistan 3 LIB - Leibniz Institute for the Analysis of Biodiversity Change, Museum Koenig, Adenauerallee 127, D-53113 Bonn, North Rhine- Westphalia, Germany https://zoobank. org/6623 BB&C-693 D-4836-91 BD-E24FE3AABE2C Corresponding author: Daniel Jablonski (daniel.jablonski@uniba.sk) Academic editor: Umilaela Arifin # Received 11 December 2021 Accepted 31 May 2022 @ Published 13 July 2022 Abstract The common Asian toad Duttaphrynus melanostictus (Schneider, 1799) complex has a wide distribution ranging from western foot- hills of the Himalaya to the easternmost range of the Wallacea, with the evidence of human-mediated introductions to some other ar- eas. In the entire distribution range, the complex is formed by several evolutionary clades, distributed mostly in South-East Asia with unresolved taxonomy. In the northwestern edge of its distribution (Pakistan), the name D. melanostictus hazarensis (Khan, 2001) has been assigned to local populations but its biological basis remained, so far, understudied and unvalidated. Therefore, we re-evaluated the available genetic data (mitochondrial and nuclear) to show the relationships between Pakistani populations (including the type locality of D. m. hazarensis) and others from across the range. Our results showed that Pakistani populations are associated with one, deeply diverged, well-supported and widely distributed clade (so-called Duttaphrynus sp. 1 according to /6S, or clade B based on tRNAGIy-ND3), that has already been detected in previous studies. This clade is further distributed in India, Nepal, Bangladesh, Malaysia, Singapore, and Indonesia and is characterized by a low level of genetic variability. This further suggests that both natural, as well as potential human-mediated dispersal, might have played an important role in setting up the current phylogeographic and distribution pattern of this clade. The clade is deeply divergent from other clades of the complex and represents a taxonomically unresolved entity. We here argue that the clade Duttaphrynus sp. 1/B represents a distinct species for which the name Duttaphrynus bengalensis (Daudin, 1802) comb. nov. is applicable, while the description of D. m. hazarensis does not satisfy the rules of the In- ternational Code of Zoological Nomenclature. Key Words Bufo, Bufonidae, Indian subcontinent, invasions, phylogeography, taxonomy Introduction The distribution of the Duttaphrynus melanostictus (Schneider, 1799) complex (Dubois and Ohler 1999) more or less corresponds to the area of the Oriental zoogeographic realm. It is a widespread group of toads with a geographic range throughout South and South- East Asia, southern China, including islands of Taiwan, Hainan, Sumatra, Java, Borneo, and Sumba (Frost 2021). The group has a high potential to spread into new localities and was artificially introduced to many regions of Asia as well as out of Asia (see Reilly et al. 2017; Vences et al. 2017; Othman et al. 2020; Soorae et al. 2020). These introductions are well-monitored be- cause they can pose a Serious threat to the local native fauna (Pettit et al. 2021). Copyright Jablonski, D et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. 276 The Indian subcontinent and northern Pakistan form the western border of the D. melanostictus complex natural range where its distribution is limited to the edge of Shivalik Hills, i.e. to several districts of Khyber Pakhtunkhwa, Punjab, and Azad Jammu and Kashmir Provinces (Mertens 1969; Khan 1972, 2001; Masroor 2012). Also, the presence of the complex in Sindh Province cannot be excluded (Masroor 2012) as this toad is known from several localities in Gujarat (India) close to the Pakistani border (Sarkar 1984). Such a distribution pattern signifies that D. melanostictus has been recorded in both the Palearctic and Oriental region, which are the main zoogeographic realms of Pakistan, (see Masroor 2012: fig. 20). Surprisingly, populations in the northwestern part of the range have been described as. Duttaphrynus melanostictus hazarensis (Khan, 2001) based on morphological data. According to Khan (2001), specimens of these populations differ from other Asian populations of D. melanostictus by bearing bony ridge, kidney-shaped parotids, subarticular tubercles that are doubled under the penultimate phalanx of all fingers, the absence of rostral ridge, and a light brown dorsal coloration. However, the taxon description and comparison provided by Khan (2001) are very brief. Further, the original description neither provides a voucher nor a catalogue number for the holotype, nor information of the type collection or museum designation. Therefore, it does not fulfil the International Code of Zoological Nomenclature (Art. 16.4.2, ICZN), rendering the subspecific name D. m. hazarensis (Khan, 2001) invalid. In order to simplify the terminology in our text, we will, however, refer to the Pakistani population and related lineage as “hazarensis’’. Due to the topographic heterogeneity and habitat diversity, it might be plausible that Pakistani amphibian populations of different genera with the Oriental affiliation may hold a different evolutionary history compared to their populations from other parts of South-East Asia. The significance of the far-western part of the Himalaya and the Indo-Gangetic plain for understanding the biogeography of taxa distributed across these geomorphological domains has been recently shown in several species groups (see Allopaa, Chrysopaa, Microhyla, Nanorana, Scutiger, Sphaerotheca;, Hofmann et al. 2017, 2019, 2021; Jablonski et al. 2020, 2021). However, except for Khan’s (2001) original work, populations of D. melanostictus from Pakistan, including toads from the type locality where the “hazarensis” population has been described, have not been studied so far in a broader biogeographical context (cf. Hussain et al. 2020; Akram et al. 2021). According to Khan (2001), the type locality has a wider range and is defined as “Ooghi, Mansehra and Datta, Hazara Division, eastern NWFP, Pakistan’, 1.e., covering several villages and districts in northern Pakistan where these toads occur. Recent molecular studies focusing on the genus Duttaphrynus were based on dense sampling, covering a wide range of South-East Asia (mostly), and different DNA markers for D. me/anostictus (see Hasan et zse.pensoft.net Jablonski, D et al.: Duttaphrynus melanostictus in Pakistan al. 2014; Wogan et al. 2016; Reilly et al. 2017; Vences et al. 2017; Othman et al. 2020; Soorae et al. 2020; Bisht et al. 2021). Because the resulting DNA sequence variation showed deep genetic structuring, we reviewed and reanalysed these data together with new sequence data of D. melanostictus from different localities of Pakistan, including the locality of “hazarensis”’ (Datta, Khyber Pakhtunkhwa; 34.2962°N, 73.2649°E). Our main goals of the study are to (1) find out the genetic affiliation of the Pakistani population of D. melanostictus in the context of known molecular phylogeny of the complex, (i1) show phylogeographic relationships of the western populations to others inside the detected clade, and (iil) comment on the validity and taxonomy of the Pakistani taxon in light of molecular phylogeny and nomenclature of the D. melanostictus complex. Materials and methods Sampling and DNA extraction Tissue samples of six individuals of the Duttaphrynus melanostictus complex from northern Pakistan were col- lected from five different localities of the Khyber Pakh- tunkhwa Province during the field studies between 2018 and 2019. Individuals were documented photographically by DJ (Nikon D810, and 105 mm macro lens) and sam- pled for toe or blood before release or being deposited in scientific collections (Suppl. material 1: Table S1). Tissue samples were transferred into 96% ethanol, kept at ambi- ent temperature during the time of the fieldwork, and later stored at -20 °C. DNA was isolated using the DNeasy Tis- sue Kit (QIAGEN), following the manufacturer’s proto- cols, and stored at -20 °C. Amplification and sequencing We amplified a fragment of the nuclear (nDNA) rhodop- sin (Rho) gene (~316 bp), of the /6S rRNA (~560 bp), and approximately 469 bp of mitochondrial DNA (mtDNA) that encodes part of the cytochrome oxidase c subunit III (COITD, tRNA glycine (tRNAGIy), and part of the NADH dehydrogenase subunit 3 (VD3) using the primers and PCR conditions following previous studies (Vences et al. 2017; Hofmann et al. 2019; Bisht et al. 2021). PCR prod- ucts were Sanger-sequenced in both directions by Macro- gen Europe (Amsterdam, The Netherlands; http://www. macrogen-europe.com) and sequence data edited where necessary in Sequencher 5.4 by eye. Dataset and phylogenetic analyses We complemented our new mitochondrial and nuclear DNA sequences dataset with published and unpublished sequence data from GenBank (for details see Suppl. Zoosyst. Evol. 98 (2) 2022, 275-284 material 1: Table S1). We mostly followed datasets pro- vided by Wogan et al. (2016), Vences et al. (2017), and Bisht et al. (2021) as well as others with available se- quences for D. melanostictus. To map the phylogeograph- ic pattern of detected clades, we used published georef- erenced localities of the sequence origin. If coordinates were not available for sequences, we georeferenced them according to the provided name of the locality or speci- men voucher numbers connected to a particular museum collection. In some cases, sequences could not be georef- erenced (e.g., “India”; especially for /6S dataset) or with uncertainty only (e.g., “Uttarakhand”, “Kalimantan”; de- tails in Suppl. material 1: Table S1). For an overall view on the distribution of the complex we downloaded data for D. melanostictus (~20,000 records) available in Glob- al Biodiversity Information Facility (GBIF, accessible via www.gbif.org). New sequences were aligned in MEGA X v.11.0 (Kumar et al. 2018). We followed Vences et al. (2017) and trimmed the mitochondrial sequences, retaining 411 bp corresponding to 71 bp of tRNA-Gly and 340 bp of ND3. We excluded the following three /6S sequences of the data set of Bisht et al. (2021) due to low data quali- ty and potential sequencing errors: EU071741, U52798, EU071765. Also, one sequence (AY458592) was ex- cluded from the data set of Vences et al. (2017) due to uncertain geographic origin. The resulting data sets con- sisted of 150 sequences for /6S (560 bp), 496 sequences for tRNAGIy-ND3 (411 bp), and ten sequences for Rho (300 bp) (Suppl. material 1: Table S1). The latter sequence alignment was used only for the haplotype network con- struction. All new sequences were deposited in GenBank (see Suppl. material 1: Table S1 for accession numbers). Bayesian phylogenetic analyses were carried out with MrBayes v.3.2.6 (Ronquist et al. 2012) using the best-fit model of sequence evolution as determined for 16S (GTR+I+G) by jModeltest v.2.1.10 (Darriba et al. 2012) and tRNAGI/y-ND3 (four partitions: K80, K80+G, HKY+G, GTR+I+G) by PartitionFinder v.1.1.1 (Lanfear et al. 2012). We run MrBayes with a random starting tree for five million generations using four parallel Mar- kov chain Monte Carlo (MCMC) simulations with four chains, sampling trees every 500 generations, until reach- ing an average standard deviation of split frequencies of <0.01. Stationarity and convergence of the runs were confirmed using the software Tracer v.1.7.2 (Rambaut et al. 2018). The first 25% of the samples of each run was discarded as burn-in. For reasons of comparison, we also inferred the BI tree for /6S based on the secondary structure of the ribosom- al RNA using RNAsalsa 0.8.1 (Stocsits et al. 2009). As the initial input, we used the constraint file based on the secondary structures of Bos taurus 16S rRNA which 1s provided with the RNAsalsa package. The sequence was then aligned with our /6S data using the MUSCLE algo- rithm in MEGA X. Starting with this initial alignment and the constraint file, RNAsalsa predicts RNA secondary structure and enhances structural alignment that results Abd in a final multiple sequence alignment together with a consensus structure. Unambiguous stem pairs were then derived based on the consensus structure from RNAsalsa, specified in the MrBayes input file, and analysed under the doublet model (16x16) proposed by Schoninger and von Haeseler (1999). To assess the evolutionary distance between our sam- ples and the Duttaphrynus sister clades recovered from the phylogenetic analysis, we calculated the uncorrected p-distances for the /6S rRNA and tRNAGIy-ND3 data set separately using Mega-X based on 1,000 bootstrap repli- cations, with the pairwise deletion option, and by consid- ering both transitions and transversions. For clades that correspond with the Pakistani population of D. melanost- ictus, and for the Rho sequence data we created a medi- an-joining haplotype network using PopART (Leigh and Bryant 2015). Results The newly generated sequences of /6S, tRNAGIy-ND3, and Rho from Pakistan belong to one, well-support- ed, monophyletic group in all investigated markers (Figs 1-3). This clade is known from previous studies, particularly Wogan et al. (2016), Vences et al. (2017) and Bisht et al. (2021). For /6S, populations from Pakistan correspond with the clade called Duttaphrynus sp. 1 (sensu Bisht et al. [2021]) that includes populations from Bangladesh, India, Nepal, and surprisingly, Malaysia (continental as well as Borneo part; Figs 1A, 2A). This clade has a sister position to clade Duttaphrynus sp. 2 (sensu Bisht et al. 2021) including populations from China, Laos, Malaysia, Myanmar, and Vietnam. Overall, seven clades of the D. melanostictus complex were detected in the /6S phylogeny. Uncorrected p-distances between these clades range from 2.3 to 4.9% (Table 1). The genetic distance between Duttaphrynus sp. 2 and Duttaphrynus sp. 1 clades is 2.7%. The clade Duttaphrynus sp. 1 1s formed by 12 haplotypes in 24 sequences with low internal sequence variation (m=1.0%). Sequences from Pakistan have a closer relationship (one mutation step) to sequences from peninsular India (KU589218, EU071749, EU071752, EU071754, EU071757, KY000453-KY000456) and Nepal (MF319517; Fig. 1B). Other sequences (haplotypes) differ by two to six mutation steps. Surprisingly, two sequences from Malaysia (including Borneo; AB530648, AB331714) are separated by only two mutation steps from Pakistani sequences, or are even identical to some sequences from India or Nepal (Figs 1B, 2A). Significantly, some sequences from India could not be georeferenced and, thus, not displayed on the map in Fig. 2A (see Suppl. material 1: Table S1). However, the overall phylogeographic pattern suggests a distribution of the clade Duttaphrynus sp. 1 in the northern part of the Indian subcontinent with a potential range limit in Myanmar, where it is replaced by the sister clade Dutaphrynus sp. 2. zse.pensoft.net 278 10 samples India 16S rRNA A mtDNA 1 sample Pakistan __— India, Nepal India India, Malaysia B) 0.95/78 Bangladesh Lea India Duttaphrynus bengalensis comb. nov. (Duttaphrynus sp.1) Bangladesh India Malaysia - AB530648, AB331714 Nepal Pakistan 1.00/75 Jablonski, D et al.: Duttaphrynus melanostictus in Pakistan tRNAGIly—ND3 mtDNA C Shahdara, Margalla Hills, Pakistan Cambodia China Lao PDR Madagascar G Malaysia Myanmar Taiwan Thailand UAE Vietnam 0.75/90 Duttaphrynus sp. 2 Malaysia Thailand Laos Vietnam Myanmar China 1.00/82 0.99/84 Indonesia, Malaysia, Singapore 10 samples 0.73/62 1 sample Pakistan 0.62/25 Pakistan India (Western Ghats) Indonesia (Sumatra) (Sulawesi) (Java) (Siberut) 0.90/91 1.00/100 - Indi . cuorivet i, Duttaphrynus melanostictus ? Duttaphrynus parientalis 79 1.00/99 7 . Duttaphrynus brevirostris -/46 Duttaphrynus cf. microtympanum Duttaphrynus melanostictus §-—~ China F Lao PDR > 1.00/97 Myanmar -/22 Vietnam ff D Myanmar \— Malaysia -/98 E Myanmar 1.00/76 C Myanmar 1.00/100 Indonesia Malaysia Pakistan 0.99/93 1.00/89 1.00/100 Duttaphrynus bengalensis ‘Singapore comb. nov. 1.00/100 A Indonesia 0.1 Figure 1. The molecular phylogenetic reconstruction of the Duttaphrynus melanostictus complex on 16S rRNA (A) sensu terminol- ogy of Bisht et al. (2021) and tRNAGIy—-ND3 (C) sensu Vences et al. (2017) together with median-joining networks (B, D) calculated from clades that include sequences of the “hazarensis” toads from Pakistan, now D. bengalensis comb. nov. (Daudin, 1802). The numbers on branches represent Bayesian Posterior Probabilities/bootstraps showing branch support. Sequence details are presented in Suppl. material 1: Table $1. The same colors of clades used for particular mtDNA trees reflect clade affiliation. The pictured individual represents a population from Shahdara, Margalla Hills, Islamabad, Pakistan. Populations of this ‘sp. 1’ clade in Western Malaysia and Borneo (Sarawak) might have been artificially introduced (see Discussion below and Fig. 2A). For tRNAGIy-ND3, Pakistani sequences correspond to the clade B (sensu Vences et al. [2017]) that includes pop- ulations from Indonesia (Kalimantan), West Malaysia, and Singapore (Figs 1C, 2B). This clade has a sister rela- tionship with the clade C distributed in Myanmar. Over- all, we confirmed seven clades of the D. melanostictus complex in the phylogeny of tRNAG/y-ND3. Uncorrected genetic distances between clades range from 2.7 to 9.3% (Table 2). The genetic distance between closely related clade B and C is 5.6%. The clade B includes three haplo- types. Two of them are found in Pakistan and differenti- zse.pensoft.net ated by one mutation step from each other. The third hap- lotype differs by four mutation steps from the Pakistani sequences and can be found in Indonesia (KU183404), Malaysia (KU183466, KU183468-KU 183470), and Sin- gapore (KU183475, KU183476, KU183484). Internal sequence variation in clade B 1s very low (a = 0.3%). Insufficient sampling across the Indian subcontinent did not allow for assessing range limits of the clade B in rela- tion to other clades in the resulting phylogeny. However, detection of the clade B in Singapore, West Malaysia, and Indonesia, (Borneo; Kalimantan, Fig. 2B) suggests an ar- tificial, probably human-mediated introduction. The results on rhodopsin (Rho) confirm mitochon- drial data. It shows that D. melanostictus from Pakistan 2t9 Zoosyst. Evol. 98 (2) 2022, 275-284 Table 1. Uncorrected genetic distances (%) between investigated /6S clades of the Duttaphrynus melanostictus complex as shown in the Fig. 1A. Given are p-distance values (lower-left) and standard error estimates (upper-right). Clade p-distances (%) cf. microtympanum brevirostris parientalis melanostictus 2 melanostictus Duttaphrynus D. bengalensis cf. microtympanum — 0.7 0.7 brevirostris 2.3 - 0.8 parientalis Ae a) = melanostictus ? 2.9 3.4 2S melanostictus 2.8 3.0 29 Duttaphrynus sp. 2 4] 4.3 4.7 D. bengalensis comb. 4.2 3:7 42 sp. 2 comb. nov. (Duttaphrynus sp. 1) 0.8 0.6 0.9 0.9 0.9 0.7 1.0 0.9 0.6 On 0.9 0.9 — 0.6 0.8 0.8 2 - 1.0 0.8 43 4.9 - 0.6 Su 39. De - nov. (Duttaphrynus sp. 1) Table 2. Uncorrected genetic distances (%) between investi- gated tRNAGIy-ND3 clades of the Duttaphrynus melanostictus complex as shown in the Fig. 1C. Given are p-distance values (lower-left) and standard error estimates (upper-right). Clade p-distances (%) A B Cc D E F G A = BH 3A 23 “2g 32) (29 B 93 = 24 29 Qe 3:0; 247 C 93- 360 — =26" 23. 2%. 24 D GES” apie . OPOay | = 1 i hey E 74 49 354 42 - Is? = 6 F Oh Ee) UO PAE BG he G Saw Ot oe CA aS eh represents an independent group (Fig. 3). This group is formed by sequences from Pakistan and differs by three mutation steps from populations from Vietnam (DQ283967; based on the geographic position it proba- bly represents the mtDNA clade Duttaphrynus sp. 2 in 16S, or clade G in tRNAGIy-ND3, respectively) and by two to six mutation steps from D. melanostictus from India (AF249097; without precise locality). The Indian sequence cannot be related to any known mitochondrial phylogenetic clade due to insufficient sampling and un- known geographic origin. The uncorrected genetic dis- tances between Rho sequences of the D. melanostictus complex range between 0.7 and 1.3%. Discussion The most prominent results of our study are the following: (i) the “hazarensis” population from Pakistan belongs to a widely distributed and taxonomically so far unresolved clade (see below), (11) Pakistani mitochondrial sequences represent unique haplotypes with a close relationship to other populations from northern parts of the Indian sub- continent, (111) comparing the overall phylogeographic pattern, the clade Duttaphrynus sp. 1/clade B has been most likely artificially translocated to South-East Asia, probably from the Indian subcontinent, (iv) the name for the clade Duttaphrynus sp. 1/clade B is available and, thus, should be applied (see Taxonomy). These results are thus important for both, evolution and taxonomy of the D. melanostictus complex. However, they still underline previous conclusions that the complex needs comprehen- sive biogeographical research with taxonomic revision based on dense DNA sampling, especially from the In- dian subcontinent (see Wogan et al. 2016; Vences et al. 2017; Bisht et al. 2021). Genetic affiliation and origin of Pakistani populations By re-analysing /6S sequences of the D. melanostictus complex we showed a pattern of three main, geograph- ically widely distributed clades (called as Duttaphrynus melanostictus, Duttaphrynus sp. 1, and sp. 2 sensu Bisht et al. 2021) with potential distribution boundary between them in western Myanmar and Malay Peninsula/Sumatra (Fig. 2). Populations from Pakistan belong to the Dut- taphrynus sp. 1 clade that is distributed in the northern part of the Indian subcontinent (Fig. 1; Bisht et al. 2021). Unfortunately, limited sampling in central India does not allow delimitating the range of this clade in the south- ern parts of India. The Indian subcontinent 1s, however, absolutely crucial to understand the range of this clade as southern India is inhabited by different evolutionary clades of the complex, 1.e. D. brevirostris, D. microtym- panum, and D. melanostictus (sensu Bisht et al. 2021; Figs 1, 2). Re-evaluation of tRNAG/y-ND3 sequence data of the complex revealed similar results as for /6S, show- ing two widely distributed clades (A and G) and four geo- graphically disjunct clades assigned to the northern and western parts of Indochina (C, E, D, F; Wogan et al. 2016; Vences et al. 2017; Figs 1, 2). Clade B (sensu Vences et al. 2017 as a part of the “mainland clade” sensu Wogan et al. 2016), where Pakistani toads belong to, was previous- ly reported from Indonesia (Kalimantan, Borneo), Ma- laysia, and Singapore (Wogan et al. 2016; Vences et al. 2017). However, since these former studies did not focus on this specific clade, no conclusions were drawn about its distribution area. Despite limited data from the Indian subcontinent, the genetic clustering of Pakistani samples within clade B, and the similar phylogenetic placement observed in the /6S tree suggest a wide distribution range of this clade. The low level of intra-clade sequence differ- entiation and star-like pattern of the haplotype networks zse.pensoft.net 280 Jablonski, D et al.: Duttaphrynus melanostictus in Pakistan [165 rRNA Datta, Pakist O D. melanostictus O D. bengalensis comb. nov. (Duttaphrynus sp. 1) O Duttaphrynus sp. 2 O D. brevirostris @ D. microtympanum — Introduced populations? [(RNAGIy-ND3] Le Se ae Be igh <> ze @ cladeA @elade B @cladec O clade D @cladeE Oclade F QOclade G ——> Introduced populations? Figure 2. Phylogeographic reconstruction and distribution of major clades of the Duttaphrynus melanostictus complex across the species range based on /6S rRNA (A) and tRNAG/Iy—ND3 (B) sequences available in GenBank and from own data (see Suppl. ma- terial 1: Table S1). The clade colors follow those of Fig. 1, clade terminology is derived from Bisht et al. (2021; 76S rRNA), and Vences et al. (2017; tRNAGIy—-ND3). The pictured individual represents the population of Duttaphrynus bengalensis comb. nov. (Daudin, 1802) from Datta, Mansehra, Pakistan. Coloration reflect the main mtDNA clades according to Fig. 1. Distribution records (white circles) have been downloaded from GBIF (gbif.org). The red line (panel A) depicts the area considered as “Bengal”, the type locality of Bufo bengalensis. implies a very close relationship between geographically For the clade Duttaphrynus sp. 1/clade B, we revealed a distant populations and potentially rapid expansions to _phylogeographic structure resembling patterns detected in wider areas, supplemented by recent human-mediated in- other widely distributed amphibians of the Indo-Gangetic troductions (see below). plain. These species colonized Pakistan from the Oriental zse.pensoft.net Zoosyst. Evol. 98 (2) 2022, 275-284 Duttaphrynus stomaticus Pakistan KJ532280 Bufotes viridis DQ283940 Duttaphrynus melanostictus India AF249097 281 Datta, Pakistan Duttaphrynus sp. 2 Vietnam DQ283967 Duttaphrynus bengalensis comb. nov. Pakistan ON320395 ON320396 ON320397 ON320398 ON320399 ON320400 Figure 3. Median-joining alleles network of sequences from the Duttaphrynus melanostictus complex and closely related members of western Palearctic Bufonidae based on rhodopsin (Rho) gene available in GenBank and from own data. The colors correspond to those used in Fig. 1. Note, the origin of sequence AF249097 from India cannot be assigned due to unknown geographic details. The pictured individual represents a population of Duttaphrynus bengalensis comb. nov. (Daudin, 1802) from Datta, Mansehra, Pakistan. region, namely Microhyla nilphamariensis Howlader, Nair, Gopalan, Merila, 2015 and Sphaerotheca maskeyi (Schle- ich & Anders, 1998) (Jablonski et al. 2020, 2021). Both Species are wide-ranging, occupying central and north- ern parts of the Indian subcontinent, from eastern India, through Bangladesh, Nepal, central, northern, and west- ern India to Pakistan. Their distribution more or less fol- lows the Indo-Gangetic plain including the Brahmaputra, Yamuna, Ganga, and Indus rivers. Thus, we assume that the clade Duttaphrynus sp. 1/clade B colonized new areas along these riverine corridors (probably from the northern India and Nepal area; see the structure of the haplotype net- work Fig. 1B) and that this clade can also be found tn cen- tral and western India, eventually also southern Pakistan (see Sarkar 1984). However, as shown by the model of the potential distribution of the complex in Asia (see supple- mentary fig. S2 in Vences et al. 2017), the distribution of these toads might be scattered in central India and denser to the region of the Shivalik Hills (see also the georefer- enced records from the GBIF database Fig. 2). The Shiva- liks (also known as the Churia Hills) may provide optimal environmental conditions for these toads (or this specific evolutionary clade). It is a mountain range of the outer Hi- malayas that stretches from Pakistan about 2,400 km east- wards (10-50 km wide) close to the Brahmaputra River, spanning across the northern parts of the Indian subconti- nent. As the western edge of these hills is currently inhabit- ed by the clade Duttaphrynus sp. 1/clade B we can assume that such an environmental corridor (the Siwalik Corridor sensu Khan 2006) also affects the distribution of this spe- cies complex, enabling migration from the east to the west. Recently, two basic biogeographic scenarios have been hypothesized for the origin of different amphibian taxa distributed in northern Pakistan: 1) dispersal of ancestral lineages during the Miocene from East or Southeast China into the far northwestern part of the Himalaya-Tibet orogen along a warm temperate corridor that existed in significant parts of Paleo-Tibet’s interior (Hofmann et al. 2021), this corridor played probably an important role in the initial diversification process of Himalayan faunal organisms (Schmidt et al. 2012), (11) dispersion on widely distributed species characterized by the low level of genetic divergence that followed Shivalik Hills and the Indo-Gangetic plain and reached the territory of Pakistan (Jablonski et al. 2020, 2021) as a result of the Pleistocene or Holocene colonisations (see Othman et al. 2020). Although these are presumable biogeographic models, the second scenario could be now applied for D. melanostictus from Pakistan. It is well known that the D. melanostictus complex is considered an invasive organism, with several evolutionary clades that were introduced out of its native range (Reilly et al. 2017; Vences et al. 2017; Othman et al. 2020; Soorae et al. 2020; Pettit et al. 2021). Evidence for potentially human-mediated introductions is visible from our anal- ysed data, showing large geographic distances between the source and apparently introduced populations (Fig. 2). Consistent with previous studies (Reilly et al. 2017; Vences et al. 2017; Othman et al. 2020; Soorae et al. 2020), proba- bly clades A and G were sources of introduced populations (sensu Vences et al. 2017, Fig. 1). These clades are widely distributed and covered by dense sampling. Moreover, these clades are distributed also in and around populated cities and ports of South-East Asia from where the artificial, hu- man-meditated introduction to other areas is highly evident (even to Madagascar or the United Arab Emirates; Vences et al. 2017; Soorae et al. 2020). On the other hand, limited sampling of the clade B (sensu Vences et al. 2017 which is equal to clade Duttaphrynus sp. 1 sensu Bisht et al. [2021] zse.pensoft.net 282 in the /6S dataset) did not allow to evaluate the origin of this clade in the Malay Peninsula or Borneo. Based on the combination of /6S and tRNAGIy-ND3 data and the dis- tribution pattern of both detected clades (Duttaphrynus sp. 1/clade B) we can hypothesize that Malaysian, Singapore, and Indonesian populations are results of the expansion of the range during ice-ages followed by recent introductions due to human activities (Othman et al. 2020) from a hith- erto unknown region of the Indian subcontinent (Fig. 2). Taxonomy All the hereby presented markers and overall phylogenetic pattern of the clade to which sequences of D. melanostic- tus from Pakistan belong to (see Figs 1—3) support previous studies discussing the unsettled taxonomy of this toad com- plex in its entire distribution (Wogan et al. 2016; Bisht et al. 2021). The taxonomy of this group is difficult to resolve due to missing data from different regions, especially from the western part of the range, and similar morphology. In our study, we provide, for the first time, the evidence that the “hazarensis’ population from Pakistan belongs to the wide- ly distributed clade specified as Duttaphrynus sp. 1 (16S, sensu Bisht et al. 2021), or clade B (tRNAGIy-ND3; sensu Vences et al. 2017), respectively. It implies that Khan (2001) correctly observed morphological differences on these toads (see Introduction) but the lack of broader comparative ma- terial probably did not allow to compare specimens from Pakistan with other specimens, especially from the Indian subcontinent. The description, moreover, is formally inval- id, since the holotype lacks information about the name and location of the repository where it has been deposited (Art. 16.4, ICZN 1999). Khan (2001) compared the Pakistani population with island populations from Taiwan (which, ge- netically, probably belong to Duttaphrynus sp. 2 clade in the 16S tree, and clade G in the tRNAGIy-ND3 tree, respective- ly) and Sri Lanka from where genetic data are lacking. Following the recommendation of Dubois and Ohler (1999), 1.e. “Jf further studies demonstrate that a distinct species, close to Bufo melanostictus, occurs in Bengal, the name Bufo bengalensis Daudin, 1802 is available for it’, and based on the /6S phylogeny, the geographic distribution of Duttaphrynus sp. | clade (sensu Bisht et al. 2021), and the Principle of priority rule (ICZN 1999), we suggest that the name Duttaphrynus bengalensis (Daudin, 1802) comb. nov. (type locality “Bengale”, neotype MNHN 4967; for description see Dubois and Ohler 1999) is applicable for this clade. “Bengale” is equal to the area that is divided between Bangladesh and several states in neighbouring India (states West Bengal, Jharkhand, Bihar and the Cachar region of Assam matching the distribution of the Duttaphrynus sp. 1 clade (Fig. 2A). According to our data, the clade is native to Bangladesh, India, Nepal and Pakistan, and apparently introduced to Indonesia, Malaysia, and Singapore. (see above and Fig. 2). Morphologically, the descriptions of D. m. hazarensis (Khan, 2001) and D. bengalensis (Daudin, 1802) comb. nov. are almost identical in several aspects: the interorbital space is concave in both “hazarensis” and zse.pensoft.net Jablonski, D et al.: Duttaphrynus melanostictus in Pakistan D. bengalensis comb. nov., the space is larger than the width of the upper eyelid and internarial distance. Tympanic oblique, its diameter tn both these species 1s about two thirds of the eye diameter; the tympanic-eye distance 1s about less than one half of its distance. Moreover, the larger size of “hazarensis”’ (max. SVL of female up to 111.0 mm) and D. bengalensis comb. nov. (max. SVL of female neotype up to 133.2 mm) differentiates this species from D. melanosctictus (Schneider, 1799) specimens from India whose maximum SVL is 71.6 mm (lectotype ZMB 3462; Dubois and Ohler 1999). The snout length in D. bengalensis comb. nov. and “hazarensis’ is longer than the horizontal diameter of eye (vs. shorter in D. melanostictus). The canthus rostralis in D. bengalensis comb. nov. and “hazarensis” is sharp, with a ridge in comparison to D. melanostictus whose canthus rostralis is angular. The tympanum diameter in D. bengalensis comb. nov. and “hazarensis” is oval, oblique, its diameter about two-thirds of the eye diameter, whereas it is vertical and its diameter is more than half of the eye diameter in D. melanostictus. The tympanum-eye distance in D. bengalensis comb. nov. and “hazarensis” is about less than one half of its diameter, whereas it is less than one sixth of its diameter in D. melanostictus. Thus we recommend to assign the name D. bengalensis (Daudin, 1802) for the clade Duttaphrynus sp. 1/clade B for future studies. Acknowledgements We would like to thank our colleagues, friends and local people for their help during the field work in Pakistan. For technical and laboratory support we thank Martina Law- son and Jana Polakova. For recommendations regarding the taxonomy we thank Frank Tillack. We are also grate- ful to S. R. Chandramouli and one anonymous reviewer for their useful comments on the manuscript. This work was supported to DJ by the Slovak Research and Devel- opment Agency under the contract APVV-19-0076, and by the German Research Foundation (DFG, grant no. HO 3792/8-1) to SH. References Akram A, Rais M, Lopez-Hervas K, Tarvin RD, Saeed M, Bolnick DI, Cannatella DC (2021) An insight into molecular taxonomy of bufo- nids, microhylids, and dicroglossid frogs: First genetic records from Pakistan. Ecology and Evolution 11(20): 14175-14216. https://do1. org/10.1002/ece3.8134 Bisht K, Garg S, Akalabya Sarmah AND, Sengupta S, Biju SD (2021) Lost, forgotten, and overlooked: Systematic reassessment of two lesser-known toad species (Anura, Bufonidae) from Peninsular In- dia and another wide-ranging northern species. Zoosystematics and Evolution 97(2): 451-470. https://doi.org/10.3897/zse.97.61770 Darriba D, Taboada GL, Doallo R, Posada D (2012) jModelTest 2: More models, new heuristics and parallel computing. Nature Methods 9(8): 772-772. https://doi.org/10.1038/nmeth.2109 Dubois A, Ohler A (1999) Asian and Oriental toads of the Bufo melanostictus, Bufo scaber and Bufo stejnegeri groups (Amphibia, Zoosyst. Evol. 98 (2) 2022, 275-284 Anura): A list of available and valid names and redescription of some name-bearing types. Journal of South Asian Natural History, Colombo 4: 133-180. Frost DR (2021) Amphibian Species of the World: an Online Reference. Version 6.1 (23 November 2021). Electronic Database accessible at https://amphibiansoftheworld.amnh.org/index.php. American Mu- seum of Natural History, New York, USA. GBIF.org (4 December 2021) GBIF Occurrence Download: https://do1. org/10.15468/dl.tdbg7m Hasan M, Islam MM, Khan MMR, Igawa T, Alam MS, Djong HT, Kura- bayashi A, Joshy H, Sen YH, Belabut DM, Kurabayashi A, Kuramoto M, Sumida M (2014) Genetic divergences of South and Southeast Asian frogs: a case study of several taxa based on 16S ribosomal RNA gene data with notes on the generic name Fejervarya. Turkish Journal of Zoology 38(4): 389-411. https://doi.org/10.3906/zo0- 1308-36 Hofmann S, Stock M, Zheng Y, Ficetola FG, Li JT, Scheidt U, Schmidt J (2017) Molecular phylogenies indicate a paleo-Tibetan origin of Himalayan lazy toads (Scutiger). Scientific Reports 7: 3308. https://doi.org/10.1038/s41598-017-03395-4 Hofmann S, Baniya CB, Litvinchuk SN, Miehe G, Li JT, Schmidt J (2019) Phylogeny of spiny frogs Nanorana (Anura: Dicroglossidae) supports a Tibetan origin of a Himalayan species group. Ecology and Evolution 9(24): 14498-14511. https://doi.org/10.1002/ece3.5909 Hofmann S, Jablonski D, Litvinchuk SN, Masroor R, Schmidt J (2021) Relict groups of spiny frogs indicate Late Paleogene-Early Neo- gene trans-Tibet dispersal of thermophile faunal elements. PeerJ 9: e11793. https://doi.org/10.7717/peerj.11793 Hussain S, Bukhari SM, Javid A, Hussain A, Rashid M, Ali W (2020) Molecular identification of genus Duttaphrynus from Punjab, Pakistan. Mitochondrial DNA. Part B, Resources 5(3): 3218-3220. https://doi.org/10.1080/23802359.2020.1810143 ICZN [International Commission on Zoological Nomenclature] (1999) International Code of Zoological Nomenclature, 4° Edn. London, International Trust for Zoological Nomenclature, [xxix +] 306 pp. Jablonski D, Khan MA, Masroor R (2020) The genus Microhyla (Anura: Microhylidae) in Pakistan: species status and origins. Zootaxa 4845(2): 293-296. https://doi.org/10.11646/zootaxa.4845.2.11 Jablonski D, Masroor R, Hofmann S (2021) Revisited molecular phylogeny of the genus Sphaerotheca (Anura: Dicroglossidae): the biogeograph- ic status of northernmost populations and further taxonomic changes. Diversity (Basel) 13(5): e216. https://doi.org/10.3390/d13050216 Khan MS (1972) The “commonest toad” of West Pakistan and a note on Bufo melanostictus Schneider. Biologia (Lahore) 18: 131-133. Khan MS (2001) Notes on cranial-ridged toads of Pakistan and descrip- tion of a new subspecies (Amphibia: Bufonidae). Pakistan Journal of Zoology 33: 293-298. Khan MS (2006) Amphibians and reptiles of Pakistan. Krieger Publish- ing Company, Malabar, Florida, 311 pp. Kumar S, Stecher G, Li M, Knyaz C, Tamura K (2018) MEGA X: Mo- lecular evolutionary genetics analysis across computing platforms. Molecular Biology and Evolution 35(6): 1547-1549. https://doi. org/10.1093/molbev/msy096 Lanfear R, Calcott B, Ho SY, Guindon S (2012) Partitionfinder: Com- bined selection of partitioning schemes and substitution models for phylogenetic analyses. Molecular Biology and Evolution 29(6): 1695-1701. https://doi.org/10.1093/molbev/mss020 Leigh JW, Bryant D (2015) PopART: Full-feature software for haplo- type network construction. Methods in Ecology and Evolution 6(9): 1110-1116. https://doi.org/10.1111/2041-210X.12410 283 Masroor R (2012) A contribution to the herpetology of northern Paki- stan: The amphibians and reptiles of Margalla Hills National Park and surrounding regions. Society for the Study of Amphibians and Reptiles (SSAR), Ithaca, New York and Chimaira Buchhandels- gesellschaft mbH, Edition, Frankfurt, 217 pp. Mertens R (1969) Die Amphibien und Reptilien West-Pakistans. Stutt- garter Beitrage zur Naturkunde 197: 1-96. Othman SN, Chen Y-H, Chuang M-F, Andersen D, Jang Y, Borzée A (2020) Impact of the Mid-Pleistocene revolution and anthropogenic factors on the dispersion of Asian Black-spined toads (Duttaphrynus Animals (Basel) 10(7): e1157. https://doi. org/10.3390/an110071157 Pettit L, Somaweera R, Kaiser S, Ward-Fear G, Shine R (2021) The im- pact of invasive toads (Bufonidae) on monitor lizards (Varanidae): melanostictus). An overview and prospectus. The Quarterly Review of Biology 96(2): 106-125. https://doi.org/10.1086/714483 Rambaut A, Drummond AJ, Xie D, Baele G, Suchard MA (2018) Posterior summarization in Bayesian phylogenetics using Tracer 1.7. Systematic Biology 67(5): 901-904. https://doi.org/10.1093/sysbio/syy032 Reilly SB, Wogan GOU, Stubbs AL, Arida E, Iskandar DT, McGuire JA (2017) Toxic toad invasion of Wallacea: A biodiversity hotspot characterized by extraordinary endemism. Global Change Biology 23(12): 5029-5031. https://doi.org/10.1111/gcb. 13877 Ronquist F, Teslenko M, Van der Mark P, Ayres DL, Darling A, Hohna S, Larget B, Liu L, Suchard MA, Huelsenbeck JP (2012) MrBayes 3.2: Efficient Bayesian phylogenetic inference and model choice across a large model space. Systematic Biology 61(3): 539-542. https://doi.org/10.1093/sysbio/sys029 Sarkar AK (1984) Ecological studies on the amphibians of Gujarat. Bul- letin Zoological Society India 6: 87-93. Schmidt J, Opgenoorth L, Holl S, Bastrop R (2012) Into the Hima- layan exile: The phylogeography of the ground beetle Ethira clade supports the Tibetan origin of forest-dwelling Himalayan species groups. PLoS ONE 7(9): e45482. https://doi.org/10.1371/journal. pone.0045482 Schoninger M, von Haeseler A (1994) A stochastic model and the evolu- tion of autocorrelated DNA sequences. Molecular Phylogenetics and Evolution 3: 2—240. https://doi.org/10.1006/mpev.1994.1026 Soorae PS, Frankham GJ, Mohamed AA (2020) The first record of the Asian common toad Duttaphrynus melanostictus Schneider, 1799 in Abu Dhabi, United Arab Emirates. BioInvasions Records 9(2): 434-443. https://doi.org/10.3391/bir.2020.9.2.29 Stocsits RR, Letsch H, Hertel J, Misof B, Stadler PF (2009) Accurate and efficient reconstruction of deep phylogenies from structured RNAs. Nucleic Acids Research 37(18): 6184-6193. https://doi. org/10.1093/nar/gkp600 Vences M, Brown JL, Lathrop A, Rosa GM, Cameron A, Crottini A, Dolch R, Edmonds D, Freeman KLM, Glaw F, Grismer LL, Litvinchuk S, Milne MG, Moore M, Solofo JF, Noél J, Nguyen TQ, Ohler A, Randrianantoandro C, Raselimanana AP, van Leeuwen P, Wogan GOU, Ziegler T, Andreone F, Murphy RW (2017) Tracing a toad invasion: Lack of mitochondrial DNA variation, haplotype origins, and potential distribution of introduced Duttaphrynus melanostictus in Madagascar. Amphibia-Reptilia 38(2): 197-207. https://doi.org/10.1163/15685381-00003 104 Wogan GO, Stuart BL, Iskandar DT, McGuire JA (2016) Deep ge- netic structure and ecological divergence in a widespread human commensal toad. Biology Letters 12(1): e20150807. https://doi. org/10.1098/rsbl.2015.0807 zse.pensoft.net 284 Supplementary material | Table SI Authors: Daniel Jablonski, Rafaqat Masroor, Sylvia Hofmann Data type: excel file Explanation note: Dataset of sequences used for the study. Copyright notice: This dataset is made available under the Open Database License (http://opendatacommons. org/licenses/odbl/1.0/). The Open Database License (ODbL) is a license agreement intended to allow us- ers 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. Link: https://do1.org/10.3897/zse.98.79213.suppl1 zse.pensoft.net Jablonski, D et al.: Duttaphrynus melanostictus in Pakistan