Geo-statistical modelling of the spatial distribution of acorn weevil, Curculio glandium Marsham (Col., Curculionidae) in the oak forest (Ilam province, Iran)

Introduction

Forests are one of the most important renewable resources that meet many human needs. The level and quality of these valuable resources have been declining in recent years due to population growth in human societies and the lack of implementation of scientific, comprehensive management methods (Marvi-Mohadjer, 2005). The forests of western and northwestern Iran are the primary habitat for oak trees and account for about 41% of the country's total forest area. In recent decades, the forests of the Zagros have faced numerous risks, including drought, fire, forest degradation, tree cutting, excessive livestock grazing, and, especially, plant pests and diseases. All of these factors have ultimately led to reductions in forests, vegetation density in natural areas, and the reproduction of important plant species (Hamzehpour et al., 2011). Several factors contribute to the destruction of forests, especially in the Zagros forests, which have been subject to numerous disturbances in recent decades. Most tree species in the Zagros forests belong to the genus Quercus L. and are considered important forest resources in Iran (Fattahi, 2002). The genus Quercus is widely distributed worldwide and is commonly found in most forest ecosystems (Jimenez-Pino et al., 2011).

The dominant species of these forests is the oak Quercus brantii Lindl, and the reproduction of trees is done both sexually (seeds) and asexually. Reproduction by seeds is the natural method of reproduction in all the forests of the world (Marvi-Mohadjer, 2005). The presence of seeds is essential for the establishment of sexual reproduction in the forest and is considered an important food source for wildlife species (Goodrum et al., 1971; Auchmoody et al., 1993; Greenberg, 2000). Due to the abundance of nutrients, beetles, moths, and other insects attack oak seeds both before and after they fall (Leiva & Alés, 2005). Unfortunately, in the forests of Zagros, sexual reproduction is associated with low success because seeds are subjected to heavy insect attack (Zargaran et al., 2018). Insects that recently damage and infect the seeds of oak trees include two pests, the acorn weevil and the acorn moth, belonging to the orders of Coleoptera and Lepidoptera.

Meanwhile, the acorn weevil (Curculio glandium Marsham) belongs to the Curculionidae family, which is a serious pest of these trees. Natural reproduction is considered one of the most important factors affecting the survival and sustainability of natural forests. Many factors affect the success of reproduction that must be known and properly controlled. Due to their high nutrient content, oak seeds are heavily attacked by beetles, lepidopterans, and other insects before and after they fall (Leiva & Alés, 2005). Unfortunately, sexual reproduction is less successful in the Zagros forests. Because the seeds of these trees are highly susceptible to attack by beetles, lepidopterans, and other insects before and after they fall, due to their high nutrient content (Zargaran et al., 2018).

Curculio glandium is an important forest pest of oak trees in many countries, including North America and Europe (Udaka & Sinclair, 2014) and Iran (Mirzaei & Mirab-balou, 2015). The pest disturbs the regeneration of host trees by feeding on their acorn. Adult females lay their eggs in developing acorns in May, and the larvae cause damage (Venner et al., 2011). The larvae are short and cylindrical in shape and move by means of ridges on the underside of the body. Until autumn, it lives inside the acorn; later on, it leaves the acorn through an exit round hole and overwinters in forest litter and 10 to 30 cm soil depth. C. glandium can pose a large risk to acorn-bearing tree populations. They are highly effective at infecting acorns, potentially rendering 70–90% of seeds non-viable (Oltean & Stana, 2007). Geostatistics is a branch of statistical science in which the coordinates of the data related to the society under investigation and the spatial structure of the relevant data are studied. Geostatistics studies variables that exhibit spatial structure; in other words, the presence or absence of spatial structure in the data must first be examined, and, if present, the data must be analyzed. Variables that have this kind of spatial structure are called regional variables. From the perspective of geostatistics, each sample has a spatial relationship with its surrounding samples up to a certain distance, in which case the possibility of similarity between the values ​​of closer samples is greater. In geostatistics, attempts are made to improve the accuracy of estimates and to determine the appropriate time and place for measuring and collecting information. Geostatistical estimates are called Kriging and include various methods, each with its own capabilities and limitations (Hassani Pak, 1997; Madani, 1994).

The use of geostatistics as an important method in sampling arthropods has expanded significantly. Using this method can reduce the number of required samples to the minimum possible (Liebhold et al., 1993). Geostatistical estimation consists of two stages: the first is the identification and modeling of the spatial structure, which is examined using a variogram (semivariogram). The second stage is the geostatistical estimation by the Kriging method, which depends on the characteristics of the variogram model fitted in the first stage (Hassani Pak, 1997). A semivariogram is a basic tool in geostatistics; it is a function of the distance between samples and quantifies the spatial correlation between measured points as a function of distance and direction. The kriging method is also used to estimate the value of the regional variable at unsampled points in a field or forest (Madani, 1994). In geostatistics, the greater the spatial dependence between samples taken from an insect population, the more accurate the population estimate at unsampled points using the kriging method (Isaaks & Srivastava, 1989). Kriging is a powerful method of spatial interpolation that uses complex mathematical formulas to estimate unknown values based on values at known points (Scheeres, 2016). In general, geo-statistical estimation is a process during which one can estimate the value of a quantity at points with known coordinates using the value of the same quantity at other points with known coordinates and determine the density of the pest population in different parts of the field or forest (Liebhold et al., 1993).

Ilam province, with an area of ​​20,027 km², has 641,667 ha of forests, which are considered dry and semi-arid forests of the Zagros Mountains. Shao et al. (2021) investigated the spatial distribution of herbivorous insects in a forest landscape containing oak trees of the species Q. variabilis Blume. and Q. aliena Blume. using geostatistics and showed that the best-fit model was exponential-linear. Also, the values ​​of C₀/(C₀+C) ranged from 0.389 to 0.685 (Shao et al., 2021). Lasmar et al. (2012) used geostatistics to determine the spatial distribution and infestation of leaf-cutting ants (Hymenoptera: Formicidae) nests in eucalyptus plantations. They concluded that the spatial distribution of this pest is aggregate (Lasmar et al., 2012). Hohn et al. (1993) investigated a geostatistical model to predict the spatial dynamics of defoliation caused by the gypsy moth (Lepidoptera: Lymantriidae) Lymantria dispar. The results showed that this approach performs relatively well in predicting the spatial pattern of damage. Overall, this model can be useful as a support tool for monitoring, management, and pest control planning in forests, but further research is needed to improve predictions of severity and the timing of spread for operational use. (Hohn et al., 1993). Hlásny et al. (2009) studied the geostatistical simulation of bark beetle infestation for forest conservation purposes in Slovakia. They used variogram modeling techniques, turning bands simulation, and selected post-processing simulation techniques. They produced maps that show the probability of exceeding the critical volume of timber harvest at a given location. They suggested that this system could be used effectively to improve pest control measures and develop strategies to reduce the rate of spread (Hlásny et al., 2009). Jimenez-Pino et al. (2011) investigated the spatial distribution of the oak seed borer Cydia fagiglandana (Zeller) in holm oak (Q. ilex L.) forests in Spain. Their results showed that the spatial distribution of this pest, as determined using geostatistics, is cumulative. They also concluded that interspecific competition with the weevil C. elephas was the most important factor affecting the spatial distribution of the oak weevil over four years (Jimenez-Pino et al., 2011). The spatial distribution of C. glandium is expected to exhibit significant aggregation rather than a random or uniform pattern, driven by intrinsic biological factors. Such aggregation typically arises from patchy host-plant availability and the species' dispersal limitations, which confine individuals to specific microhabitats. Consequently, employing spatial statistical tools is essential to move beyond simple density measurements and to rigorously quantify whether these observed patterns represent true biological clustering or merely stochastic spatial noise.

By using information on the spatial distribution of the pest population in each region, pest management can be carried out in the shortest time and at a lower cost than with other sampling methods. In addition, by preparing a Kriging map, we can obtain a map of the pest's distribution in fields/Garden or forests. Therefore, the purpose of this study was to analyze the spatial distribution of the acorn weevil in the Dalab forest area of Ilam province (western Iran) using geostatistics over two years, to obtain information useful for implementing a sustainable management strategy for this forest pest. Despite extensive studies of forest insect fauna in Iran and international research on the ecology of oak weevils, our knowledge of the spatial dynamics of C. glandium in Iranian forest ecosystems remains extremely limited. Most domestic studies have been limited to checklists or simple biological surveys. So far, no geostatistical model has been presented to explain the patchy distribution of this pest in the country's oak stands. This study represents the first geostatistical characterization of the spatial distribution of Curculio glandium in Zagros oak forests, providing kriging-based prediction maps that can be used as a basis for integrated pest management (IPM) to reduce pesticide use and focus on pollution hotspots in oak forests.

Materials and methods

Study site

This study was conducted in the Dalab Strait (N 33° 41', E 46° 25'), which is located in Chavar city, Ilam province, west of Iran (Fig. 1). The area is about 850 ha, of which 370 ha was studied for this study according to the conditions of the region. The climate of the studied area is semi-arid and cold; the average annual rainfall is 525 mm, and the average annual temperature is 16.9^ºC. The altitude range is 1100 to 2000 m above sea level. In terms of vegetation, the area is mostly forested, with Iranian oak (Q. brantii Lindl.) as the dominant species.

Sampling method

The present study was conducted in two seasons, autumn and winter, from October to March in the years 2022 and 2023. Sampling was done monthly. The area of this strait is about 850 hectares; for this study, due to regional conditions, nearly 370 ha were examined. In this research, a systematic sampling pattern was used to uniformly cover the entire area. so that the whole area is divided into six altitude classes (1300-1360 - Oct, 1360-1420- Nov, 1420-1480- Dec, 1480-1540- Jan, 1540-1600 – Feb, and 1600-1660 – Mar) meters above sea level, and in each class an area of ​​about 60 ha was determined and 30 trees were examined. The 60-hectare area was divided into a grid of cells. 30 cells were randomly selected and an oak tree sample was taken from the center of each selected cell. In each sampling, four oak seeds were randomly collected from four geographic directions and transferred to the laboratory of the Department of Plant Protection at Ilam University. In the laboratory, the number of seeds infected with the pest was counted. Each of the four acorns in the four main geographical directions in the shade of each oak tree was determined as a sampling unit.

Fig. 1. The location of the Dalab Strait (Chavar County, Ilam Province, Iran), showing forest at the top

At each visit, oak fruits were collected directly by the first author. In the laboratory, the number of seeds infected with the acorn weevil was counted. Oak fruit borer damage can usually be identified by the shape of the exit hole created by its larvae. This means that the weevil larvae open their exit hole in a circular shape, while the butterfly larvae's exit hole is oval and smaller than the weevil's exit hole. In addition, the butterfly larvae's excrement is in the form of larger grains, while the weevil larvae's excrement is fine-grained or powdery. All four oak fruits in the four main geographical directions (North, South, East, and West) in the shade of each oak tree were determined as a sampling unit. During each visit, 120 oak fruits were collected from 30 oak trees. Each oak tree represents a single spatial observation point, and the four acorns collected from each tree are combined into a single quantity.

Spatial distribution analysis

The spatial distribution of acorn weevil larvae, based on the number of infested fruits, was studied using classical statistical methods and by developing population density maps. In general, geostatistics has two main parts: variogram and kriging (Cressie, 1993). The semivariogram is used to describe the spatial relationship between values of a regional variable (for example, the population density of a pest) across different parts of the study area. As a result, it is a basic tool in geostatistics. The Semivariogram shows the spatial correlation between measured points as a function of distance and direction. The relationship of the semivariogram (𝛾) is as follows (Madani, 1994).

In which N(h) is the even number of samples that are located at a horizontal distance of h from each other, z(x_i) and z(x_i+h) are the values ​​of the regional variable at two points at the coordinates (x_i) and (x_i+h). Each of the calculated values ​​of the variable or semivariogram, together with their equivalent value h, specifies a point in the coordinates h (Madani, 1994). The formula for each of the models is as follows:

In which γ(h) is the semivariogram for the distance h, C₀ is the segment effect, C is the structural variance, and R is the range. Heroot-mean-square error statistics were used to evaluate the performance of different interpolation models. In variogram model fitting, the residual sum of squares (RSS) is defined as the sum of the squared differences between the experimental se][\-mivariogram and the values predicted by the] fitted model. It can be expressed as:

Where  𝛾(hi) denotes the experimental semivariogram at lag , is the fitted model value at the same lag, and nnn is the number of lag classes. A smaller RSS indicates a better fit between the model and the observed variogram data (Cressie, 1993).

Kriging interpolation formula

The kriging method was used for spatial interpolation and spatial map preparation. Kriging is the most common geostatistical estimation method, widely used for its ability to minimize error variance and provide unbiased estimates (Polhmann, 1993). In this study, the ordinary kriging method was used. Kriging is a weighted moving average, and in fact, it is a technique of interpolating a variable in unsampled areas with the help of variable values ​​in adjacent points and the weights determined by the change model (λi), and its formula is as follows (Webster & Oliver, 2000):

 In which  Z_x0^*  is the estimated value at the point  and Z(_xi) is the actual value of the sample at points x_i. N the number of observations located in the neighborhood of the point we want to estimate, and λ_i is the statistical weight of the sample Z(_xi) it is found in the neighborhood of the characteristic point . To evaluate the effectiveness of interpolation methods, the mean squared error of estimation (RMSE) was used to assess the models (Mohammadi, 2006).

In which Z^*(_xi) is the estimated value at point x_i, Z(_xi) is the measured value at point x_i, and n is the number of points. Theoretically, when this criterion becomes zero, it means that the estimated value is exactly equal to the actual quantity. Cross-validation was used to evaluate the efficiency of interpolation methods. The Root Mean Square Error (RMSE) was used to assess the models. Spatial variability was also examined using semivariograms. The spatial dependence was classified based on the nugget-to-total semivariance ratio (expressed as a percentage), where a ratio of less than 25% indicated strong spatial dependence, between 25% and 75% indicated moderate dependence, and more than 75% indicated weak spatial dependence (Cambardella et al., 1994).

Data analysis

Data collected from sampling in both years of study were entered into Excel version 2010 software (Microsoft 2010). The normality of the data was assessed using SPSS version 23 (SPSS 2015) and the Kolmogorov-Smirnov test. If needed, they were normalized using the data transformation method (Box-Cox Transformation). Finally, the data were analyzed using GS⁺ version 5.1 software (Robertson, 2008), and variogram curves and population kriging maps were prepared.

Results

The results of this study showed that in both years of study, the infestation rate of the acorn weevil was higher than that of the acorn moth (Figs. 2 and 3). The summary of the geo-statistical characteristics of data on oak damaged by acorn weevil (in 2022 and 2023 in the Dalab Strait) is shown in Tables 1 and 2, along with the kriging maps derived from these variograms (Figs 4 and 5). In this study, data-fitting results indicated that among the twelve sampling dates, nine cases were consistent with the spherical model, while three were consistent with the exponential model. The ratio representing the proportion of the structured component, C₀/(C₀+C) exceeded 50% in all instances. This signifies that the structured component accounts for more than half of the total variance, suggesting that a significant portion of the acorn weevil distribution in the oak forest is governed by spatially correlated processes rather than purely random fluctuations. The performance of the fitted variogram models varied significantly with the sampling period. While most models showed satisfactory coefficients of determination (R²), others showed relatively lower values, indicating that predictive power varied across data sets. Based on the results, the effective range of this pest was relatively high in the exponential models, reaching up to 564.50 meters, whereas it was lower in the spherical models.

Additionally, except for two instances, the explanatory coefficients were greater than 0.506. In all cases, the nugget effect (C₀) was lower than the sill C₀/(C₀+C), which enhances the precision of the kriging interpolation. It is important to note that while a high C₀/(C₀+C) ratio indicates the dominance of the structured component over the random component, it should be interpreted as the relative contribution to the total variance. In the present study, the majority of the cases yielded ratios exceeding 0.501, indicating that the structured component represents the primary driver of the observed spatial variability. Semivariogram analysis demonstrated significant spatial structuring in all samples. The nugget-to-semivariance ratio indicated that 50% of the samples possessed moderate spatial dependence, reflecting suitable correlation between nearby points, while the other 50% showed weak dependence. These results support the effective use of spatial interpolation techniques, such as kriging, to map these features in the study area.

Fig. 2. The percentage of oak fruit-eating insect infestation in the oak forests of Ilam province, 2022 year

Fig. 3. The percentage of oak fruit-eating insect infestation in the oak forests of Ilam province, 2023 year

Fig. 4. Semivariogram and kriging maps of the infestation of oak fruit by the Curculio glandium in Dalab Strait, year 2022-2023

The kriging maps of acorn weevil distribution based on data obtained from pest distribution variography, based on the kriging interpolation method, were used to estimate variable values ​​in unsampled locations (Figs 4 and 5). The spatial distribution and continuous variations in the population density of acorn weevils across the study region are illustrated by the kriging maps. As shown in Figures 4 and 5, the areas with the highest infestation levels are indicated by the red zones in the legend. In both studied years, a higher density of the pest was associated with the 1540–1600 m elevation class (Mean ± SE: 2.10 ± 0.18 in 2022 and 1.67 ± 0.19 in 2023), which showed significant differences compared to other elevation classes (F: 7.51, p < 0.05 in 2022; F: 2.63, p < 0.05 in 2023, as determined by One-Way ANOVA). Due to the temporal-spatial nature of the sampling design, these results reflect temporal–spatial associations rather than a direct causal effect of elevation. The performance of the ordinary kriging model was evaluated using cross-validation. The R² values, which were above 0.506 in most cases, indicate the model's ability to explain the spatial variance. In contrast, the low RMSE values confirm the accuracy of the predictions relative to the observed data. Overall, the cross-validation results demonstrate the efficiency and reliability of the ordinary kriging model for estimating pest density in the study area.

Fig. 5. Semivariogram and kriging maps of the infestation of oak fruit by the Curculio glandium in Dalab Strait, year 2023-2024

Discussion

Libehold et al. (1993), while investigating the distribution of egg masses of the silkworm moth Lymantria dispar (L.) in American forests across years, showed that spatial dependence values varied from year to year and from region to region. It is consistent with the results of this study, indicating that the degree of spatial dependence differs somewhat between the two years under study. Jimenez-Pino et al. (2011) investigated the spatial distribution of the acorn moth, Cydia fagiglandana (Zeller), in holm oak (Q. ilex L.) forests in Spain. Their results showed that the spatial distribution of this pest, as determined using geostatistics, exhibits a spatially dependent pattern, which is consistent with the present study (Jimenez-Pino et al., 2011). The reasons for this similarity can be attributed to the method of sampling oak fruit from the ground and to the collection period, which runs from late November to early December, when most of the fruit is on the ground. In addition, the samples were collected from the lower parts of oak tree crowns, consistent with the sampling method used in this study. In addition, Lasmar et al. (2012), using geostatistics, showed how to determine the spatial distribution and the level of infestation of leaf ant nests (Hymenoptera: Formicidae) in eucalyptus plantations. They concluded that the spatial distribution of this pest is spatially dependent, and their results are consistent with the present study. Cocco et al. (2010) investigated the spatial distribution analysis of the gypsy moth (L. dispar) population (Lepidoptera: Lymantriidae) in oak forests of Sardinia, Italy, using geo-statistical and climate models. Their results showed that the spherical model is the best fit to the data, which is similar to the present research.

Also, the nugget values ​​were close to zero. The C₀/(C₀+C) ratio values were greater than 0.75, which is somewhat close to this research, and the differences in some cases may be attributed to the type of insect studied, weather conditions, and the types of trees mentioned. In the study by Cocco et al. (2010), a total of 40 trees were selected at each site, and the number of eggs was counted in four main directions. Shao et al. (2021) investigated the spatial patterns of herbivorous insects in a forest landscape that included oak trees of Q. variabilis Blume. and Q. aliena Blume. showed that the best model to fit the data was exponential and linear. Also, the values of C₀/(C₀+C) ranged from 0.389 to 0.685, which are lower than the present results and indicate moderate spatial dependence. One of the reasons for the difference with the present study is the research method: in Shao et al. (2021), the study area covered ​​200 ha, which was divided into 10 x 10-meter grids, and sampling was then conducted.

Previous studies, such as Zargaran et al. (2018), have reported higher infestation rates of acorn weevils and acorn moths in higher-altitude classes than in lower ones in the oak forests of Ilam. In the present study, the generated maps indicated that population density was higher in southern-exposed locations and lowest in the northernmost areas. Although microclimatic parameters (e.g., solar radiation and ambient temperature) were not directly measured, the higher pest density in southern-exposed areas is consistent with existing literature suggesting that increased thermal energy may accelerate pest life cycles (Hanna et al., 1982). Consequently, the observed spatial pattern could reflect the influence of microclimatic variations. Furthermore, the maps showed a trend of increasing pest density with elevation during the sampling periods. This observed association between acorn weevil density and higher elevations may reflect a complex combination of environmental and ecological factors, rather than a direct linear effect of elevation itself. While higher elevations are typically associated with cooler temperatures and shorter growing seasons, they may create favorable conditions for this particular pest species depending on acorn quality, habitat structure, and altered pressure from natural enemies. For instance, acorns at higher elevations might exhibit reduced resistance to weevil invasion due to physiological stresses such as lower temperatures, higher wind speeds, or nutritional limitations in the host trees.

Additionally, elevation-related changes in the timing of acorn ripening and dispersal could concentrate weevil densities in areas with greater host availability. The low nugget-to-sill ratio observed in 50% of the samples indicates a predominant spatial structure. According to the geostatistical framework established by Cressie (1993), the observed variability is primarily driven by spatially correlated processes rather than random noise, thereby validating the reliability of spatial interpolation via Kriging.

There is also a change in biological interactions along the altitudinal gradient. Natural enemies of the oak weevil, such as predatory insects and pathogens, may be less abundant or less active at higher elevations, leading to poorer control of the pest population. Also, competition from other insect species or rodents that feed on oak seeds may be less at higher elevations, allowing the oak weevil to increase its population. Microclimatic conditions, including humidity and temperature fluctuations, may also affect survival, larval development, and the ability of the weevil to complete its life cycle in oak cones. Thus, the observed pattern likely results from a complex interaction among the availability and quality of oak seeds, environmental stress on host trees, and changes in the community dynamics of natural enemies and competitors across elevations. In geo-statistical research, the greater the effective range, the greater its utility in pest management. So that in sampling, the distance between samples can be increased up to 75% of the effective range. Therefore, by increasing the effective range within a given field area to estimate a pest's population density, fewer samples will be required (Hassani Pak, 1997). The comparison of this study's results with other studies shows that the geo-statistical method is effective at simulating the population distribution of the acorn weevil, as with other pests.

Nowadays, the use of geostatistics as an important sampling method for arthropod populations is increasing. Geostatistics reduces the number of samples needed to the minimum possible (Liebhold et al., 1993). The generated Kriging maps offer potential for precision management by identifying areas of high pest density. Utilizing these maps could facilitate targeted interventions, potentially reducing insecticide consumption and supporting natural enemy populations. Nonetheless, these findings should be viewed as preliminary indicators, and further studies are needed to validate their effectiveness under practical management conditions and to define specific treatment thresholds. The low R² value in October 2023 indicates that the fitted model only moderately captures the data's empirical spatial structure. This relatively low fit may reflect the strong local variability of the study population, the limited number of sampling points for this period, and the presence of short-range spatial dependence that is difficult to model with a simple theoretical variogram. Also, spatial predictions for this sampling period should be interpreted with caution, primarily as an indication of broad spatial trends rather than as highly accurate local estimates.

Some sampling periods exhibited extremely low nugget values (e.g., C₀ = 0.001). Such values may indicate a highly structured spatial pattern; however, they may also reflect characteristics of the model fitting process. The low nugget values are interpreted as evidence of strong spatial continuity and homogeneity in the distribution of the studied population during those specific sampling periods. This suggests that the factors influencing the spatial distribution were operating consistently at the scales captured by our sampling. Beyond the methodological utility of geostatistical tools, the spatial patterns revealed in this study provide significant ecological insights into the distribution of C. glandium habitats. The spatially dependent pattern of the pest distribution suggests that this species is not randomly distributed, but is highly sensitive to oak trees, especially oak fruit. This spatial structure also suggests that C. glandium population dynamics are closely linked to elevational changes in the region, providing a clearer understanding of how this pest occupies and uses its habitats in the study area.

Conclusions

This study demonstrates the potential of geostatistical modeling to characterize the spatial dynamics of C. glandium within the oak forests of Dalab Strait, Ilam Province. Our analysis indicates that the spatial distribution of the pest can be effectively described using a spherical model, which provided a better fit for the observed spatial structure. The generated Kriging maps identified areas of relatively higher pest density (e.g., in the southern part), providing high-resolution spatial information to support forest management decisions. These findings suggest a possible transition from broad-scale monitoring toward more targeted, site-specific management approaches. By prioritizing monitoring and intervention efforts in identified high-risk zones, resource managers may optimize resource allocation and reduce operational costs. However, it is important to note that these maps should be used as decision-support tools rather than absolute indicators of infestation. Future research should optimize sampling grids to determine the minimum number of samples per hectare required to maintain acceptable mapping accuracy. Furthermore, implementing co-kriging models that incorporate environmental covariates, such as topography, microclimate, or canopy density, could provide deeper insights into the drivers of the observed spatial patterns. Finally, conducting field-based validation of these models under practical management conditions is essential for establishing economic thresholds and evaluating the efficacy of site-specific interventions.

Author's Contributions

Pegah Nazari: investigation, sampling, laboratory experiments; Zahra Rafiei-Karahrodi: methodology, analysis; Majid Mirab-balou: identification, draft preparation, final review and editing; Elham Sanatgar: analysis.

Author's Information

Funding

The Islamic Azad University of Arak, Iran, financially supported this study.

Data Availability Statement

All data generated or analyzed during this study, including detailed methodology, are available from the corresponding author upon reasonable request.

Acknowledgments

This paper is a part of the doctoral thesis of the senior author, which was financially supported by the Islamic Azad University of Arak branch. We want to thank Dr. Behzad Miri, who helped us with data analysis. Sincere thanks to the anonymous reviewers for their valuable comments on the manuscript.

Ethics Approval and Consent to Participate

Insects were used in this study. All applicable international, national, and institutional guidelines for the care and use of animals were followed. This article does not contain any studies with human participants performed by the authors.

Conflict of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Generative AI statement

The authors declare that no Gen AI was used in the creation of this manuscript.