The Oriental Tiny Frog of the Genus Microhyla Tschudi, 1839 (Amphibia: Anura: Microhylidae) Revealed across Geographical Barriers of the Wallace Line

ABSTRACT The frog genus Microhyla was considered as the South, East, and Southeast Asian frog species. Microhyla orientalis was described in 2013, distributed in Java and Bali, Indonesia. Thenceforth, it was known as the easternmost distribution of this genus within the oriental region, but recently this species was recorded from the Timor Island and Sulawesi on the Wallace regions. We applied molecular analysis to evaluate the taxonomic status and the origin of the Wallacean population. Phylogenetic analysis using the partial 16S mitochondrial gene demonstrated that the Java, Timor and Sulawesi populations were not significantly different from the Bali population. This Wallacean population of M. orientalis was originated from Java and possibly it is accidentally distributed by humans through the expansion of agricultural activity.


INTRODUCTION
Microhyla Tschudi, 1838, a genus of tiny narrow-mouthed frogs that consisted of 46 species which were widely spread from Japan (Ryukyu Islands), southern China across India to Sri Lanka, and through Southeast Asia region into the Indonesian archipelago (Matsui et al. 2005;Kurabayashi et al. 2011; and Sulawesi (identified as Microhyla sp. in Wiantoro et al. 2019) based on morphological characters.
The disjunct distribution of Microhyla orientalis from the western and eastern parts of Wallace line raised a question about the origin of this species. By using a molecular technique, we evaluated the species identification of the Sulawesi and Timor populations then reconstructed its phylogenetic relationships. The estimation on the phylogenetic relationship using molecular technique has succeeded to prove the origin of certain vertebrates (e.g. Moritz et al. 1993;Suzuki et al. 2011;Suzuki et al. 2014;Johnson et al. 2011, Hadi et al. 2020 as well as in the frogs (Kuraishi et al. 2009;Wogan et al. 2016;Reilly et al. 2017;Vences et al. 2017).
Here, we provide the distribution history encompassing the entire populations of M. orientalis in Indonesia. The presumption of M. orientalis genetic variation in the Wallace region, which is possibly associated with natural or human-mediated distribution, was appraised.

Materials
We examined specimens of Microhyla orientalis stored in Museum Zoologicum Bogoriense (MZB), Research Center for Biology, Indonesian Institute of Sciences. A number of 25 partial 16S mtDNA sequences of the M. orientalis were analysed (Table 1). We compared several populations viz. Java, Bali, Timor (specimen from Reilly et al. (2020)), and Sulawesi Island (specimens from Wiantoro et al. (2019) previously identified as Microhyla sp.) ( Figure 1 and Table 1).

Figure 1.
Indonesian Archipelago map of the localities of Microhyla orientalis specimens in Sulawesi (purple), Timor Islands (cyan), Java (green), and type locality, Bali (orange). Outgroups are shown by the dark blue circle. The specimen numbers represent those listed in Table 1. (Map modified from ArcMap 10.7.1, February 5 th , 2021).

Methods
Fragments of the 16S mtDNA (ca. 462) were employing a method as delineated by Matsui et al. ( , 2013. DNA sequences obtained in this study were checked and edited using the ChromasPro software (Technelysium Pty Ltd., Tewantin, Queensland, Australia). The newly M. orientalis sequences were deposited in GenBank with accession numbers MW683205-MW683215 together with those from GenBank (Table 1) aligned applying Clustal W in MEGA X (Kumar et al. 2018). Phylogenetic trees were reconstructed using Neighbor-Joining (NJ), Maximum Likelihood (ML), and Bayesian Inference (BI) analyses. The NJ tree was administered in MEGA X using p-distances with 1000 bootstrap replicates. The Akaike Information Criterion (AIC) performed using Kakusan 4 were used to identify the models of rate evolution for ML and BI analyses (Tanabe 2011). ML analysis was performed using Treefinder ver. March 2011 (Jobb et al. 2004) with general time-reversible (GTR) and a gamma shape parameter (G). with general time-reversible (GTR) and a gamma shape parameter (G). Analysis of the Markov Chain Monte Carlo (MCMC) for the dataset was run for 5 million generations and every 100 cycles, trees were sampled. The convergence of the runs was determined by a split frequency of < 0.01 standard deviations and potential scale reduction factors of ~1.0. We discarded the first 20% of the sampled trees as burn-in and generate a majority-rule consensus tree using the remaining samples. Strong supports of tree nodes were considered when possessing bootstrap values of 70% or more for ML and NJ analyses (Hillis & Bull 1993). The genetic distances of the 16S mtDNA gene were computed using uncorrected p-distances with MEGA X. We regarded tree nodes with BI posterior probabilities values over 0.95 as strongly supported; values between 0.90 and 0.95 were considered as moderately supported; while the lower values were regarded to have no nodal support (Huelsenbeck & Hillis 1993).

Analysis on the structure of populations
Molecular variance (AMOVA) and fixation index (F ST ) (Wright 1951) analysis was performed for the group within intraspecific populations of M. orientalis using the Arlequin v.3.5 programs (set up, 1000 permutations; significance level threshold, α = 0.05). The analyses allowed the approximation of the overall extent of the genetic variation and differentiation level within M. orientalis population. Furthermore, population differentiation and its significance between sampling locations were also calculated using pair-wise estimates (Weir & Cockerham 1984;Excoffier et al. 1992;Weir 1996).

Genetic haplotypes analysis
The median-joining method was generated to build the haplotype network using Network v10.2.0.0 program (Bandelt et al. 1999). Haplotypes distribution for each location was presented in the informative map ( Figure  1) to show the recent genetic connectivity and distributions among populations.

RESULTS AND DISCUSSION
Morphologically, Timor and Sulawesi populations of M. orientalis showed high similarity in appearance with Java and Bali populations. The diagnostic character such as bright orange to pinkish brown on the upper forearm, a vertebral line at the dorsal of the body, and the web formulae of toes confirmed the species identification as M. orientalis (Figure 2). We obtained 462 bp 16S mtDNA fragments of 25 samples including outgroups, 371 nucleotide sites were conserved, 87 site variables, and 28 sites parsimoniously informative. A topology with the highest log likelihood -1223.5300, gamma shape parameter 0.1407, and frequencies of the nucleotides: A=0.333, T=0.251, G=0.197, and C=0.217 resulted from the ML analysis. BI analysis with nucleotide frequencies: A=0.335, T=0.248, G=0.202, C=0.215, and a gamma shape parameter 0.3701. All analyses produced similar topologies, which differed only in several low-supported nodes. Thus, only the ML tree is shown (Figure 3).

Genetic Diversity
Sequences (n= 22) from four localities (Java, Timor, Bali, and Sulawesi) resulting in 462 bp 16S mtDNA was gained. All sequences, including samples from Bali as type locality of M. orientalis resulted in 14 haplotype variations with 17 polymorphic sites ( Table 2).
The contrast of the genetic diversity among person tests of M. orientalis uncovered the nearness of variety (Table 3). Among the samples, it was obtained that the haplotype diversity (Hd) ranged from 0.222 to 1.0 and the nucleotide diversity (π) ranged from 0.000 to 0.015. Samples from Java were observed as the highest genetic diversity (Hd = 1.0; π = 0.015), the second is Timor, which shows a single haplotype and nucleotide diversity (Hd = 1.0; π = 0.000). The M. orientalis population from Sulawesi revealed a fairly high genetic diversity (Hd = 0.933; π = 0.008) and the lowest was from Bali (Hd = 0.222 and π = 0.011).  Table 3. Genetic diversity samples of Microhyla orientalis derived from sample size (n), haplotype number (Hn), Haplotype diversity (Hd), and nucleotide diversity (π).

Populations of Microhyla orientalis
The fixation index (F ST ) and P-values analysis between and within M. orientalis populations (Java, Bali, Timor, and Sulawesi) are presented in Table 4. The comprehensive F ST value within four populations were highly significant (F ST = 0.605; P < 0.001). This was predicted due to multiple subdivisions. Values of the pair-wise F ST between the populations of M. orientalis are provided in Table 5. Significant differentiation between Bali and the other three populations (Java, Sulawesi, and Timor) was detected through the overall pair-wise analysis using the distance method.

Genetic Connectivity of M. orientalis
The haplotype network analysis recognized a complex ancestor or origin of Wallacean M. orientalis population (Timor and Sulawesi) (Fig 4). Two major groups of haplotypes are named clade Bali and Java-Timor-Sulawesi. Bali as the type locality of M. orientalis shows two different haplotypes H1 and H2. Timor (one haplotype) and Sulawesi (five haplotypes) have mixed connectivity with each other, several probabilities of genetic connection of various ancestor points evolved from the Java population (six haplotypes).
Two confined haplotype groups between Bali and Java-Timor-Sulawesi M. orientalis were differentiated by 6 discrete nucleotide bases due to the mutations. Concerning the significant value of F ST among populations (F ST = 0.605; P < 0.001) as well as haplotype network ( Figure 4) and phylogenetic tree (Figure 3).

Discussion
Here, we confirmed the species identification of Yudha et al. (2019), Wiantoro et al. (2019) and Reilly et al. (2020) as Microhyla orientalis genetically . Genetic variation has not been significantly detected among the sundaic and wallacean population of M. orientalis. This assumed that the Wallacean M. orientalis were genetically mixed, feasibly limited in number, and accidentally introduced, as stated earlier by Reilly et al. (2020). It is predicted that oceanic dispersal to adjacent islands must also have occurred in association with human activities, such as the condition on some frogs species in Japan (Ota et al. 2004;Kuraishi et al. 2009 (Wogan et al. 2016;Vences et al. 2017) and Wallacean region (Reilly et al. 2017). Madagascar D. melanostictus was found to be a single origin from the Southeast Asian lineage. The Southeast Asia mainland population shows high haplotypes diversity than the island's population, 51 haplotypes in the mainland, two and four in the island's population (Wogan et al. 2016). As well as, Wallacean D. melanostictus were shared single haplotype originated from Sundaic region of Sumatra and Java (Reilly et al. 2017). Microhyla orientalis from Java showed higher haplotypes diversity than Bali, Timor and Sulawesi. The haplotype network of M. orientalis showed two major haplotype groups, Bali and Java-Timor-Sulawesi, with low genetic diversity. Since the Wallacean population of M. orientalis is closest to Java than Bali, we assumed that this Wallacean population is originated from Java similar to Reilly et al. (2017) finding.
M. orientalis were suspected as an accidentally introduced species by the colonization event back to The Dutch East Indies era in 1905 (Kementerian Desa, Pengembangan Daerah Tertinggal dan Transmigrasi RI 2015a). One of the major sources of people who move to other areas of Indonesia is Java. It was mentioned that Timor and Central Sulawesi were the destinations of the colonization called "transmigrasi" program of the Indonesian government (Kementerian Desa, Pengembangan Daerah Tertinggal dan Transmigrasi RI 2015a; Kementerian Desa, Pengembangan Daerah Tertinggal dan Transmigrasi RI 2015b). Colonization is usually also involved with agriculture and fisheries for living. Possibly , people from Java also brought their cultural agriculture and fisheries to fulfill their basic needs in the new land. It is likely that the species was recently introduced by human activity, and genetic analysis has proven the human interfere of M. orientalis distribution. Moreover, the paddy field is the suitable habitat of M. orientalis (Matsui et al. 2013), the expansion of agricultural activity in the Wallace region might play an important role in distributing this species out from Java.