Alternative DNA Markers to Detect Guam-Specific CRB-G (Clade I) Oryctes rhinoceros (Coleoptera: Scarabaeidae) Indicate That the Beetle Did Not Disperse from Guam to the Solomon Islands or Palau.

1. Introduction

The mitochondrial DNA (mtDNA) genome (mitogenome) is largely maternally inherited and generally consists of 13 protein coding genes (PCG’s), 22 tRNA genes, 2 rRNA genes, and an AT-rich region that is low in nucleotide complexity (e.g., Crozier 1990). Due to the maternal inheritance nature, the mitogenome in general does not undergo recombination (e.g., Saville et al. 1998; Leducq et al. 2017). Hebert et al. (2003) demonstrated the use of the partial mtDNA cytochrome oxidase subunit I (mtCOI) gene sequence to aid in species diagnostics, and this system helped transform understanding of species diversity (Hebert et al. 2004; Rugman-Jones et al. 2013; Tay et al. 2016). Multiple partial mtCOI databases (e.g., BOLD; GenBank) subsequently provided considerable contribution to disentangling species status (e.g., De Barro et al. 2011; Behere et al. 2007; Tay et al. 2017a). However, the use of partial mtCOI is not without its limitations; association with endosymbionts, effect of selective sweep, and impact of pseudogenes, amongst other factors, can all lead to inaccurate interpretations (Hurst and Jiggins, 2005; Tay et al. 2017b). Analysis of population-wide partial mtCOI gene diversity has also found that the 5’ gene region typically has low nucleotide diversity (i.e., conserved nucleotide sequence) in some arthropod groups such as Hemiptera, Lepidoptera, and Coleoptera (e.g., Tay et al. 2022a; Tay et al. 2022b; Behere et al. 2007). Subsequently, an over reliance on partial mtCOI gene sequences has resulted in some misidentification of species status (e.g., Tay et al. 2017b), and some misunderstandings of population dynamics (e.g., Goergen et al. 2016 vs. Tay et al. 2022a; see also review by Tay et al. 2023). Examples of confounded interpretations include the invasion history of Spodoptera frugiperda (Cock et al. 2017; Nagoshi et al. 2018, Nagoshi et al. 2019) and population expansion patterns of Helicoverpa species (Leite et al. 2014 vs. Arnemann et al. 2018).

As technology advances and costs decrease, it is now easier and cheaper than ever before to obtain greater genetic data from specimens to provide richer information content for analysis and interpretation. Combining sequence data from multiple mtDNA genes, full mitochondrial DNA genomes (mitogenomes) and/or whole genomes is now regularly being shown to provide superior analysis to partial mtCOI genes alone for all applications. Examples include: identifications of species (e.g., Walsh et al. 2019), sub-species (e.g., Anderson et al. 2016; Elfekih et al. 2021; Zhang et al. 2022), hybrids (e.g., Anderson et al. 2016; Elfekih et al. 2021; Valencia-Montoya et al. 2020), populations (Tay et al. 2022a; Rane et al. 2023), patterns like demographic expansion (e.g., Pearce et al. 2017; Lannucci et al. 2021), and pest incursion histories (e.g., Tay 2016; Otim et al. 2018; Benelli et al. 2023; Li et al. 2022).

Here, we examine draft full mitogenomes of the coconut rhinoceros beetle (CRB; Oryctes rhinoceros), a pest that causes economic yield losses to coconut and oil palm (Bedford, 1980; Indriyanti et al. 2019). These mitogenomes were generated through whole genome sequencing (WGS) from multiple geographically distinct locations to develop additional molecular markers for tracking and monitoring genetically distinct populations of this species. Particular attention is given to the CRB-G (clade I) group determined using partial mtCOI gene assessment (Marshall et al. 2017), because of the current biosecurity focus on this mitocnondrial haplotype variant due to its reported resistance to known isolates of the Oryctes rhinoceros Nudivirus (OrNV) biological control agent (Marshall et al. 2017), and new incursions within the Pacific region (Paudel et al. 2023; Marshall et al. 2023). Specifically, we test whether or not the Guam CRB-G (clade I) was the source population for the CRB-G (clade I) in Solomon Islands and in Palau as suggested in some publications (e.g., Datt 2020, Caasi 2023). We do this by identifying two partial mitochondrial gene regions that more confidently differentiated CRB (clade I) individuals that invaded Guam from other CRB individuals, including other individuals classed as CRB-G (clade I) using the partial mtCOI gene, but which were collected on other Pacific islands (e.g., Solomon Islands, Palau). We discuss the benefits of mitogenomes as resources for developing alternative diagnostic markers, and assess efficacies of the partial mtCOI gene as the current preferred standard diagnostic DNA marker to distinguish CRB populations.

2. Material and Methods

2.1. Samples

CRB samples collected between 2019 and 2022 were sourced from Guam, Palau, Indonesia, Malaysia, Hawaii, and Philippines (Table 1). The gut of each specimen was dissected, preserved in 100% ethanol, and stored at -20 °C until needed for DNA extraction. Additionally, a historic specimen collected from Guam (04-Or5; circa 2014) was used as a reference to enable matching of more recently collected CRB individuals from Guam that classed within the CRB-G (clade I) haplotype grouping sensu Marshall et al. (2017).

2.2. Whole Genome Sequencing (WGS)

We used the Qiagen Blood and Tissue DNA extraction kit (Duesseldorf, Germany) and the manufacturer’s protocol to extract genomic DNA. All extracted DNA was eluted in 200 µL EB and kept frozen until used for WGS. We assessed the quality of the extracted DNA using Qubit 2.0 prior to sequencing. WGS was carried out by the Australian Genome Resource Facility (AGRF) in Melbourne, Australia, or by AZENTA Life Sciences in China. The WGS data returned an average of 25x coverage, 150 bp paired-end reads/sample, assuming a genome size of approximately 350 Mbp.

2.3. Mitogenome Assembly and Annotation

We assembled all mitogenomes by importing the raw sequence reads into Geneious Prime 2022.2.2 (Build 2022-08-18 14:34) (Biomatters Ltd., Auckland), and used the published mitogenome (MT457815, Filipovic et al. 2021) as the reference sequence. We used Geneious Mapper with ‘Low Sensitivity / Fastest’ option and selecting no fine tuning (i.e., None (fast / read mapping)) during the mitogenome assembling process. Although we received pair-ended reads for all samples, mitogenomes were assembled using forward reads only due to the high genome coverage for each sample. All assembled mitogenomes were initially annotated using the MITOS program and selecting invertebrate mitochondrial genetic code (Bernt et al. 2013). As a final quality assessment, the annotated CRB mitogenomes were visually inspected. The assembled and annotated mtCOI, ATP6, and COIII genes used in this study are available from the CSIRO data repository (Tay et al. 2024).

2.4. Mitogenome Identity Assessment

The non-recombination nature of the mitogenome implies that CRB individuals classified as CRB-G (clade I) based on the partial mtCOI gene assessment method of Marshall et al. (2017) (e.g., Solomon Islands MT457815, Taiwan MW632131) would share mitogenome identity with our reference Guam specimen (i.e., 04-Or5; Table 1), if a single source of invasion entered into Guam and subsequently spread from to other locations. To assess this, randomly selected CRB specimens from Guam (i.e., NZ-20-738; Guam-01_GDoA, Guam-02_GDoA, Guam-09_GDoA, Guam-13_GDoA, Guam-17_GDoA; Table 1) that were collected in more recent times (2020 and 2021) were compared with the representative historical Guam individual (04-Or5) to visually assess and confirm mitogenome identity. This was then followed by comparison with all other CRB individuals including CRB-G type individuals collected from elsewhere (Table 1). Individuals were compared based on the partial mtCOI sequence analysis (described in Marshall et al. 2017) as well as sequence similarity of other mitochondrial genes assessed by this work.

2.5. Alternative CRB Marker Development to Identify the Original CRB Population that Invaded Guam

To identify candidate mitochondrial genes as alternative DNA markers specific to individuals from Guam, all mitogenomes generated from this study, as well as publicly available CRB mitogenomes from GenBank, were downloaded and aligned within GenBank using MAFFT V7.490 (Katoh and Standley 2013; Katoh et al. 2002) with default setting options (i.e., algorithm: FFT-NS-2; Scoring matrix: 200PAM / k = 2; Gap open penalty: 1.53; Offset value: 0.123). We visually identified candidate polymorphic sites unique to individuals from Guam (i.e., 04-Or5, NZ-20-738, Guam-01_GDoA, Guam-02_GDoA, Guam-09_GDoA, Guam-13_GDoA, Guam-17_GDoA), but absent in all other CRB individuals (Table 1). SNPs identified were analysed for potential restriction endonucleases to develop polymerase-chain reaction (PCR) restriction fragment length polymorphism (RFLP) solutions (PCR-RFLP), for a simple and easy-to-use approach to confidently differentiate CRB-G (clade I) that invaded Guam from all other genetically distinct CRB, including CRB from elsewhere classed as CRB-G (clade I) by partial COI assessment. Design and analysis of PCR-primers were through the Primer Analysis Software version 7.60 (Molecular Biology Insights, Inc., Cascade, CO, USA). Primers were optimised for minimal false primer annealing sites, minimal primer dimer and duplex formation, and minimal/no hairpin structure, with a Ta (theoretical annealing temperature) of ≥60 °C (calculated as Ta = 4(G + C) + 2(A + T)), and an optimal amplicon length of between 500-600 bp to facilitate ease of Sanger sequencing. The candidate restriction endonuclease was initially selected for a single cut site with in silico analysis of all different mitochondrial DNA haplotypes in Enzyme X version 3.3 <http://nucleobytes.com/enzymex/>. We reconfirmed primer efficacies and RFLP conditions by randomly selecting and analysing DNA from specimens collected from Guam and elsewhere, as well as by PCR-Sanger sequencing to confirm primer amplification accuracy. We used the restriction digestion conditions as recommended by the manufacturer of the BmpI restriction enzyme (New England BioLabs). Visualisation of the RFLP was on a 1.5% 1x TAE agarose gel.

2.6. Mitogenome Analysis

The mitogenomes from the GenBank database and those generated from this study were aligned to estimate pairwise nucleotide identity and distances (p-dist) between: (i) full mtCOI gene vs. full ATP6 gene, and (ii) full mtCOI vs. full COIII genes. The related O. nasicornis mitogenome (0K484312) was included to provide comparison of inter-specific nucleotide distance with CRB. We also inferred phylogenies of the CRB individuals based on the widely used partial mtCOI gene region (676 bp) versus our proposed alternative mitochondrial ATP6 and COIII partial gene regions (excluding nucleotides at primer annealing sites, see Marshall et al. 2017). The APT6 and COIII partial gene sequences were concatenated before phylogeny inference. We used IQ-Tree (Trifinopoulos et al. 2016) and selected the ‘Auto’ option for estimating substitution models, and 1,000 bootstrap alignments to estimate branch support using the ultrafast bootstrap approximation (UFBoot) (Hoang et al. 2018) algorithm. We used Dendroscope 3 (Huson et al. 2007) for visualisation and post analysis editing for both COI and ATP6+COIII phylogenies.

3. Results

3.1. Mitochondrial Genome Analysis

Mitochondrial genomes were assembled and annotated from an average of 1,472 fragments (mean standard deviation 997 fragments) per sample. Across all the mitochondrial COI, ATP6, and COIII gene sequences, nucleotide differences between CRB individuals were low (<2% difference) suggesting that all were the same species (i.e., O. rhinoceros) (Table S1). The assembled and annotated mitogenome from the Guam 04-Or5 specimen (collected in 2014) provided evidence that all Guam individuals examined here shared the same mitogenome (Table 1). Nucleotide differences within the mitochondrial ATP6 and COIII genes were identified in Guam individuals, but these nucleotide polymorphisms were absent from specimens collected elsewhere, including those classed as CRB-G (clade I) by partial COI assessment. For the ATP6_T4430C and COIII_C5390T SNPs identified from the partial ATP6 for COIII genes (see Table 1), there were, on average, 1,876 and 1,740 reads at each of these nucleotide sites to confirm differentiation of CRB-G (clade I) that invaded Guam from other CRB (i.e., equivalent to the diagnostic SNP for ATP6 and COIII being independently confirmed an average of 1,876 and 1,740 times, respectively).

Pairwise nucleotide analysis of the complete mtCOI gene sequence vs. complete ATP6 gene sequence, and also the complete mtCOI gene sequence vs. the complete COIII gene sequence, showed that the seven Guam CRB, one Solomon Islands CRB (MT457815), one Taiwan (MW632131), two Palau CRBs (Palau-03, Palau-04), and four Philippines CRB (Phil-01. -02, -03, -04) specimens analysed shared 100% identity across the complete mtCOI gene sequence. However, when the comparison included the full ATP6 and full COIII gene sequences, only the Guam individuals remained 100% identical to each other. CRB from Solomon Islands (MT457815), Taiwan (MW632131), Palau (Palau-03, Palau-04), and Philippines (Phil-01. -02, -03, -04) all had polymorphisms in these two alternative mitochondrial marker genes (Table 1).

3.2. Alternative Primers to Identify the Original Invasive CRB Population Present in Guam

Two alternative sets of primers were developed (Table 2) to distinguish CRB-G (clade I) that invaded Guam from other CRB, including those collected elsewhere classed as CRB-G (clade 1) by partial COI assessment. One primer amplifies a partial ATP6 gene region of 494 bp length, and the other amplifies a partial COIII gene region of 469 bp length. Optimal PCR annealing temperature for both ATP6 and COIII was 52˚ C with 1.0 µM primer concentration for both ATP6 and COIII primer pairs, and 0.5 mM dNTPs concentration and 1 unit of DNA polymerase in a 50 µL PCR reaction volume.

3.3. Phylogeny

Phylogenetic analysis was carried out using specimens (see Table 1) with available mitogenome DNA sequence data to allow comparison using the mitochondrial gene regions from COI, COIII, and ATP6. Based on the widely used mtCOI partial gene, Figure 1a three clades could be recognised. One clade (red branches) included seven Guam CRB-G (clade I) specimens, four Philippines (green circles) and two individuals from Palau (yellow circles). The other two clades (blue branches) included six individuals from Indonesia, eight from Malaysia, and two from Palau B. Phylogenetic analysis based on concatenation of both partial ATP6 and partial COIII genes returned a different population clustering pattern (Figure 1b). Together, the use of the ATF6 and COIII gene regions showed that Guam CRB-G (clade I) individuals (red branches) clustered by themselves, whereas the Philippines, Malaysian, and Indonesian individuals clustered largely according to their geographical distributions. CRB from Palau (yellow circles) appeared to have multiple origins, clustering with specimens collected from both the Philippines and Indonesia. However, branch node confidence values for Indonesia (54-78) and Philippines (46-48) were low, suggesting longer sequence lengths from both mitochondrial and inclusion of nuclear genes, as well as more samples, are required for confident assessment. Notably, CRB populations in Malaysia appeared to be consisted of two diverse evolutionary lineages based on both COI and the concatenated ATP6-COIII partial genes, with the unique ON764799 and ON764801 individuals originating from both coconut palm and oil palm hosts from the state of Johor (Anggraini et al. 2023).

4. Discussion

In this study, we characterised and reanalysed the draft mitogenomes of CRB individuals from both the native (i.e., Indonesia, Malaysia, Philippines, Taiwan) and exotic (i.e., Guam, Palau, Solomon Islands) ranges. This is also the first time the mitogenome of all recently collected Guam CRB individuals analysed in this study were found to shared sequence identity with specimens historically collected from Guam by possessing the same mitogenome sequence across all 13 protein coding genes (results not presented), and specifically to the ATP6 and COIII genes that exhibited nucleotide differences with CRB from other locations. This resulted in the ATP6 and COIII protein coding genes being used as alternative DNA markers for differentiating Guam-specific CRB-G (clade I) from the other tested CRB individuals. The remaining individuals from elsewhere, however, including those designated as CRB-G (clade I) (based on the partial mtCOI assessment approach), did not share the same maternal lineage as the Guam CRB-G (clade I) individuals. In other words, the multi-gene assessment (albeit with a limited number of specimens), provided strong supporting evidence that the CRB invasion into Guam was distinct from the CRB invasions detected in from Solomon Islands and in Palau, and therefore Guam was not the source of the CRB that invaded these other locations. Our finding based on full mitogenome and especially the ATP6 and COIII genes vs. the widely used partial COI gene in the CRB studies reported to-date (e.g., Marshall et al. 2017; Reil et al. 2018; Etebari et al. 2020) indicated that the reported mitogenome of CRB-G from Solomon Islands identified from the partial COI gene signature (Filipović et al. 2021) represented a mistaken identity. Increasing sampling of CRB from Guam, Palau, and Solomon Islands is needed to further increase confidence of the specificity of the ATP6 and COIII alternative markers to differentiate CRB-G (clade I) from other CRB populations.

For the PCR-RFLP primers focused on the partial ATP6 gene sequence, separation is based on a BpmI restriction site. Individuals classed as Guam CRB-G (clade I) produced two fragments (i.e., 271 bp and 223 bp), whereas all other CRB remained uncut (i.e., 494 bp) (Table 2). A second primer set was developed based on the COIII gene that can also differentiate between Guam CRB-G (clade I) from other CRB; however, this diagnostic method requires sequence analysis (such as through Sanger sequencing) to detect the presence of a ‘C’ or a ‘T’ base at nucleotide position 5,390 (see Table 1). Although these new markers improve the differentiation between CRB that invaded Guam and other CRB populations, assessment of more CRB individuals from native populations (e.g., Malaysia, Singapore, Sri Lanka, India, Bangladesh, Myanmar, Cambodia, Laos, Vietnam, southern China, Indonesia, Philippines, Taiwan, Thailand), will be needed to provide a more robust confirmation of CRB invasion histories. Also, for all work using molecular diagnostics of CRB, use of either the mtCOI or mtCOIII genes is recommended as an initial approach to first confirm that samples are O. rhinoceros. For example, while the T4430C SNP site within ATP6 from Guam specimens was a T, it was also a ‘T’ in O. narsicornis (see Table 1). Therefore, a direct PCR-RFLP without first confirming species status could lead to misidentification of O. rhinoceros among other Oryctes spp.

The CRB is a hitchhiker pest (Hoffmann et al. in press) and is continuing to disperse to new locations, being recently reported in the Marshall Islands (The Marshall Islands Journal 2023) and multiple Hawaiian islands (HDOA 2024). Notably, our results found that Palau CRB appear to have multiple origins (Figure 1b). The node confidence support estimates in Figure 1b displayed a range of values, with some of the individuals (e.g., from Palau, Indonesia) appearing low (less than 60), which limited the power of inference for better understanding the invasion history of this pest across its distributional ranges. It is likely that future detailed genetic assessments of CRB will provide the resolving power required to further elucidate CRB invasion histories.

Increasingly, WGS and multigene approaches have provided greater analytical power than partial genome assessments, and are therefore rapidly becoming more widely adopted for the interrogation of demographic history and evolutionary relationships of some of the world’s most significant transboundary invasive plant pests (Anderson et al. 2016; Anderson et al. 2018; Zhang et al. 2022; Elfekih et al. 2018; Elfekih et al. 2021; Tay et al. 2022a; Rane et al. 2023), including CRB (Reil et al. 2018; Etebari et al. 2021; Filipović et al. 2021). Given that the WGS/multigene approaches can provide more comprehensive evidence than single gene analyses (e.g., partial mtCOI), and we have found exactly this result with this analysis, we suggest that a detailed study using these more detailed genetic assessments is needed to further improve current understanding of CRB invasion biology.