Morphological and Molecular Identification of Bark and Ambrosia Beetles and Their Associated Predators in Honduran Coniferous Forests

1. Introduction

In the last decades, the faunistic composition of forest ecosystems has experienced significant changes because of bark beetle attacks [1]. It has been estimated that bark beetle pests affect 30 million hectares of forest annually [2]. Consequently, they represent a substantial burden in terms of both economics and ecology due to their negative impact on conifer forests worldwide. Bark beetles (Coleoptera: Curculionidae, Scolytinae) comprise a large and diverse group of insects with more than 6,000 species [3,4] distributed in all continents except Antarctica [3]. Although few genera are recognized as important disturbance agents in conifer forests, invasive species represent a major threat to native biodiversity [5,6], which may lead to tree death and ultimately affect the productivity dynamic of forest ecosystems [7]. Further, several biological and mobility drivers can favor the risk of introduction and establishment of invasive species, increasing their geographic expansion [8]. In addition, warmer temperatures and long-lasting droughts have been linked to tree stress, but also bark beetle survival and population growth [7,9], contributing to changes in tree physiological processes (e.g., leaf water potential, water content, osmotic potential) and greater coniferous forest infestation [10].

Several bark beetle species have been the subject of research in Central America, particularly those in the genus Dendroctonus due to their aggressive behavior in invading and destroying pine forests [11,12,13]. Honduras, located in the geographic center of the Central American isthmus, is characterized by mountainous regions, from which coniferous forests cover approximately 1,951,977 ha, representing 30.91% of forest ecosystems [14,15]. These forests consist primarily of Pinus oocarpa and P. caribaea [16], playing a key role in many environmental and economic services, including habitat provision, genetic resources, water filtration, climate regulation, carbon sequestration, and erosion control, among others. However, some studies have also reported the detrimental effects of beetle outbreaks in Honduran coniferous forests [16,17], leading to incalculable ecological and economic impacts. For instance, the Honduras Forest Institute of Conservation (ICF) reported that coniferous-dominated forests were extremely affected by Dendroctonus frontalis between 2014-2016, a period characterized by prolonged droughts and low precipitation. Moreover, the continuous emergence of Ips beetle attacks has triggered ecological succession in Honduran pine forests, with potential cascading effects on forest ecosystems [18]. Given that bark beetle assemblages have profound ecological effects, it is critical to assess their taxonomic and spatial distribution to adequately mitigate beetle propagation, overall forest health, and resiliency.

In North and Central America, Wood’s monograph [19] has been widely used for insect identification based on morphological features. However, species identification based on dichotomous keys often does not allow a precise description at a lower taxonomic level due to a lack of information for new genera or species. Additionally, DNA barcoding analyses have been largely implemented to validate morphological identification in insect taxonomic inventories [20,21,22]. As a result, integrating morphological and molecular methods enables accurate identification in taxonomic surveys and allows the assessment of population structure. Despite its broad relevance in taxonomic studies, the Honduran insect fauna has rarely been investigated applying molecular techniques, highlighting some previous reports in thrips [23], flies [24,25], and mosquitoes [26]. Moreover, bark and ambrosia beetles, as well as their biological predators have never been evaluated using DNA barcoding analysis in Honduras. Consequently, our analyses comprise the first attempt to apply DNA barcoding to identify insect species associated with pine forests. Here, we used both morphological and molecular methods to evaluate the Scolytinae community composition across different regions in Honduras, as well as other families of beetles. Additionally, due to secondary bark beetles of the genus Ips cause significant damage in coniferous forests [27] and the actual status distribution is poorly known in our country, we also determined the spatial distribution and genetic diversity of Ips apache beetle assemblages by analyzing intra-specific variation in the mitochondrial gene COI (cytochrome C oxidase subunit I).

2. Materials and Methods

2.1. Study Area and Collection of the Specimens

During 2019-2023 a total of 1,131 insects of the order Coleoptera were collected from 12 departments of Honduras (e.g., Comayagua, Copán, Cortés, El Paraíso, Francisco Morazán, Intibucá, La Paz, Lempira, Olancho, Santa Bárbara, Yoro, and Atlántida) (Figure 1). The data on the individuals and localities are shown in the Supplementary Table S1. The specimens were trapped with Lindgren multiple funnel traps, placed approximately 0.5-1.0 m above the ground, and baited with aggregation pheromone frontalin and tree-emitted compound α-pinene. Each trap was filled with coolant to preserve captured insects. In most cases, traps were inspected every two weeks to collect the specimens and then transported to the “Laboratorio de Diagnóstico Sanitario Forestal” at the ICF, Tegucigalpa, Honduras. Specimens that were initially screened as beetles were selected, and those specimens identified as Dendroctonus spp. were separated from the rest of the individuals for a separate analysis whose results are not shown in this study. Consequently, this study only evaluated insect fauna considered secondary bark and ambrosia beetles, as well as natural predators. The selected specimens were transferred to glass vials (16 mm diameter × 50 mm height) containing 70% ethanol and transported to the Genetics Research Center at the National Autonomous University of Honduras (UNAH).

To magnify morphological features, specimens were observed under stereomicroscopes (Zeiss Stemi DV4 and UNITRON Z10 series) and identified using dichotomous keys [3,19,28,29]. Identification accuracy was assisted using bark and ambrosia beetles’ websites (https://www.barkbeetles.info/, https://ambrosiasymbiosis.org/), and validated by experts.

The species were classified according to the feeding habits of the immature stages in the following categories: (a) Myelophages or with feeding in the medulla; (b) phyllophagous to those that feed in the inner cortex; (c) feeding of the roots; (d) spermatophagous or seed feeding; (e) xylomycetophagous or fungal growers, and (f) xylophages to those that feed on xylem tissues. In a simplified way, all the genera that feed on phloem tissues were considered as bark beetles, and those that cultivate fungi were considered as ambrosial. Those that feed on tissues outside the trunk (such as seeds or roots) or internal tissues such as the pith were not considered.

2.2. Alpha Diversity Analysis

Descriptive statistics were used to calculate the number of families, subfamilies, genera, and species collected per department. An assessment was conducted to determine the departments with the highest diversity. For this purpose, the mathematical framework of Hill numbers was used. Hill numbers differ by the parameter q, which determines the sensitivity of the measure to the relative abundance [30]. Using the Hilldiv2 R package v.2.0.2 [31], we computed the orders q = 0 which corresponds to the species richness, and q = 1 which considers species according to their relative abundances and is equivalent to the Shannon’s entropy index.

2.3. DNA Extraction and Barcoding

After morphological identification, individuals were stored separately in 1.5 mL tubes with 70% ethanol. Then, a small portion of soft tissue was retrieved from the abdomen to perform DNA barcoding analysis. A subset of 351 samples (31%) were molecularly evaluated. Those genera with less than three individuals were not sequenced since specimens were fully preserved to corroborate their identification and as vouchers. Detailed information about the sequenced samples can be found in the Supplementary Table S2.

Genomic DNA was isolated using the Extracta DNA Prep for PCR kit^® (QuantaBio, Beverly, MA, USA) following the manufacturer’s instructions. Briefly, soft tissues were macerated and incubated at 95°C in 100 μL of extraction reagent for 25 minutes and then cooled to room temperature for 5 minutes. Then, 100 μL of stabilization buffer was added and the final volume (200 μL) was transferred to a new 1.5 mL vial. A fragment of the mitochondrial gene COI was amplified using one of three primer pairs (Table 1).

The PCR reactions were conducted in a 50 μL reaction mixture, containing 25 μL of KOD One^TM PCR master Mix (Toyobo Co, Ltd. Tokyo, Japan), 2 µL of each primer (100 µM), 2 µL of acetylated albumin - BSA (10 mg/mL), 15 μL of nuclease-free water, and 4 μL of DNA template. PCR reactions were performed under the following conditions: initial denaturation at 98°C for 30 s, followed by 37 cycles of 10 s at 98°C, annealing at 50°C for 10 s, elongation at 68°C for 2 min, and a final extension at 68°C for 2 min. For each set of reactions, a negative control was included. The PCR efficiency was then visualized through a 1% agarose gel electrophoresis with ethidium bromide. The resulting PCR products were purified and sequenced by Psomagen^® (Rockville, MD, USA) (https://www.psomagen.com).

2.4. Bioinformatic Analyses

The raw sequences were edited and assembled using Geneious^® prime software (Dotmatics, Boston, MA) [36]. The consensus sequences were deposited in GenBank, and accession numbers were assigned. Each sequence was subjected to a BLAST analysis on the NCBI platform (https://blast.ncbi.nlm.nih.gov/Blast.cgi) to validate the taxonomic identification made by the experts. If the BLAST result revealed a percentage of identity of less than 95% with respect to the result with the greatest similarity, morphology-based identification was employed to identify the deposited sequences. In addition, an open-access project titled “Project - COLEH Biodiversity of Bark beetles and other Coleoptera in Honduras” was created in the Barcode of Life Data Systems database (http://www.boldsystems.org).

To analyze the intraspecific genetic diversity of I. apache, 52 assembled sequences were aligned, and the resulting alignments were utilized to create a phylogenetic tree using the Tamura-Nei genetic distance model and the Neighbor-Joining tree construction method and performing 1000 iterations of Bootstrap. A homologous sequence of Enoclerus sp. was included as an outgroup in the analysis. The collecting site was considered to ascertain the presence of a population structure for I. apache.

The number of polymorphisms among the I. apache sequences, the number of segregating sites (S), the average nucleotide differences (k), the number of haplotypes (H), the haplotype diversity (Hd), and the nucleotide diversity (π) were calculated using DnaSP Software (Version 6.12.03). Additionally, Tajima’s D test [37] was performed in DnaSP to evaluate the neutrality theory of evolution. Haplotype networks were generated using the Median Joining algorithm in Network (Version 10.2).

3. Results

3.1. Morphological Identification

This study evaluated the beetle fauna associated with coniferous woods in Honduras. Using morphological traits, we identified 27 genera of Coleoptera collected in 12 departments of Honduras. The community composition was grouped into four families and eight subfamilies (Table 2). Our analyses described three ecological groups: bark and ambrosia beetles, and natural predators. Most genera were classified as bark and ambrosia beetles and only about 10% as predators (Table 2). A photographic record of the genera described in this study is shown in Supplementary file S1.

In terms of relative abundance, our results revealed that Ips (54%), Xyleborus (16.5%), Temnoscheila (10%), Tomolips (5%), Euplatypus (3%), Hypothenemus (2%) and Gnathotrichus (2%) comprised more than 93% of the total collected individuals, whereas the other 20 genera accounted 7% of the collection. In addition, according to the number of elytral declivity spines, two species of Ips were identified: I. apache and I. cribricollis, most of them collected in Francisco Morazán, Yoro, and Olancho.

The largest number of specimens were collected in El Paraíso (n = 209; 18.5%), followed by Olancho (n = 186; 16.4%), Yoro (n = 178; 15.7%), Copán (n = 167; 14.8%), Francisco Morazán (n = 121; 10.7%). The remaining departments contributed less than 7% of the total (Figure 2). To our knowledge, this is the first report of the presence of Xylomeira and Stephanopachys in Honduran pine forests.

3.2. Diversity Analysis

Alpha diversity was computed using the Hill numbers by parameterizing the q value. For instance, diversity analysis at order q = 0 revealed that four departments have the greatest richness of genera: Olancho and El Paraíso have effective numbers of 17, Copán = 10, and Francisco Morazán = 9; whereas taxonomic diversity at order q = 1 increased in El Paraíso (n = 6.5), Olancho (n = 5.8), Santa Bárbara (n = 5), La Paz (n = 4.5) and Copán (n = 4.4) (Figure 3). While these departments displayed greater diversity, the low diversity for the remaining sites could be related to the small sample sizes and therefore affected the statistical power of our analysis.

3.3. Barcoding

DNA was isolated from 376 individuals, and successful amplification was achieved for 212 (56.4%) using at least one of the PCR procedures, from which 124 sequences (55.5%) yielded a high-quality score. A total of 56 sequences were identified as I. apache (44.8%), 21 as I. cribricollis (16.8%), 12 as Xyleborus spp. (10.4%), 6 as Euplatypus spp. (4.8%), and 5 as Temnoscheila spp. (4.0%). Additionally, five sequences were from Gnathotrichus spp. and Hypothenemus spp., and three sequences were from Cryptocarenus lepidus. Moreover, four sequences were from Enoclerus spp., two sequences each from Pityophtorus spp. and Tomolips spp., and one sequence each from Hylastes spp. and Xyleborinus spp. GenBank accession numbers and BOLD identification numbers (BINs) are listed in Supplementary Table S2. Most of the sequences obtained in this study are the first reported in GenBank for the identified species. Furthermore, one of the sequences was identified as a mite of the taxon Trichouropoda, with 85% percent identity in NCBI. Three sequences were identified as the intracellular bacterium Rickettsia belli with 97% identity. Moreover, six sequences were identified as the nematode Deladenus spp. with identity percentages between 86-89%.

3.4. Intraspecific Diversity of Ips apache

Regarding the analysis of intraspecific diversity of I. apache, it was observed that among 51 sequences with a length of up to 630 bp, the number of identical sites was 568 (81.5%), and the Pairwise % identity coefficient was 98.6 %. A population structure related to the collection site was not demonstrated for I. apache (Figure 4).

To analyze the number of haplotypes in the I. apache population, 49 sequences of 527 bp were aligned. The number of segregating sites (S) was 93. The nucleotide diversity (π) was equal to 0.0142, and the average nucleotide difference (k) was 7.156. Twenty-three haplotypes (H) were found (Figure 5), with a diversity index (Hd) of 0.9073 (SD = 0.027). Haplotypes 5 and 1 were the most frequent, with 12 and 8 sequences respectively. 18 of 23 (78.3%) haplotypes included a single sequence. Haplotype 11 was formed with 3 sequences, while haplotypes 2 and 12 were formed with 4 sequences each. The haplotype networks were constructed considering the department of origin of the specimen (Figure 5). No type of correlation was observed between the haplotypes and the collection site.

To explore the influence of natural selection on population structure, Tajima’s D test was performed on the COI sequences and calculated a D value of –2.45039 (Statistical significance: **, p < 0.02). Fu and Li’s D* and F* tests resulted in values of –5.47644 and –5.19970 respectively (Statistical significance: **, p < 0.02), suggesting an excess of rare variation, consistent with population growth, or positive selection.

4. Discussion

Few studies about the biodiversity of beetle fauna related to coniferous forests in Honduras have been done, limiting our knowledge of their ecological and functional roles. Therefore, a better understanding of insect inventories and spatial distribution is required to increase awareness of beetle outbreaks and their negative impacts on ecosystem services at regional scales. In this study, by incorporating morphological and molecular methods, we examined the community composition of bark and ambrosia beetles, and their natural predators in pine forests. Overall, our integrative analysis provided a more complete and suitable assessment of the forest insect community at lower taxonomic levels (e.g., genus and species). In line with our findings, recent studies have also evidenced the suitability of integrated morphological and molecular approaches to taxonomic analysis in insect communities [20,21,22,43].

4.1. Beetle Diversity and Spatial Distribution

Our collection revealed a high diversity of beetles, composed of 27 genera, with the majority found in the subfamily Scolytinae (60%). We observed that Ips, Temnoscheila, and Xyleborus displayed the broadest range of distribution in Honduras. Recent studies have suggested that these beetles play a crucial role in forest ecosystem dynamics, enabling ecological balance, natural renewal process, forest structure, and succession [44,45]. In addition, the large-scale geographical distribution and colonization of these genera can be related to their high dispersal capacity over short distances within a forest stand or long distances above the forest canopy [46,47]. Particularly, it has been argued that Ips spp. can fly above the canopy and move up to 55 km for several hours [47], and thus can colonize a wide variety of forest ecosystems in Honduras. The extensive distribution of Ips spp. (secondary predator), in addition to the previously reported prevalence of the southern pine beetle Dendroctonus frontalis (primary predator) [14,16] is particularly noteworthy due to the potential resource conflict that may occur between them [48]. However, further studies are needed to properly address their interspecific competition at temporal and spatial scales.

Moreover, according to the statistical analysis at the genus level, Olancho and El Paraíso harbored a higher alpha diversity. Strikingly, in the last decades, El Paraíso has experienced the greatest attacks by bark beetles (Dendroctonus spp. and Ips spp.), losing at least 265 ha of healthy forest [15]. The long-term consequences of bark beetle attacks have resulted in extensive tree mortality, alteration of forest structure, and low productivity. Despite this continuous and persistent infestation, the high degree of ecological integrity and resilience of pine forests seems to maintain a great diversity of beetles.

4.2. Predator-Prey Interactions and Ecological Relationships

Our findings revealed that the secondary bark beetle species of the genus Ips and their natural enemies, i.e., Temnoscheila (Trogossitidae) and Enoclerus (Cleridae), occur simultaneously in most sampled sites, indicating an overlapping distribution. Both genera have been extensively identified as associated predators of harmful forest pests such as Ips and Dendroctonus species [49]. While Temnoscheila and Enoclerus showed a greater diversity and broad distribution [50,51], their diversity, natural history traits, population dynamics, and interaction with primary/secondary bark beetles have been little explored in Honduras. Additionally, the results showed a positive response of bark beetles and their associated predators to Lindgren funnel traps baited with frontalin and α-pinene. The use of targeted semiochemicals in Scolytinae trapping studies provides a unique opportunity to assess their associated insects and improve biological control strategies. For instance, Aukema & Raffa [52] reported that Thanasimus dubius (Cleridae) was strongly attracted to frontalin and α-pinene pheromones. However, further research on selective attraction as well as testing new combinations of semiochemicals will improve forest management strategies during bark beetle infestation.

4.3. Genetic Diversity and Phylogenetic Analysis

DNA barcoding has become a powerful tool for estimating intraspecific genetic diversity in Coleoptera [53,54,55,56]. However, the effectiveness of DNA barcoding depends on comprehensive reference libraries. If the reference database lacks sequences for certain species, accurate identification becomes challenging. For most of the genera/species described in this study, there were no COI sequences available in the databases. As a result, it was not possible to compare these sequences with those acquired from specimens collected in other countries. Therefore, this study presents, for the first time, partial sequences of the COI gene from many different species and genera of beetles that reside in coniferous forests in Honduras.

On the other hand, the phylogenetic analysis and the haplotype network show intraspecific genetic differences among individuals of I. apache, however the sequences clustered without any clear structure. Moreover, a π value of 0.0142 suggests a moderate level of genetic diversity whilst a haplotype diversity (Hd) of 0.9073, indicates a high level of genetic variability. These results suggest a substantial level of genetic diversity, which is beneficial for the population’s adaptability and resilience to environmental changes. This lack of genetic differentiation among populations coupled with moderate to high levels of genetic diversity, may suggest gene flow among individuals inhabiting different coniferous forests in Honduras. In addition, the predominance of 2 haplotypes implies that certain haplotypes are more common within the population, which could be due to selective advantages or historical demographic events.

Furthermore, neutrality tests suggest an excess of rare variation, consistent with population growth, or positive selection, which implies that I. apache populations are potentially well adapted to environmental conditions in Honduran coniferous forests, increasing their persistence and spread within these habitats. Ips apache is a well-known species for its host specialization and feeding primarily on pine forests [57]. Hence, a level of substantial genetic diversity and allelic richness, coupled with limited population differentiation, might be attributed in part to the lack of physical barriers and the very favorable environmental conditions that promote population expansion. In addition, given that Ips bark beetles are a priority issue for local governments, predicting future trajectories will provide a solid scientific basis to estimate their impacts on forest ecosystems. However, more genetic and ecological studies are needed to assess the genetic divergence and environmental adaptations of Ips engraver beetles in Honduran coniferous forests

5. Conclusions

As far as we know, this research is one of the most comprehensive investigations conducted on the diversity of bark beetles in Honduras. The methodology employed in this research combines morphological and molecular techniques to offer a comprehensive understanding of the diversity and composition of beetle communities that reside in coniferous forests. Furthermore, the mitochondrial COI gene reveal significant genetic diversity among I. apache populations but no structure, suggesting gene flow among individuals from different localities. Additionally, due to the greater impact of bark beetles in Honduran pine forests, our findings are important for pest management strategies and phytosanitary measures to prevent the risk of introduction and spread of harmful beetle species. Lastly, future studies should focus on understanding the trophic relationships between bark beetles and natural predators, providing valuable insights into fundamental ecological processes in coniferous forests.