Fifteen into Three Does Go: Morphology, Genetics and Genitalia Confirm Taxonomic Inflation of New Zealand Beetles (Chrysomelidae: Eucolaspis)

Eucolaspis Sharp 1886 is a New Zealand native leaf beetle genus (Coleoptera: Chrysomelidae: Eumolpinae) with poorly described species and a complex taxonomy. Many economically important fruit crops are severely damaged by these beetles. Uncertain species taxonomy of Eucolaspis is leaving any biological research, as well as pest management, tenuous. We used morphometrics, mitochondrial DNA and male genitalia to study phylogenetic and geographic diversity of Eucolaspis in New Zealand. Freshly collected beetles from several locations across their distribution range, as well as identified voucher specimens from major museum collections were examined to test the current classification. We also considered phylogenetic relationships among New Zealand and global Eumolpinae (Coleoptera: Chyrosomelidae). We demonstrate that most of the morphological information used previously to define New Zealand Eucolaspis species is insufficient. At the same time, we show that a combination of morphological and genetic evidence supports the existence of just 3 mainland Eucolaspis lineages (putative species), and not 5 or 15, as previously reported. In addition, there may be another closely related lineage (putative species) on an offshore location (Three Kings Islands, NZ). The cladistic structure among the lineages, conferred through mitochondrial DNA data, was well supported by differences in male genitalia. We found that only a single species (lineage) infests fruit orchards in Hawke’s Bay region of New Zealand. Species-host plant associations vary among different regions.


Introduction
It has been estimated that about 86% of extant species on Earth are yet to be described even after 250 years of taxonomic classification [1], creating impediment to many areas of biological research. It is also increasingly apparent that adherence to strictly defined species concepts is Insects Adult Eucolaspis beetles were collected from various locations and host plants throughout New Zealand (Fig 1 and S2 Table) and preserved in 95% ethanol. In addition, we examined representative named Eucolaspis specimens in the New Zealand Arthropod Collection (NZAC, Landcare Research Ltd., Auckland) including available type material. Specimens of related taxa, Atrichatus ochraceus, A. aenicollis, Peniticus sp. and Pilacolaspis sp. in the Entomology Research Museum at Lincoln University, Lincoln (LUNZ) were also examined.

External morphology
The main characters Broun [10][11][12][13] used to delineate species were body size, body colour, pronotum shape, pronotum size, density of pronotal punctures, and density of elytral punctures. Shaw [15], whose study was based almost entirely on reexamination of Broun's specimens at the British Museum of Natural History (BMNH), primarily used external shape and pronotal puncture density. However, neither Broun nor Shaw quantified the variation in morphology within and between putative species, but rather used relative estimates to indicate variation. It is clear that some of the characters, such as body colour, shape and size of pronotal punctures can vary greatly within populations (P. Doddala, pers. obs.).
In this study, we quantified a set of the external morphological characters previously used (Fig 2A). External morphological characters were recorded from randomly selected individuals from each locality / sample using a digital camera (Moticam 2000 2.0 MP USB 2.0; Motic Group Co., Ltd.) fitted to a dissecting microscope (Zeiss Stemi 2000-c; Carl Zeiss, Inc.). Motic Images Plus v.2.0 (Motic Group Co., Ltd.) was used to record measurements from those images.

Morphology of male genitalia
The morphology (structure, shape and size) of internal genitalia was studied in a subset of randomly selected male beetles among fresh specimens used for external morphology. Male beetles were soaked in cold 10% potassium hydroxide for 12 hours, rinsed thoroughly in 70% ethanol followed by rinsing in dH 2 O, and then soaked for an hour in hydrogen peroxide [23]. The clearing procedure was repeated as necessary, and the cleared genitalia ( Fig 2B) examined under a dissecting microscope with measurements taken using Motic Images Plus software.

DNA extraction
Total genomic DNA was extracted from selected beetles using either a salting-out extraction method [24] with excised legs or the QIAGEN DNeasy blood and tissue kit (QIAGEN N.V.) with whole body samples. Extraction of DNA from dry museum specimens was carried out in a dedicated Ancient DNA laboratory at Massey University, Palmerston North. DNA extractions were checked for quantity and quality by gel-electrophoresis and spectrophotometry (NanoDrop; Thermo Fisher Scientific Inc.). One percent agarose gels with SYBR Safe DNA gel stain (Life Technologies Corp.) in TAE buffer (Tris-HCl, glacial acetic acid, EDTA and H 2 O) were used for electrophoreses.

Data analysis
Morphometric data were tested for distributional normality using multivariate procedures. Morphometric differences between the two sexes were assessed using a t-test. Stepwise discriminant analysis was used to identify which variables (characters) contributed significantly to delineation of sample classes. Subsequently, canonical discriminant analysis was performed using the variables identified by stepwise discriminant analysis, to verify similarity / diversity of samples grouped according to ecological region, host plant, genetic lineage or genitalia shape. Sample locations were assigned to recognised New Zealand ecological regions to test for association of taxa and environment [28]. The "Axial ranges" ecological region was represented by a single sample locality and so was excluded from analyses. A 95% level of confidence was used as a significance level for all the statistical analyses. All analyses were performed using SAS v.9.2 (SAS Institute, 1992).

Morphometric analysis
Morphological variation in fresh samples. Body length (BL) in fresh beetles varied from 2.69 mm to 4.45 mm (mean 3.56 mm) whereas body width (BW) varied from 1.54 mm to 3.16 mm (mean 2.14 mm) (n = 135). Punctures were denser on the pronotum than on elytra or head in all insects. Punctures on pronotum (PPD) varied in density from 160 to 810 per mm 2 , whereas punctures on head (HPD) varied from 20 to 320 per mm 2 . Elytra were less densely punctured at 50 to 180 punctures per mm 2 . There was noticeable sexual dimorphism in body shape with male beetles significantly smaller and more slender than female beetles, and having longer antennae (Table 1). There was, however, no significant difference between male and female beetles in the density of punctures on head, pronotum and elytra.
Punctures on pronotum (PPD) and head (HPD) were the two characters contributing to the separation among beetles of different ecological regions of New Zealand. Beetles from Northern North Island and Central Volcanic Plateau were morphometrically similar to each other (p = .197), while Leeward districts samples were morphologically distant from other regions (p < .001) ( Fig 3A). Canonical variable 1 (Can1) explained about 91% of variation among the ecological regions ( Fig 3A). PPD contributed most to the Can1, while HPD contributed most to canoncial variable 2 (Can2). PPD, pronotum length (PL), body length (BL) and anterior elytral puncture density (AEPD) were the characters that significantly separated beetles from different host plants. Samples from apple and blackberry appeared to cluster together, while samples from manuka were very diverse with no clustering. PPD contributed more to the separation than the other three characters (F (7, 123) = 12.3, p < .001). Can1 and Can2 together explained about 89% variation among the beetles from different host plants.
Museum samples. Morphometric analysis of identified specimens in museum collections showed no distinct clusters. Overlap of data from different species, such as E. vittiger, E. colorata and E. brunnea suggested poor phenotypic separation of current species (Fig 3B). Individuals assigned to E. picticornis were morphologically highly variable and did not cluster together. The four E. montana paratypes from Broun's collection varied considerably, highlighting Table 1. Sexual dimorphism in New Zealand Eucolaspis beetles. Data from representative individuals among fresh beetle samples collected throughout New Zealand. BL Body length, BW body width, EL elytra length, EW elytra width, AL antennae length, PL pronotum length, HPD head puncture density, PPD pronotal puncture density, AEPD anterior elytral puncture density, PEPD posterior elytral puncture density. instability of the existing taxonomy ( Fig 3B). Among the ten morphological characters measured, only elytral width (EW), puncture density on head (HPD) and posterior elytral region (PEPD) contributed to significant variation among "species". PEPD contributed the most  variation (F (10, 58) = 11.77, p < .001). Can1 explained about 73% of variation among "species", whereas Can2 explained about 19% of variation ( Fig 3B).

Genetic, ecological and geographic associations
The 117 aligned mtDNA COI sequences (617 bp) from mainland New Zealand Eucolaspis comprised 39 haplotypes, with an additional haplotype identified from Three Kings Islands specimens. Haplotype diversity (Hd ± S.D.) was 0.97 ± 0.01, and nucleotide diversity per site (P ± S.D.) was 0.0634 ± 0.0046. The alignment contained 129 variable positions and 97 parsimony informative sites. Reconstructed phylogeny of these haplotypes, using Three Kings Islands haplotype (HapTK-NZAC) as an outgroup, showed three well-supported lineages in the mainland ingroup (Lineage 1, 2 and 3) (Fig 4). Phylogenetic inference using different methods (Minimum Evolution, Maximum Likelihood and Bayesian inference) yielded near identical topologies, with minor variation in placement of haplotypes within lineages. Species delimitation analysis using Geneious conducted on a Bayesian inference tree, confirmed the monophyly of the three lineages with the probability of correct identification of an unknown specimen by the sequence tree ranging from 0.87 to 0.98 ( Table 2). The overall, mean genetic distance (p-distance ± standard error) among haplotypes was 0.068 ± 0.006 (measured using MEGA6). Lineage 2 had the highest within-group mean genetic distance (0.018 ± 0.004) compared to the other two lineages (Lineage 1: 0.012 ± 0.002; Lineage 3: 0.007 ± 0.002). Intra-lineage pairwise genetic distances ranged from 0.1% to 3% whereas inter-lineage pairwise genetic distances ranged from 8% to 12.7% (Fig 5). Mean net genetic distance (P-distance) measured as sequence divergence between lineages varied from 7.3% (lineages 1 and 2) to 10% (lineages 1 and 3) ( Table 3). Lineages 1 and 2 were genetically more similar to each other than to Lineage 3 at this locus ( Table 2). Phylogenetic analysis using entire CO1 fragment (~1400 bp) for respresentative individuals from these three lineages also conferred similar evolutionary relationships among the lineages (S1 Fig).
The spatial distribution of the three Eucolaspis lineages showed clear demarcation between East and West, with lineage 1 occupying mostly Eastern areas (Fig 6). Where multiple beetle specimens were available we found that two lineages were sometimes represented at the same locality, whereas one locality (Torere-Bay Of Plenty) had representatives of all three lineages (n = 5). Lineage 2, which had more within-lineage genetic diversity than the other two, had wide geographic distribution and occured frequently in sympatry with one or both lineages 1 and 3. Among the host plants sampled in mainland New Zealand, blackberry (Rubus fruticosus) and manuka (Leptospermum scoparium)were used by all three Eucolaspis lineages, apple (Malus domesticus) was used only by lineages 1 and 2, and kanuka (Kunzea ericoides) was used by lineages 2 and 3.

Male genitalia, morphology and molecular data
Three forms of male genitalic appendage or aedeagei, which differed primarily in the shape of the tip of aedeagus proper (type 1 -apiculate and tapered apically, type 2 -apiculate and broad apically, and type 3 -not apiculate and subacute apically) were found in a sample of 60 male beetles (Fig 4 and S2 Fig). Aedeagus type 3 was most different from the other two types, lacking a well-defined beak (tip). Individuals that belonged to mtDNA genetic lineage 1 possessed type 1 aedeagei, individuals of lineage 2 possessed type 2 and individuals of lineage 3 possessed type 3 aedeagei. An exception was in two males from a single locality in Nelson (collected on apples) that had type 1 aedeagei and belonged to lineage 3. Lineages 1 and 2 were genetically more similar to one another, and their aedeagei appeared to be relatively similar (Fig 4 and S2 Fig). In morphometric analysis, elytra length (EL), pronotum length (PL) and pronotal punctures (PPD) were the only characters that differed among the males with different aedeagus types (Fig 3C). Can1 explained 71% variation, while Can2 explained 29% variation (Fig 3C). Males with aedeagus type 2 differed from the other two groups in having longer elytra and pronotum and lesser density of punctures on pronotum (squared Mahalanobis distance: between types 1 and 2 = 27.12, p < .001; between types 2 and 3 = 24.10, p < .001). Males with aedeagei type 1 and 3 mainly differed from one another in terms of puncture density on pronotum (squared Mahalanobis distance = 5.98, p < .001).
PPD, HPD, AEPD and EW differed significantly among beetles belonging to different haplotype lineages (Fig 3D). Can1 explained 86.4% variation among the lineages, while Can2 explained 13.6% (Fig 3D). Along Can1 (X-axis) PPD, HPD and AEPD were higher among Table 2. Species delimitation analysis confirms monophyly of the thee mainland New Zealand lineages of Eucolaspis. Inter Dist closest = mean pairwise tree distance between the members of the focal species and members of the next closest species; P ID(strict) = mean probability of correctly identifying an unknown specimen of the focal species using placement on a tree sequence; Av (MRCA) = mean distance between the most recent common ancestor of a species and its members; P (randomly distinct) is the probability that a lineage has the observed degree of distinctiveness due to random coalescent processes. Input tree was constructed by Bayesian inference method using GTR+G+I model.

Relationship with other Eumolpinae genera
In addition to our data, COI and 18S rDNA sequences corresponding to 9 different global Eumolpinae genera were obtained from GenBank. In the resulting mtDNA COI phylogeny constructed using Bayesian inference (Fig 7), Eucolaspis and Atrichatus formed a monophyletic lineage. New Zealand Peniticus was clearly more distant to Eucolaspis than Atrichatus, and indeed the current data yield a polytomy comprising Eucolaspis and Atrichatus. Other New Zealand genera in the subfamily Eumolpinae appear to be closely related to the Eucolaspis lineages. Analysis of rRNA 18S sequences constructed using the Maximum Likelihood criterion placed the New Zealand Eucolaspis lineages in a well-supported lineage along with another unidentified Eumolpinae taxon from New Caledonia (S3 Fig). This topology suggests that New Zealand Eucolaspis is loosely separated from other pacific Eumolpinae genera.

Discussion
Mitochondrial DNA sequences, male genitalia and morphometric data provide strong evidence for just three mainland New Zealand species of Eucolaspis (lineages 1, 2 and 3) and a probably fourth on the Three Kings islands (HapTK-NZAC). Phylogenetic analysis showed well-supported lineages, sufficiently distinct to be consistent with different species. There was no support for a larger number of Eucolapsis species in mainland New Zealand. The smallest inter-lineage genetic distance in Eucolaspis at this locus was about 8%, whereas mean interlineage genetic distance was about 10%; the inter-specific (8-12.7%) and intra-specific (0.1-3%) genetic distances did not overlap. There was a prominent gap between the intralineage and interlineage pairwise distances (Fig 6), indicating that taxonomic division is not being arbitrarily imposed on a continuous distribution of diversity. A similar gap was reported in Crioceris (Coleoptera: Chrysomelidae), where the maximum intraspecific genetic distance was about 2.5% and the interspecific distances ranged from 16.9 to 20.3% [39]. Similarly, in Arsipoda (Coleoptera: Chrysomelidae) the interspecific genetic distances varied from 8.1% to 14.4% while intraspecific genetic distances were much smaller (0.3-0.6%) [40]. Although this single locus evidence is insufficient on its own for taxonomic distinction [41,42], it is notable that just three mainland lineages are indicated rather than 15 [11][12][13] or 5 [15]. Analyses of phylogenetic relationships among New Zealand Eumolpinae genera suggest that the genera Eucolaspis Sharp and Atrichatus Sharp are more closely related to each other, than either is to Peniticus Sharp and this confirms the suggestions of Broun [11] and Shaw [15]. A fourth genus, Pilacolaspis Sharp, could not be included in this study as the only available specimens were old and did not yield amplifiable DNA.
Analysis of morphological characters from previously identified museum voucher specimens was not consistent with existing classification. Body size (length and width) of the beetles, the main characters that Broun [10][11][12][13] used in addition to body colour to describe many of his 13 Eucolaspis species, did not differ significantly among the randomly selected sample of different species named voucher specimens. Instead, other characters such as the width of elytra and puncture density on head and posterior elytra partitioned the species into clusters ( Fig  3B). This provided a good impartial test of existing Eucolaspis taxonomy. Overlap of morphology of specimens supposedly representing different described species (such as E. vittiger, E. colorata and E. brunnea) supports in part the synonymy proposed by Shaw [15], although he based his inference on a different set of characters. Our examination of beetles from different genetic lineages in regards to the shape of the punctures, the main character Shaw (1957) used to delineate species, suggested that this character is highly inconsistent. The shape of the punctures varied among individuals within a population, and differences among individuals of different populations (and lineages) showed no consistency. However, puncture density on pronotum, head and elytra, characters also used by Shaw and Broun in species descriptions, displayed consistent differences among the genetic lineages.
Male genitalic shape and morphometric data coincide with genetic data, reiterating three mainland New Zealand lineages. Genitalic shape in males was consistent with a shared common ancestor of lineages 1 and 2. Variation in the shape of male genitalia also indicated that these reproductive structures are under evolutionary selection and this may reflect reproductive isolation, especially in sympatric populations. Reproductive isolation mechanisms such as variation in size and shape of cerci of male grasshoppers (Parapodisma setouchiensis and P. subastris) [43] and difference in cuticular hydrocarbon profiles that act as sex pheromones in leaf beetles Chrysochus auratus and C. cobaltinus [44] have been reported in sympatric populations. However, there is no information on reproductive isolation and / or incompatibility between Eucolaspis "species".
Our results showed that only one lineage (putative species)-Eucolaspis lineage 1 infests apple orchards in Hawke's Bay, New Zealand, while apples elsewhere in the country (e.g. Nelson) are infested by beetles of a different lineage. Eucolaspis feed on many different native and exotic plant species in New Zealand [20], and the wide range of host plants contributed to our sample of the three mainland New Zealand lineages suggests they are polyphagous and all could infest exotic fruit crops.
The three Eucolapsis lineages (putative species) were partitioned into North-West and South-East populations, and this was especially apparent in the Leeward Districts ecological region of New Zealand, which was occupied by lineage 1 (Fig 6). The Leeward districts, which is the driest region included in our sampling, are separated geographically from the rest of the country by the Ruahine and Tararua axial ranges in the North Island and the Southern Alps in the South Island. These ranges may act as a physical barrier, limiting mobility, however, Palmerston North (Windward districts) beetles were genetically similar to Hawke's Bay (Leeward districts) populations, suggesting that contact and dispersal between regions is possible. We do not know if this dispersal is due to discontinuity in the ranges, or anthropogenic, or refelcts intermediate environmental conditions in this area. A similar East-West partitioning of distribution has been suggested in other New Zealand invertebrates including Onychophora [45], Paryphanta snails and corophiid amphipods (in [46]).
The genetic, genitalic and morphometric data utilized in our study complement each other but are not mutually exclusive, and therefore, integrative taxonomy is possible in this genus. Such an integration of different characters provides reliable taxonomic decisions [5], that reflect evolution [2]. Congruence of different types of data has been reported in many recent studies of Coleoptera (e.g., [47,48]). We conclude that there are only three putative species in mainland New Zealand unless others are very scarce or isolated. This is unlikely as our pattern of sampling encompassed the areas used to provide specimens for most of the earlier descriptions [10][11][12][13], which used material from just a few isolated locations. We also sampled through the North Island and in the Nelson-Marlborough and Canterbury regions of the South Island. Eucolaspis beetles are scarce or absent south of Canterbury region of New Zealand [49].
We therefore propose three mainland Eucolaspis taxa, distinguished by haplotype lineage, aedeagus shape, puncture density (on pronotum, head and anterior elytra) and elytra width. Beetles that belong to lineage 1 are distinguished morphologically by having denser puncturation (on pronotum, head and anterior elytra) and narrower elytra than the other two lineages. Given the available data, we propose the following names for Eucolaspis lineages as being appropriate: lineage 1 -Eucolaspis puncticollis (Broun 1880), based on resemblance of aedeagus tip shape with that described by Shaw [15]; lineage 2 -E. picticornis Broun 1893, based on comparison of 18S rDNA data with that of BMNH voucher 69636 GenBank accession DQ337133; lineage 3 -E. jucunda (Broun 1880), based on analogous aedeagus tip shape in Shaw [15] and congruence of 18S rDNA with that of BMNH voucher 696321 GenBank accession DQ337120.
Shaw [15] suggested E. picticornis as a junior synonym of E. brunnea (Fabricius, 1781), however, we feel that E picticornis is more appropriate to use in application to our data. Although E. brunnea (F., 1781) is the earliest name, we did not see the type material and cannot confirm that E. brunnea is consistent with the data we have gathered. In addition, there is a long-standing homonymy between New Zealand E. brunnea (originally described by Fabricius as Chrysomela brunnea, and later moved to Colaspis by White [9]), and North American "grape colaspis" Colaspis brunnea (Fabricuis, 1798) (originally described as Galleruca brunnea, and moved to Colaspis in 1801). To add more confusion, both species are horticultural pests; the North American grape colaspis has been sometimes referred to as the "bronzed beetle" (due to brown colour), and has been indexed in the American Review of Applied Entomology under an incorrect name of Eucolaspis (see discussion in Barber [50]). Although NZ Eucolaspis brunnea is the senior homonym, the name Colaspis brunnea is widely used for the North American species and the homonymy remains unresolved. The higher order evolutionary relationships of Eucolaspis (inter-generic and intra-subfamilial) within New Zealand, Pacific, and the world need to be investigated further.   Broun (1880Broun ( , 1893Broun ( , 1903 and Shaw (1957). (PDF) S4 Table. GenBank (NCBI, USA) accession numbers (GI and Version) for global Eumolpinae taxa 18S rDNA sequences used in the current study. (PDF)