**4. Discussion**

In this research we aimed to find a reliable pattern of differences related to resistance type using genetic and geometric morphometric analyses. For population structure analysis, we used DArTseq SNP markers. One of the questions we were interested in was whether resistant WCR populations differ at the genetic level. We found no significant evidence of high genetic diversity in any of the assumed populations. However, the estimated values were congruen<sup>t</sup> with moderate genetic diversity across the genotyped beetles. The STRUCTURE revealed three genetic clusters. This classification was also supported by the VanRaden kinship algorithm, where Cry3Bb1\_Cry34/35Ab1 individuals and Cry34/35Ab1 were separated from Cry3Bb1 adapted to crop rotation and non-resistant individuals, although some non-resistant individuals mixed between Cry34/35Ab1, which could be due to the normal evolutionary process. The fact that Cry3Bb1 non-resistant and adapted to crop rotation populations are mixed suggests that they are genetically similar (Figure 4). The neighbor-joining tree separated the individuals adapted to crop rotation, which is to be expected given that the first evolved resistance (not including insecticides) was to crop rotation [12]. Afterwards, all other resistance evolved, and we can see that clearly in this result. The fact that the non-resistant population is not separated could be due to an evolutionary process, as we mentioned earlier.

High-throughput sequencing has provided deeper insight into the molecular mechanisms of resistance [74]. It allowed us to find that many point mutations are found in different genes, suggesting that these mechanisms can occur simultaneously, making it more difficult to understand which one is really responsible for the resistance phenotype [75,76]. In our research, we focused on resistant populations, and we determined that there was some variability between them, but there was no exact pattern. Recent molecular studies show us that different sets of genes are involved in resistance [76–79], which makes it unlikely that universal markers of resistance can be developed to accurately determine the likelihood of a population becoming resistant to a particular compound [75,77,79]. A different number of genes may be involved in resistance, and individuals within a population exhibit different evolutionary patterns of resistance evolution. Therefore, resistance can be found throughout the whole genome, but it is not conditioned by the differences. However, certain shifts could be a warning that some changes in the genome have occurred. Through estimates of genetic diversity, population structuring, and genetic relatedness between individuals, information on the effectiveness of control strategies can be obtained, and recommendations to improve the efficacy of control programs may be possible.

The actual sample size of each site does not need to be large when using SNPs. SNP markers provide the power, not the sample size, as SNPs have genome-wide coverage and there end up being many thousands of SNPs by the time genotyping is complete [80]. The paper by Trask et al. [81] states, "Given that each SNP marker has an individual

evolutionary history, we calculated that the most complete and unbiased representation of genetic diversity present in the individual can be achieved by including at least 10 individuals in the discovery sample set to ensure the discovery of both common and rare polymorphisms." The second paper by Li et al. [82], who also worked with beetles from the order Coleoptera, found that "a minimum sample size of 3–8 individuals is sufficient to dissect the population architecture of the harlequin ladybird, *Harmonia axyridis*, a biological control agen<sup>t</sup> and invasive alien species." They also estimated the optimal sample size for accurately estimating genetic diversity within and between populations of *Harmonia axyridis*. They determined that six individuals are the minimum sample size required.

Wing morphology (size and shape) is the most important trait of an insect's dispersal capacity. For this reason, the integration of different techniques to understand the plasticity and variation of this trait is vital to understanding how they adapt to new environments and to coordinating strategic planning ahead of possible new invasions [3]. Different types of wing morphotypes have been studied to determine the dispersal capabilities of flying insects [83–85]. Le et al. [86] found that narrowed wings in beetles are more efficient for flapping low-level flights. Additionally, for *D. v. virgifera*, wing shape has been identified as a very good trait to measure in different agronomic studies, including studies of life history (sexual dimorphism) and interspecific and intraspecific shape variation [47–49], and wing shape has also been a useful variable when combined with other monitoring tools (genetics (e.g., microsatellites) and traditional traps (e.g., pheromones)) [3].

Mikac et al. [46] showed that beetles adapted to crop rotation had broader wings (cf. susceptible beetle). Mikac et al. [45] expanded the use of differences in hindwing size and shape to examine changes in WCR associated with the development of resistance, specifically to examine potential differences between (*Bt*)-resistant, non-resistant (or susceptible), and adapted to crop rotation populations in the US. In general, the hindwings of non-resistant beetles were significantly more elongated in shape and narrower in width (chord length) compared with beetles resistant to *Bt* maize or crop rotation. This result was confirmed by our study. Mikac et al. (2019) did not separate the *Bt*-resistant populations in their study, but considered them as one population. Therefore, in our study, we separated all *Bt*-resistant populations to see the differences between them. Cry3Bb1\_Cry34/35Ab1 individuals had the broader shape and a more robust wing with an expansion of landmark 14 and a contraction of landmark 9. Cry3Bb1 individuals had the narrower wings, while individuals resistant to Cry34/35Ab1 had similar but smaller wings, distinguished by the expansion of landmarks 3 and 4. The more stable and elongated wing shape was that of the population adapted to crop rotation, in which there was an extension to landmarks 1 and 2 to the left and an elongation to landmark 9 to the right. The non-resistant population is also slightly wider than the population of Cry3Bb1-Cry34/35Ab1, with the movement of landmarks 14 and 2 also slightly to the right and the wider shape that is also produced by the movement of landmark 7 to the upper left. Elongated wings are more aerodynamic and are considered to be involved in migratory movement [46]. Mikac et al. [46] also suggested that this could be a useful invasive dispersal strategy for mated females. In our research, individuals adapted to crop rotation had more stable and elongated wings, suggesting that these individuals could fly long distances. Such differences may impact upon the dispersal or long-distance movement of resistant and non-resistant beetles. Understanding which beetle morphotype is the superior flyer and spreader has implications for managing WCR through integrated resistance strategies. These findings confirmed GM as a reliable technique for resistance detection. In this study, we aimed to confirm the results from SNPs markers with GM. We found that geometric morphometric tools could provide important clues to differentiate resistant and non-resistant populations. One of the principal results was the similarity of the hindwing shape variation between the population after the STRUCTURE analysis, where using both monitoring techniques showed that the more differentiated population was the resistant Cry34/35Ab1.

Here we describe a possibility that combining genetic and geometric morphometrics could be a reliable technique that can be used to reveal differences among WCR populations. Hence, geometric morphometrics can be used as a biomarker for resistance detection as part of a larger integrated resistance managemen<sup>t</sup> strategy for western corn rootworm.

In Croatia, WCR have been investigated in detail (traditional monitoring, genetic monitoring, and GM monitoring), and knowledge about dispersal and adaptive abilities of these invasive insects is well known [3,47,87,88]. Our future work will focus on populations collected in intensive maize-growing areas in Croatia, where WCR populations have become established since their introduction 30 years ago. We will use the comparative techniques presented in this paper to determine whether Croatian populations are potentially resistant and which US WCR population was the source population for Croatia and Europe. This knowledge would help to detect resistant individuals that might invade geographical areas where they are not ye<sup>t</sup> present (e.g., beetles adapted to crop rotation invading Europe where such variants are not present). Such information is very important for biosecurity measures, resistance management, and future control strategies for this pest worldwide.

**Author Contributions:** Conceptualization, M.K.B. and D.L.; data curation, M.K.B., K.M.M., H.A.B. and R.B.; formal analysis, M.K.B., H.A.B. and D.L.; funding acquisition, R.B.; investigation, M.K.B., K.M.M., H.A.B., R.B. and D.L.; methodology, M.K.B., K.M.M., H.A.B., R.B. and D.L.; project administration, R.B.; resources, R.B.; software, M.K.B. and H.A.B.; supervision, K.M.M., R.B. and D.L.; validation, R.B. and D.L.; visualization, M.K.B. and H.A.B.; writing—original draft, M.K.B. and H.A.B.; writing—review and editing, K.M.M., R.B. and D.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This study was supported by the Croatian Science Foundation through the project Monitoring of Insect Pest Resistance: Novel Approach for Detection, and Effective Resistance Management Strategies (MONPERES) (IP-2016-06-7458) and the young researchers' career development project training of new doctoral students (DOK-01-2018).

**Institutional Review Board Statement:** Not applicable. Western corn rootworm is an established pest of maize in the USA and Southern Europe. No special permission was needed for its collection in this study.

**Acknowledgments:** The authors thank Wade French and Chad Nielson for Bt-maize-resistant rootworm colonies from South Dakota and Joseph Spencer for providing field-collected beetles adapted to crop rotation from Illinois. The authors are very grateful to Reza Talebi and João Paulo Gomes Viana for their help and advice with data analysis. The authors also thank colleague Zrinka Drmi´c and student Patricija Majceni´c for help in preparing the WCR wings.

**Conflicts of Interest:** The authors declare no conflict of interest.
