Kikihia (Dugdale 1971) is a monophyletic endemic New Zealand genus of cicadas (Buckley et al. 2002, Arensburger et al. 2004a). Currently, there are 14 described species and 16 yet-to-be-described species and subspecies (Dugdale 1971; Fleming 1984; Marshall et al. 2008, 2011; Larivière et al. 2010). Species of this cicada genus are found on the three main islands of New Zealand including North, South, and Stewart islands as well as smaller near-by islands and the more distantly located Kermadec and Norfolk islands. Kikihia is widespread with species found in most terrestrial habitats including grass, shrub, and forest ecosystems from lowland to subalpine regions of New Zealand.

Phylogenetic and phylogeographic studies using mitochondrial and nuclear gene sequences have provided insight into the evolutionary history of the genus (Arensburger et al. 2004a,b; Marshall et al. 2008, 2009, 2011). Extensive field study and molecular investigation into this genus identified 20 contact zones between different species that are possible sites of hybridization (Marshall et al. 2011). Individuals found in these proposed hybrid zones possess morphological and acoustic traits that are intermediate between parental species found away from the hybrid zones. Although species in the genus are estimated to have diverged more than 6 million years ago, hybridization seems confined to species that diverged 2 million years ago or less (Marshall et al. 2008).

Here, we present six newly discovered microsatellite markers to investigate the population dynamics, levels of gene flow, and impact of introgression on species boundaries throughout this diverse cicada genus. Because processes related to speciation are difficult to observe after the fact, the study of the interaction at contact zones between incipient species pairs or between recently diverged species provides insight into the speciation process (Barton and Hewitt 1989, Harrison 1993, Mallet 2007, Hewitt 2011). We have developed six novel highly polymorphic microsatellites that amplify in seven species of Kikihia; Kikihia “northwestlandica,” Kikihia “southwestlandica,” Kikihia muta (Fabricius, 1775), Kikihia angusta (Walker, 1850 ), Kikihia “tuta,” Kikihia “nelsonensis,” and Kikihia “murihikua.” Data for pure populations (as far from hybrid zones as possible) are used to test the microsatellite markers' ability to differentiate these species.

Materials and Methods

Specimen Collection. Cicada specimens were identified in the field based on morphology and song traits. Taxonomic descriptions are pending for many Kikihia species as denoted by informal names in quotes. Populations located as far from putative hybrid zones as possible were chosen as species-typical or “parental” populations. Subsequently, 24–36 individuals per species from 1 to 5 largely species-typical populations were investigated to test these new microsatellite markers (Table 1). Whole body specimens were placed in 95% ethanol or three legs were removed and stored in 95% ethanol and the bodies were pinned. All ethanol specimens are stored at -20°C.

Microsatellite Development. Nine species of Kikihia from throughout the genus were used to develop microsatellites: K. “southwestlandica,” K. muta, K. angusta, K. “nelsonensis,” Kikihia cutora cutora, Kikihia scutellaris, Kikihia “peninsularis,” Kikihia horologium, and Kikihia “aotea western.” Genomic DNAwas extracted from cicada legs using the Qiagen DNeasy Blood & Tissue kit (Qiagen, Boston, MA). RNAwas removed with 20 µg RNAase A (New England Biolabs, Ipswich, MA) per sample during DNA purification. The 454 GS FLX Rapid Library Preparation kit (454 Life Sciences, Branford, CT) with Rapid Library multiplex identifiers (MIDs) were used for one specimen from each of the nine species. MIDs are unique sequences incorporated in the 454 primers for the purpose of identification. Sequencing was done unidirectionally. Sample preparation and sequencing were performed according to the manufacturer's instructions (454 Life Sciences, Branford, CT). One library (K. muta) was sequenced on 1/16th of a GS FLX Titanium PicoTiter Plate. The eight remaining libraries were pooled and sequenced on one region of a two-region GS FLX Titanium PicoTiter Plate. GS De Novo Assembler v2.3 was used to align reads with the heterogeneous and large genome options. This resulted in 3,683 contigs with an average size of 743 bp and a range of 100–5,909 bp. MsatCommander (Faircloth 2008, Abdelkrim et al. 2009), a program that searches for di-tri-, and tetranucleotide repeats with sufficient flanking sequence for primer design, was used to scan contigs. Primer3 (Untergrasser et al. 2012) was used to design primers. In total, 30 primer pairs were screened. Six primer pairs were identified representing variable microsatellite loci that amplified in all species. These primer pairs were tested in multiple Kikihia populations.

Table 1.

Kikihia species and populations used for microsatellite discovery

Fragment Analysis. Polymerase chain reaction (PCR) primer pairs were tested with the forward primer tagged with a M13 sequence on the 5'-end and the reverse primers for each locus plus a 6-FAM-labeled M13 primer (Schuelke 2000). PCR amplifications were performed in 15 µl reactions containing 0.4 U Ex Taq DNA polymerase (TaKaRa Biomedical Inc., Otsu, Shiga, Japan), 1× PCR buffer, 0.3mM dNTP mix, 3 pmol forward primer, 12 pmol reverse primer, 12 pmol M13 primer with a 5' 6-FAM fluorescent label reactions, and up to 100 ng template DNA. PCR conditions were as follows: an initial 94°C for 5 min followed by 30 cycles of 94°C for 30 s, 57°C for 45 s, 72°C for 45 s, followed by 8 cycles of 94°C for 30 s, 53°C for 45 s, 72°C for 45 s, followed by a final extension of 72°C for 10 min. PCR was diluted 1:9 in HiDi Formamide (Invitrogen Life Technologies, Grand Island, NY) and run on an ABI 3130xl DNA sequencer (Applied Biosystems, Grand Island, NY). Alleles were designated according to amplicon size relative to LIZ 600 size standard (Invitrogen Life Technologies, Grand Island, NY). PCR forward primers that amplified specimens from across the Kikihia phylogeny and revealed variable microsatellites were then fluorescently labeled with the G5 (Invitrogen Life Technologies, Grand Island, NY) labels (Table 2). Testing of Kikihia populations for the six described primer pairs were multiplex amplified using the Qiagen Type-it microsatellite PCR kit according to the manufacturer's instructions (Qiagen, Boston, MA), including an annealing temperature of 57°C for all multiplexed PCRs.

Microsatellite Analysis. GeneMarker v2.2.0 (SoftGenetics, LLC, Kalamazoo, MI) was used to visualize and score alleles.

Files were converted between different formats using CONVERT v1.31 (Glaubitz 2004). Arlequin v3.5.1.2 (Excoffier et al. 2005) was used to assess the number of alleles per locus, the allele frequencies, and expected and observed heterozygosity. Deviations from Hardy-Weinberg equilibrium (HWE) for each locus in each population were tested with Genepop v4.2.1 (Rousset 2008) using theMarkov chain method. Microchecker v2.2.3 (van Oosterhout et al. 2004) was used to test for the presence of null alleles and large allele dropout.Microsatellite neutrality was tested using the Fst-outlier method implemented in LOSITAN (Beaumont and Nichols 1996, Antao et al. 2008). LOSITAN was run for 100 simulations using the neutral mean Fst and force mean Fst settings for a 0.99 confidence interval and infinite allelesmodel.

Results

Primer sequences, microsatellite repeat motifs, PCR amplicon sizes, and the number of alleles are provided in Table 2. Only one marker, A553, showed stutter patterns that made scoring of some individuals difficult. The largest peak was used, and individuals that were ambiguous because of stutter were marked as missing data. The other five markers had minimal stuttering that did not influence peak calling. The number of alleles per locus varied from 14 to 47 in a total of 213 individuals from seven species. Table 3 summarizes the genetic diversity estimated for each locus in each population for seven species. Observed heterozygosity ranged from 0 to 1 and the expected heterozygosity ranged from 0.129 to 0.945. There was a total of 62 private alleles with all but one population, OL.NLW having at least one private allele. The HWE test, which indicates deviations between observed and expected heterozygosity showed significant deviations in three loci for one to two populations each after a Bonferroni correction (a=0.00055). Three loci were in HWE in all populations. Null alleles were detected in two of the six loci, A1267 and M2333, which had overlapping amplicon ranges and were multiplexed in the same reaction. Null alleles in these two markers are likely due to some interference between them in the electropherograms. Using the Fst-oulier method implemented in LOSITAN, one locus was identified as being a candidate of balancing selection, M2333. This may be due because of the presence of null alleles, which would likely reduce heterozygosity and make this locus. The A1267 locus was unusual because over 50% of the specimens tested were homozygous for one very common allele, whereas the other loci were much more heterogeneous. Some populations failed to amplify one of the multiplex reactions.

Table 2.

Characterization of the microsatellite markers

Table 3.

Characterization of microsatellite markers for each population

Discussion

Next-generation sequencing using 454 pyrosequencing was employed in the development of novel microsatellite markers for the New Zealand cicada genus Kikihia. This technology is becoming more commonly used in microsatellite discovery from a range of taxa (Abdelkrim et al. 2009, Ekblom and Galindo 2010, Gardner et al. 2011). This is the first published discovery of microsatellite markers for any New Zealand cicadas. Unlike most studies that are focused on intraspecific diversity, we are focused on intra- and interspecific diversity in a genus that is approximately 10–12 million years old with the earliest extant species split more than 6 million years ago (Marshall et al. 2008). Microsatellite markers that were variable and cross-amplified in all seven species of Kikihia are presented here.

Three of the six loci were found to violate the HWE in no more than 2 of the 15 populations tested. Heterozygote deficiency increases due to factors such as inbreeding, population stratification, null alleles, and genotyping errors. It is not surprising that species of Kikihia display high levels of inbreeding or population stratification since field observations suggest low dispersal rates. However, these were not explicitly tested here.

The large number of potential hybrid zones between multiple species of Kikihia makes this genus uniquely suited for the study of hybrid zones. Many of these potential hybrid zones are between nonsister species that vary in their divergence times from less than 1 million years to more than 3 million years (Marshall et al. 2008, 2011). The markers developed in this study will be useful for investigations of the evolutionary past and future of this interesting species radiation.