A Small-Scale Investigation into the Viability of Detecting Canopy Damage Caused by Acantholyda posticalis Disturbance Using High-Resolution Satellite Imagery in a Managed Pinus sylvestris Stand in Central Poland

Abstract

As the effects of climate change progressively worsen, many scientists are concerned over the expanding geographic range and impact of forest-defoliating insects. Many are currently pointing to this form of disturbance becoming a key focus of remote sensing research in the coming decades; however, the available body of research remains lacking. This study investigated the viability of detecting and quantifying damage caused to a managed Scots pine forest in central Poland by insect defoliation disturbance using high-resolution multispectral satellite imagery. Observed leaf area index (LAI) values were compared to frass observations (insect detritus) to assess the relationship between LAI and defoliating insect activity across a single life cycle of A. posticalis Mats. Across four managed plots, four vegetative indices (NDVI, GNDVI, EVI, and MSAVI2) were calculated using multispectral satellite imagery from a PlanetScope (PSB.SD instrument) satellite system. Then, 1137 point-sampled digital number (DN) values were extracted from each index, and a correlation analysis compared each to 40 ground-observed LAI data points. LAI was modeled on the basis of NDVI values. Three models were assessed for their performance in predicting LAI. They were fit using a variety of regression techniques and assessed using several goodness-of-fit measures. A relationship between observed LAI and frass observations was found to be statistically significant (p-value = 0.000303). NDVI was found to be the correlated LAI values (rho = 0.612). Model 3, which was based on concepts of the Beer–Lambert law, resulted in the most robust predictions of LAI. All parameters were found to be significant post fitting of the model using a nonlinear least squares method. Despite the success of the Beer’s law model in predicting LAI, detection of A. posticalis damage was not achieved. This was predominately due to issues of resolution and plot condition, among others. The results of this analysis address many interesting facets of remote sensing analysis and challenge the commonly held view of the impeachability of these methods.

1. Introduction

Shifts in global climate are a reality of the modern landscape; however, recognizing and characterizing the consequences these events have on dynamic environmental systems is a formidable task. One such impact of increasingly irregular seasonal weather patterns is the expansion of insect-driven disturbance events in terms of scope and amplitude in regions that previously were untouched by such occurrences.

Forest insect pest disease disturbance (FIPD) is generally considered to be the most devastating biotic disturbance to forest ecosystems at all levels and is the subject of the largest body of modern disturbance research [1]. While the gravity of the disturbance depends on the feeding behavior of the particular insect species, it is universally true that the presence of FIPD affects forest tree photosynthesis, transpiration, and nutrient cycling [2]. There is little doubt that this form of biotic disturbance has broad implications for forest health, regeneration, and carbon sequestration, among other natural processes and systems.

The great web-spinning sawfly Acantholyda posticalis Mats., a significant defoliator species common in pine-dominated forests of Asia and Europe, is a notable example of one such FIPD species. Specifically adapted to feeding on Pinus sp., mass outbreaks have been recorded to have defoliated thousands of hectares of P. sylvestris L. (Scots pine) forest, irrevocably altering local forest systems and timber production [3]. Voolma et al. noted that the range of mass A. posticalis outbreaks, once limited to central Europe, has expanded into the northern Baltic states and Finland. Though considered to be a common species in these regions, no large-scale A. posticalis outbreaks have ever been recorded [4]. In 2008, 250 ha of mid-successional pine forest were defoliated as the result of an A. posticalis outbreak in Estonia, and two years prior, a similar event was recorded in Finland—both representing the first event of its kind to have ever been recorded in these countries. Despite A. posticalis’s ubiquity within the region at large, current methods of monitoring for this insect’s presence and geographic distribution are considered lacking [5,6,7].

Remote sensing (RS) methodologies represent an increasingly popular means of addressing the mentioned knowledge gaps among FIPD researchers. Employment of these tools not only represents a reduction in the cost and time usually associated with conventional FIPD monitoring but also promise validated, accurate results.

In broad terms, this body of work relies on a wide spectrum of remote sensing platforms, from decades-old satellite systems to machine learning, in an effort to better understand the phenomena. Differentiating between the methods currently employed largely comes down to the biological differences of the insects they study and their habitat. However, due to the dynamic interaction of forest insects with forest ecosystems, including the interplay of biotic and abiotic disturbances, RS analysis results must be built on a plethora of correct assumptions and principles to be considered accurate. For instance, an explicit temporal scope as well as the scale at which the studied interactions occur are basic considerations for a robust analysis [8,9]. In the next decade, implementation of novel remote sensing methods will play a key role in improving spatial mapping and impact analysis of FIPD events via the use of more effective scanning units and processing algorithms. According to a literature review conducted by Mngadi et al. [10], the body of remote sensing research regarding FIPD disturbance is growing; however, no research has been conducted regarding either P. sylvestris defoliation or A. posticalis. Such tools will not only broaden our understanding of these destructive outbreaks; they will further facilitate the adaptation of forest management strategy at all scales. There is a consensus in the scientific community that this adaptation is crucial in the face of severe changes in global climate. However, can these methods be fully relied upon to accurately describe forest disturbance in a meaningful way? This study looks to assess this dilemma by investigating a means to detect damage caused by insect defoliation disturbance using a suite of remote sensing and statistical methods currently being investigated by disturbance researchers.

The primary aim of the study is to investigate the viability of detecting and quantifying damage caused to protect Scots pine forests from insect defoliation disturbance using high-resolution multispectral satellite imagery. Modeled leaf area index (LAI) estimates were selected as the metric by which detection was attempted. The primary impact of A. posticalis infestation is defoliation; thus, LAI makes for a logical primary investigation metric. This relationship was thoroughly investigated via robust statistical means and approached by framing the issue into two distinct avenues of investigation:

1. Establishing a clear relationship between observed LAI values and A. posticalis defoliation activity;

2. Assessing the quality of empirical methods of LAI estimation via RS imagery.

In summary, the study adopted the following generalized hypothesis h1: Modeled LAI from remotely sensed imagery of the research plots equals LAI values observed in those plots from the ground, with the null hypothesis being that those groups of LAI values are not equal.

2. Materials and Methods

2.1. Area of Interest

The geographic scope of the data is focused on the Kolumna Forest district near the central Polish city of Łódź (51°37′28.4″ N, 19°18′37.7″ E). The area is categorized as state forest and is managed by The State Forest Holding (Lasy Państwowe) on behalf of the State Treasury (Figure 1). The foundation of this study compares the results of analysis between plots described as infested and control. Infested plots are those that have observed Acantholyda posticalis activity, and the control plot is that in which no activity is detected. Plot 104a represents an area of 7.3 ha, 104b an area of 2.7 ha, 103a an area of 3.7 ha, and 236a represented an area of 1.7 ha. References to infested plots include 103a, 104a, and 104b, and the control is plot 236a, where no A. posticalis activity was observed. Divisions 103 and 104 are a part of a continuously forested area, approximately 5.5 square kilometers in size. Division 236 is roughly six kilometers south of 103 and 104 in a stand that makes up about eight square kilometers. Bisecting the two stands is an area of anthropomorphic development, including a large highway in addition to agricultural areas. The Pinus sylvestris-dominated stands that make up the area of study are characterized as humid pine forest/moderately humid mixed pine forest, respectively, for infested and control plots. The soil profile of all plots was generally classified as course-grained soil, with a high concentration of sand, and the general topography can be described as flat with slight rises and troughs, particularly in plots 236a and 103a.

2.2. Ground Observed Data

Within each of the four plots, five litter trap boxes were placed as subplot centers. These litter boxes with an opening of 2.5 m² were placed along a transect with roughly 10 m between each. From the litter traps, A. posticalis frass and LAI averages were sampled for a period of four months. The temporal scope of the project covered the months of April, May, June, and July of 2021. LAI measurements included, with the exception of April, which had one day of sampling, three collection dates throughout the study period. Frass observations spanned May through July with three measurements per month, with the exception of July, which had two. “Frass” refers loosely to the relatively solid excreta of insects and to certain other related matter and applies to excreted residues of anything that insects have eaten and similarly to other chewed or mined refuse that insects leave behind. This term is commonly used in the entomological literature, and in this study, it represents observed A. posticalis activity [11]. Frass deposits were collected (g/m²) within the litter traps and extrapolated using the volume of the traps to quantify activity throughout the plot. This provided one value per plot per collection date. These values correspond to the time frame and spatial resolution of LAI collection data. LAI is measured in leaf area per ground area (m² m⁻²) but due to the shortening of the units is a unitless measure [12]. Terrestrial observations of LAI were measured using the app “VitiCanopy” (University of Adelaide), which calculates LAI using gap analysis based on upward-facing imagery input into the application. All imagery used to measure LAI was taken on a Samsung S7 [13]. LAI measurements were made one meter above the ground adjacent to the litter traps. For each of the collection dates, 18 measurements were taken at each of the five liter traps, and their values were averaged to provide a single value per plot (n = 40). LAI values were calculated directly in the app and stored in the smartphone’s memory. They were then transferred to the computer for further analysis. Assessment of the relationship between frass and observed LAI was conducted via linear regression using the function “lm” in R package (4.3.3). However, the data were found to be non-normal, and a logarithmic transformation was applied to mitigate this discrepancy. The complete validation data included 40 averaged LAI values and 35 frass measurements. Each average was calculated based on six original measurements. Frass values were compared visually and statistically to observed LAI values to assess their relationship (Table 1).

2.3. Remote Sensing Data

All multispectral imagery was provided by Planet Labs Earth Imaging, a satellite imaging company based in San Francisco, California (Planet Labs PBC, 2022). The satellite system used in this study for analysis was the “Super Dove (PSB.SD)”. The Super Dove consists of a “PSBlue” telescope with a 47-megapixel sensor that was launched in 2017. The PSB.SD instrument is capable of imaging eight spectral bands, the most important of which to this study were the blue (465–515 nm), green (547–583 nm), red (650–680 nm), red edge (697–713 nm), and NIR (845–885 nm) bands of light. Planet Labs offers pixel resolutions up to 30 cm for some instruments; however, due to the specific temporal window in which A. posticalis activity can be observed, the study was limited to instruments capable of 1 m pixel resolution. Products were selected by dates (those which corresponded to dates of ground data collections) and the quality of the imagery. Suitable images in terms of quality were those that possessed little to no cloud cover (≥0.1%) and had no visible image artefacts in or around the area of analysis. In total, 39 images from 2021 were used: three images from April, three from May, four from June, and three from July, and the remainder were unsuitable for remote sensing analysis due to lack of ground observation data and were used as an aid for visual reference. The study aimed for an even spread of dates from the beginning, middle, and end of the month. However, due to a semi-cloudy May and an unusually cloudy July, the representation of dates for these months skewed either to the beginning or end of the desired temporal distribution. The preprocessing and postprocessing of Planet Labs products was carried out prior to downloading. This service was performed automatically for all products of the type used in this study. The coordinate reference system used to project all imagery in the study was WGS 84/UTM Zone 34N (EPSG: 32634). All imagery was converted on-the-fly by QGIS (3.36.1) software to the Polish State Coordinate System PL-1992 (EPSG: 2180).

2.4. Vegetative Indices Selection

Vegetative indices (VI) were the fundamental tool by which this study connected patterns seen in terrestrial observation to patterns indirectly observed via remote sensing methods. Four vegetation indices were compared and used in the study: the normalized difference vegetative index (NDVI), enhanced vegetative index (EVI), modified soil adjusted vegetative index (MSAVI2), and green normalized difference vegetative index (GNDVI). All four add texture or enhance the greenness of vegetation based on varying principles; leveraging these differences, the indices were compared in order to select the best basis for LAI modeling [14,15]. All indices were applied using QGIS’s SCP “bandmath” tool [16]. No corrections were applied to the output rasters.

(1)$\mathrm{p}=1-{\displaystyle \frac{6\sum {\mathrm{d}}_{\mathrm{i}}^2}{\mathrm{n}\left({\mathrm{n}}^2-1\right)}}$

The performance of vegetative indices was assessed using a Spearman’s rank correlation coefficient test (1) [17], where p is Spearman’s rank correlation coefficient, d_i is the difference between the two ranks of each observations, and n is the number of observations. Values were compared based on the significance of the correlation to the observed LAI trends and the resulting rho value (Increased rho value indicates higher level of correlation) [18]. Though vegetative index values do not directly represent LAI values, they do offer an indication of how closely any model will predict observed values.

2.5. Canopy Segmentation

Within the context of this study, canopy segmentation represented a means of increasing sampling accuracy by limiting the sampling extent to the crowns of trees located in sampling plots. All segmentations were generated in R package (4.3.3) using the lidR library. DTMs (digital terrain models) were generated using LAS (LiDAR—light detection and ranging) data supplied by the Polish Head Office of Geodesy and Cartography. All LAS files used for DTM creation were originally captured in 2011 and were clipped to match the extent of the area of interest. The DSMs (digital surface models) used in this study were captured via drone imaging (DJI Mavic 2 Pro, Shenzhen, China). The ground control points were located in the area of interest and measured by the GNSS receiver connected to the Polish network of reference stations (ASG-EUPOS). The horizontal accuracy of measurements did not exceed ±0.03 m [19]. The results of the measurements were used for generating the final orthomosaic and DSM model using Agisoft Metashape (Agisoft LLC, St. Petersburg, Russia). The estimated total control point accuracy was around ±0.05 m. Since the DSM component represents the actual current structure of the forest, it was important that the point clouds were captured concurrently with the dates of terrestrial collections, and thus, we limited the number of canopy segmentations to those dates. The data used for DSM generation were obtained on four dates in 2021 (Table 1).

2.6. Modeling LAI

Using multispectral imagery, the vegetative indices NDVI, EVI, MSAVI2, and GNDVI were calculated within QGIS (QGIS.org). Random sampling points were generated over the bounding area of the vegetative index rasters and then clipped using the canopy segmentation extent, which created 1137 sampling points. Vegetative index digital number values were sampled, averaged, and then compared to ground observation values via regression. Three models were assessed for their ability to predict LAI from NDVI data (2, 3, and 4). Two of the models calculated LAI through correlation with NDVI, and the third was a parameterized model derived from Beer’s law [20]. Model 1 was based on a study conducted in Japanese pine forests [21], model 2 was based on a study in the semi-arid grasslands of central Mongolia [22], and model 3 was derived from a study involving wheat LAI monitoring in Eastern China [23]. The coefficients used in models 1 and 2 were adapted from starting parameter estimates found when fitting model 3.

(2)$\mathrm{M}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l} 1:\mathrm{L}\mathrm{A}\mathrm{I}=-0.0964+1.722\times {\mathrm{I}}_{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}$

(3)$\mathrm{M}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l} 2:\mathrm{L}\mathrm{A}\mathrm{I}=1.53\ast \mathrm{e}\mathrm{x}\mathrm{p}\left(3.27\times {\mathrm{I}}_{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}\right)$

(4)$\mathrm{M}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l} 3:\mathrm{L}\mathrm{A}\mathrm{I}={\displaystyle \frac{-\mathrm{l}\mathrm{n}\left({{\mathrm{a}}_1\times \mathrm{I}}_{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}+{\mathrm{b}}_1\right)}{\left({{\mathrm{a}}_2\times \mathrm{I}}_{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}+{\mathrm{b}}_2\right)}}$

Models 1 and 2 were fitted using linear regression post-logarithmic transformation where a₁, a₂, b₁, and b₂ are parameters, and I represents vegetative indices. Model 3 was fitted using a non-linear least squares method. Model 3 was weighted based on number of observations per plot and treatment (infested/control). Final model selection was based on the adjusted coefficient of determination (R2), root mean square error (RMSE), and mean absolute error (MAE). Pearson residual plots were additionally generated to assess the quality of the selected modeling results. The output of the selected model was then tested for normality, and the observed LAI and predicted LAI results were compared using linear regression at the plot level. Linear regression and non-linear least-squares were calculated with the R packages “lm” and “nls”, respectively (Figure 2).

3. Results

3.1. Ground Data

Mean LAI values for individual plots (104a = 0.951, 104b = 0.833, and 103a = 0.962) differed only slightly between infested plots, while the control plot (no defoliation disturbance activity) 236a displayed the highest mean of 1.169 (Table 2).

A stark pattern of decreasing LAI with increasing frass values was observed with each infested plot and date throughout the typical temporal range of insect activity (Figure 3 and Figure 4). The relationship was found to be significant (p = 0.0003, Table 3).

Canopy segmentation was performed for each date with the available drone imagery from which a DSM could be constructed (Table 1 and Figure 5). The accuracy of exact crown delineation varied based on the relatively divergent conditions of the plots. However, centroids, which were placed in the “cells” and checked visually for accuracy against orthophotos, were highly successful in terms of landing on some section of crown and required no removal of points placed over gaps.

3.2. Vegetative Indices Results

Due to the theoretical design of the vegetative indices, the sampled pixel values varied greatly (Table 2); however, once sampled, they displayed relatively similar patterns of association (Figure 6).

While it is often reported that NDVI suffers from oversaturation, especially in forested areas [24], this issue was not observed in any of the VI-transformed imagery. This is most likely due to the forest structure present in central Poland and the semi-aridness of the environment in comparison to regions at higher latitudes. NDVI and EVI distributions, when separated by plot, displayed similarities in the behavior of each area over time, while GNDVI and MSAVI2 were more similar to each other. Both EVI and NDVI immediately displayed clear differences between plots that GNDVI and MSAVI2 were unable to differentiate (Figure 7).

The results of the Spearman’s rank correlation analysis showed significant correlation (all p-values < 0.01) between all VIs and observed LAI; however, the correlation coefficient (rho) showed a much stronger relationship between NDVI and LAI than the other indices (Table 4). As a result of the correlation analysis, NDVI was selected as the basis on which to model LAI.

3.3. LAI Modeling

As described in the Methods Section, three models were assessed for their ability to predict LAI from NDVI digital number (DN) values (Figure 8).

Results of the various regression analyses (Table 5) revealed significant results in all cases. Pearson’s residual plots were generated for each of the three models (Figure 9).

Plots for M1 and M2 displayed similar results with tight grouping around zero and semi-cone shaped distributions characteristic of heteroscedastic data. Model 3’s residual plot displayed a larger spread of values above and below the zero line than the other plots; however, it displayed grouping to the left of the plot due to several outliers.

Goodness-of-fit measures were used to determine the efficacy of the models (Table 6).

The measures used were adjusted R2, root mean squared error (RMSE), mean absolute error (MAE), and the Akaike information criterion (AIC). Despite the significance of the presented p-values, adjusted R2 displayed weak signals for models 1 and 2, with values of 0.366 and 0.403 respectively, and a stronger value of 0.754 for M3. RMSE and MAE values indicated M3 predictions were subjected to less error than M1 and M2 results (Table 6). AIC values for models 1, 2, and 3 indicated M3 was the best performing model. All goodness-of-fit results found M3 predictions were the most robust of the models assessed. Model 3 mean LAI values for individual plots can be found in Table 2. Generally, the mean values are close to that of the observation dataset; however, they do not exhibit the same dynamic between plots.

An additional regression analysis was performed comparing observed LAI and M3 predicted LAI results to evaluate the ability of the modeled results to assess LAI at the plot level. Results of this further regression indicated a significant relationship between the individual plots and observed LAI for plots 104b and 236a (Table 3). This additional investigation indicated a more complicated relationship between the modeled results, observed LAI, and the RS imagery than initially anticipated and is explored further in the following sections.

4. Discussion

4.1. The Relationship Between Acantholyda posticalis and Ground Observation Data

Results of the analysis performed on the ground-observed data were in line with observations from the literature in regard to the relationship between defoliator disturbance impact and LAI values [25,26]. As mentioned in previous sections, there should be an observable negative relationship between insect defoliation disturbance and forest LAI values; however, despite the large number of studies researching forest damage by insect pests, no studies were found in the literature that investigated this relationship using remote sensing methods for A. posticalis. Frass collection is an effective way to measure the level of activity of A. posticalis during peak disturbance, as their characteristic web-spinning and feeding activity create a large amount of visible detritus, which can be collected and estimated across a plot. Though a relatively rough measure, it could be theoretically calculated down to the individual insect in terms of frass generation in future study. Thus, frass can serve as proxy for a wide range of estimations including but not limited to population density, breeding success, and level of impact. Within the scope of this study, frass measures were used as a proxy for level of activity of A. posticalis. Therefore, according to the literature, results should show a negative relationship between frass levels and LAI values. The study results of the linear regression analysis between frass and LAI indicated a significant negative relationship (p-value = 0.0003); this was further confirmed in the distribution plotting (Figure 4). As such, LAI was used as an indicator of A. posticalis activity, with lower LAI values indicating greater insect activity and vice versa.

In plotting observed LAI distribution by date and plot, a clear pattern was observed for each of the two plot conditions (control and infested), as seen in Figure 4. The control plot represented a typical LAI curve for Pinus sylvestris in the growing season, with values steadily increasing through April and May and peaking around the end of July into August. Dips in this cycle can be attributed in large part to various other environmental factors, such as drought and storm activity. Infested plots, however, demonstrated a pattern of delayed growth indicative of forest stands affected by defoliation disturbance. The pattern is characterized by a similar initial growth period to that of the control plot followed by a flattening of or even reduction in LAI during peak insect feeding. After the insects descend to the forest floor and begin the overwintering process, in late June into July, a period of immediate regrowth can be observed. This pattern was described in de Beurs and Townsend [27] and matches the results observed in this study regarding LAI distribution.

4.2. Influence of Ground Cover on Canopy Analysis in Remote Sensing

Understanding the role ground cover plays in any remote sensing study involving forest systems is crucial. If the focus of the research is on the canopy, ground cover will be a vector of error, which must be addressed in the scope of the study. In the case of this research, the canopy segmentation was designed to mitigate the influence of understory VI values being sampled partially. One clue to the inconsistencies in the VI distributions could be the stark differences that can be observed when comparing the state of the canopies immediately surrounding the sampling area to the sampling area itself (Figure 5). As was the case with the VI analysis, the canopies of plots 103a and 104a displayed similarities that set them apart from the other two plots. The understory also displayed the same characteristics around the segmentation. In both 103a and 104b, the understories were far denser and greener than that of the other two plots (especially 104b). In the case of 236a, the canopy had a substantially higher percentage of cover than the other plots and, correspondingly, a less developed understory. Plot 104b, which demonstrated the lowest values in both observed and VI-derived data, had the sparsest coverage in both the canopy and the understory. The vegetation characteristics of plots 103a and 104b (both a sparse canopy and a thriving understory) were likely to introduce error into the remote sensing analyses, as vegetative indices often give much higher values of “greenness” than those of coniferous tree species [28]. The developing photosynthetic material of a vigorous herbaceous understory species, which has a high degree of leaf coverage in the growing season, is more likely to reflect greater levels of NIR light than that of the much sparser needle coverage of Pinus sylvestris trees, especially those in the midst of an insect defoliation event.

A further point in the consideration of error lays in the method of organizing data or the scale at which analysis is occurring. Developing a method compatible with practical application in forestry most likely necessitates monitoring at the stand level. Since the assumption is that managers would not know if or where defoliation disturbance is occurring, a full-stand RS survey would be the primary indicator if there is any such activity. Within the study’s dataset, observed LAI is organized individually for each plot by date. These assessments, however, were analyzed via Spearman’s correlation analysis at the stand level (all plots). While the organization of the analysis more closely matches the orientation of the overall study’s goal, it may cause loss of information when comparing dynamics amongst the various plots.

4.3. Study Limitations

Accurate remotely sensed LAI measurement has long been a goal of researchers involved in disciplines related to forest and forest-adjacent research. Most LAI estimates for systems like MODIS and Sentinel-2 have globe-spanning maps of LAI, sometimes going back nearly a decade (in the case of Copernicus data). However, spatial resolutions are at best 300 m (LAI300), and the values are often found to be more internally consistent than that of ground-observed LAI values [25]. The limitations evident in these established systems limit their practical application in the field, especially for studying localized insect defoliation disturbance. While this study’s attempt to improve upon LAI estimates using VI-augmented satellite imagery saw encouraging results, several deficiencies persisted as a part of the analysis.

The ground observation sample size used in the validation of the models presented in this study is one issue that must be acknowledged. In the case of LAI measurements, regardless of whether the values were averaged, the sample would still represent a relatively low level of replication, which makes observations at the individual tree level difficult. This most likely accounts for the model distribution and average performance compared to studies that make use of similar modeling designs [10,29]. The effects of this persistent issue limit the generalization of the models and the methods employed in this study. Follow-up analysis with larger ground observation datasets to confirm the trends presented is needed. The benefit of the presented research is that it represents a cost-effective solution to RS-based FIPD monitoring. High-resolution imagery is expensive and limited in terms of geographic scale, whereas large ground observation datasets, which allow for lower spatial resolution imagery and the implementation of machine learning algorithms, are also prohibitively expensive, especially in the context of smaller organizations. Other studies have indicated machine learning algorithms offer much better performance than the algorithm described in this research; however, as hinted at earlier, these algorithms require large ground observation datasets to be validly implemented [9,10].

4.4. Model Construction and Fitting

Model 3 was derived from work conducted by Tan et al. [23] and was based around integrating concepts of NDVI with the Beer–Lambert law. Beer’s law, which states light absorbance is proportional to its direction and length of travel through a sample, depending on concentration or in this case density. This is applicable to predicting LAI, as it applies to an extinction coefficient (K) and the light energy ratio (I/I0). I/I0 represents the light intensity both above and below the canopy, essentially representing light absorbed in transmission through the overstory. K, the extinction coefficient, represents a non-real value of light absorbed through a material. This is represented mathematically according to the formulas which were originally published by Tan et al. [23]:

(5)$\mathrm{I}/{\mathrm{I}}_0={\mathrm{a}}_1\times \mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}+{\mathrm{b}}_1$

(6)$\mathrm{K}={\mathrm{a}}_2\times \mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}+{\mathrm{b}}_2$

(7)$\mathrm{L}\mathrm{A}\mathrm{I}=-\ln \left({\displaystyle \frac{\mathrm{I}}{{\mathrm{I}}_0}}\right)\times {\displaystyle \frac{1}{\mathrm{K}}}$

where a₁, a₂, b₁, and b₂ are unknown parameters partially determined by LOV values; I/I0 is the light energy ratio; and K is the extinction coefficient. In the original research conducted by Tan et al. [23] the parameters, as seen in the above equations, were representative of leaf orientation value (LOV) at different K and I/I0 cutoffs. In adapting this formula to the analysis presented in this study, these parameters were taken as unknown values to be estimated using a nonlinear least squares method (NLS). In implementing this method, these values were estimated algorithmically to produce the least possible sum of squares error (residual sum of squares). In this way, LOV was simultaneously backwards-estimated through K and I/I0 to the best ability of the model. Model 3 parameter estimates (Table 5), based on assessments from de Wit [30], indicate an LOV of a species with spherically orientated leaves, of which P. sylvestris is a member. Models 1 and 2 assume constant K values, a common solution for the few NDVI to LAI models found [9], which is most likely why they failed to perform as well alongside a model that has parameterized these estimates.

The distribution of the modeled results (model 3) seen at the plot level were improved in comparison to raw VI trends but further indicated a lack of capturing all observed nuance. However, when observed at the stand level, the results improved substantially, with a far tighter distribution, reporting at an adj. R2 value of 0.75. The results seen at the stand level are about average in terms of agreement with the results reported in many of the studies attempting to model LAI from VI-adjusted imagery, and these are an indication of the model’s success [10].

This discrepancy between scales, however, is one that is particular to the nature of the study’s research objective. All conventional means of modeling LAI from such imagery do so at coarser geographic scales and usually in agricultural settings. Among these few studies, there is not a single one that has tracked defoliation-induced disturbance. While the primary goal of this study is to develop means of detecting damage caused by FIPD, a secondary goal is assessing the adequacy of modeling LAI for this purpose. This goal can be further synthesized as follows: if LAI can be modeled from vegetative indices successfully, at what scale and how fine is the resolution of the results? The body of research presented in this study indicates success at the stand level, but at the plot level, the story is not as concise. The dataset utilized in this study was built around plots that assessed LAI during a cycle of A. posticalis disturbance, with an additional control plot for referencing healthy stand values. Due to the limitation of observed LAI values, there is relatively low data replication, which lends itself to skewing. Moreover, without the addition of more observations from other similar control plots, there would be no way to adjust to a more evenly distributed dataset without some form of further data transformation. Weights were applied to model 3, which slightly improved overall distribution. This is not an uncommon issue among other studies found in the literature [25,29], as ground-observed LAI datasets are time-consuming to measure manually, and open-data sources are near non-existent.

Ultimately, increasing the number of observations of both infested and control conditions would largely mitigate the issues of outlying data points. However, one available solution to the outliers found in the dataset would be to conduct model fitting individually at the plot level. This would partially alleviate skewing by comparing values on a more normal scale. However, as discussed in the vegetative indices section above, this solution is at odds with the aim of the study. Additionally, such a reorganization would impede on the robustness with which the model estimates variety in LAI values and would severely limit sample size. Therefore, the study further analyzed observations at the plot level after model fitting.

Additional analysis conducted at the plot level using plots 104b and 236a further demonstrated the robustness of model 3’s predictions. As described in the results section, the regression analysis performed on both separated plots returned significant findings (104b p-value = 0.0128; 236a p-value < 2 × 10¹⁶). These results suggest that model 3 shows an ability to predict LAI given specific conditions of the stand. Additionally, these results coincide with the findings discussed in the above vegetative indices section regarding plot condition. The same patterns that were found between plots were continually present in the modeling data and likely the cause of the discrepancies in the model distribution. The error introduced by sampling pixels that partially overlay the understory is due to the spatial resolution of the source imagery.

4.5. Impacts of Resolution

Many studies, such as those presented in Mngadi et al. [10], either make use of higher-resolution multispectral imagery (and often at equally small geographic scale) or lower-resolution imagery, such as Sentinel-2 or Landsat, at a larger geographic scale with larger ground observation datasets. This study did not have the benefit of higher-resolution imagery or a large ground observation dataset; however, it does represent a more cost-effective solution to RS-based FIPD monitoring.

Results from the findings of the VI and LAI modeling analysis indicate that 1 m spatial resolution may be insufficient to effectively detect A. posticalis defoliation. At 1 m spatial resolution, many of the pixels in the source imagery used in this study partially cover forest canopy and the understory layer. Since pixels themselves are averages of DN values sampled by the sensor at the maximum extent of resolution, a high degree of understory coverage would have a proportionate impact on its value. Regardless of where the canopy-segmented sampling points lie, spatial resolutions seen from the source imagery would not allow for mitigation of error from the understory. This is demonstrated by Figure 10, which compares several commonly seen levels of resolution, including that of Sentinel-2 (Figure 10).

UAV-based multispectral sensor systems, as discussed previously, promise a far higher degree of spatial resolution. While UAV systems are limited geographically, they could be easily applied to small-scale analysis methods similar to the ones employed in this study. At such resolutions, issues with understory sampling would be mitigated, and it could even allow for sensing and detection of physical manifestations of A. posticalis, such as web spinning. Furthermore, the application of such platforms, if not already in use, are increasingly available for forest managers in many parts of the world, a trend that is likely to continue.

Ultimately, the results found in the course of analysis indicate deficiencies in establishing clear patterns of defoliation damage as a result of A. posticalis disturbance; the chief cause of which is a lack of spatial resolution. This finding, however, is a technical manifestation of error, and while modeled LAI was inconsistent at this level, the inability to predict insect disturbance was not universal. Plots 104a, 104b, 103a, and 236a were relatively distinct from each other and representative of several scenarios of understory and overstory structure. In this regard, a scale of difficulty can be superimposed on the modeled LAI results, upon which confidence in outputs can be assessed. In practice, plot 104b would be assessed as having sparse canopy coverage and minimum understory with likely poor soil quality. On the proposed scale, 104b would rank lower than that of 103a, with 103a being considered a difficult-to-measure plot with its sparse canopy coverage and a thriving understory. Once ranked, appropriate weights could be applied to each that would affect confidence in the assessment at this scale and indicate where overestimation is likely. These approaches represent valid means by which the limitations of spatial resolution can be overcome.

4.6. Remote Sensing Difficulties

Many studies that investigate RS-based methods of detection or estimation often report claims of statistically significant models that predict their aims with a modestly high degree of success [31,32,33,34,35,36,37,38,39]. However, in some cases, it is evident on closer inspection that these studies have made an incorrect statistical assumption or have modeled RS data based on products originally derived from RS sources in order to describe an observed terrestrial phenomena, etc. It is not the authors’ opinion that this describes even a modest portion of the available body of research or that these types of errors were made in bad faith but that the promise of using this technology is alluring.

The use of remote sensing for biotic disturbance mapping promises a faster, cheaper, and easier result. Accuracy, however, is a promise that comes with many qualifications [40]. Strict procedures and a disciplined approach to sample design, image processing, statistical operations, model fitting, and interpretation (to name a few) are vital to ensuring any degree of accuracy, and all the while, the results found are only as good as the terrestrial observation dataset they are based upon. Even if these criteria are fulfilled, it is still unlikely to find any meaningful signal that demonstrates the success of a given method. That said, in the case of such studies, something was most likely accurately observed, but it may have not been what was originally intended or otherwise wanted by the researchers involved.

RS is not a panacea to the difficulties in describing these dynamic and subtle natural occurrences. All things considered, its use presents a much more challenging prospect than is generally assumed. This study, for one or more of the reasons discussed, most likely did not observe significantly representative modeled LAI values indicating the presence of defoliation disturbance such as that caused by A. posticalis.

4.7. LAI Modeling Confidence

Although the analysis was unable to detect damage caused by A. posticalis disturbance, model 3 (Beer–Lambert model) demonstrated a reasonable ability to estimate LAI from vegetative indices. As with the issues discussed regarding resolution, modeling of LAI was found to be fairly robust at the stand level but failed to capture the much finer patterns seen at the plot level. One consideration is that the dataset employed in this study was not designed around creating a representative sample of normal LAI values but instead describes values of a stand undergoing defoliation disturbance. What the results of the analysis do indicate is that the model may actually predict normal (“healthy stand”) LAI distributions with a more reasonable degree of accuracy.

A “healthy stand”, as seen in plot 236a, with a mostly closed canopy represents a good candidate for the application of the Beer’s law model. Unlike that of plot 103a, the control plot mostly guarantees sampling of the canopy and demonstrates a significantly closer fit with observed values (Table 2). However, any application of the model presented must be assessed alongside the interplay with the dynamics of stand condition and spatial resolution, as these all play a role in the final distributions seen. In summation, while the model is deemed effective at predicting observed LAI, it has serious shortcomings in its application to A. posticalis-driven insect defoliation mapping at a 1 m spatial resolution.

4.8. Adaptation for Forest Management

To adapt a similar method to management practices, regardless of the model of LAI used or the type of vegetative index employed, the workflow would need to satisfy four specific requirements:

1. High-resolution multispectral imagery ≤ 0.5 m along with LiDAR point cloud data (though not presented in this study, a modified workflow that only uses LiDAR sensors could also be utilized);

2. Employment of vegetative indices as a means of augmenting imagery to display underlaying biological responses;

3. Development of a canopy segmentation that is capable of defining and labeling all unique individual crowns;

4. Implementation of a model that accurately predicts LAI from sampled VI values that can later be used as the basis for constructing disturbance mapping.

Evidence from this study suggests that such a workflow would likely produce accurate mapping at higher spatial resolutions and potentially could be adapted to monitor other insect defoliation disturbance species. Furthermore, if all four criteria are met and the if model in question has been sufficiently validated for the region, implementation of such a workflow could be exceedingly fast. Imagery would likely be captured via some form of UAV platform by which stand-wide sampling could be conducted in a matter of days.

4.9. Future Studies

Future studies looking to improve RS-based insect defoliation detection should concentrate primarily on aerial laser scanners (ALS) and drone-imaging systems, as they provide significantly higher spatial resolutions and are appropriate for the required geographic scale. This increased spatial resolution opens up a suite of novel analysis techniques that were precluded by the scope of this study. Statistical approaches to modeling LAI should investigate the application of mixed-effect models, which better allow for assessing effects from categorical variables such as plot and treatment. Additionally, machine-learning and CNN (convolutional neural networks) techniques will likely play an important role in future research on this subject [26]. High-resolution VI and LAI sampling present the possibility of large datasets, for which the application of ML and CNN techniques are uniquely suited, especially considering larger geographic areas than that studied herein.

5. Conclusions

• LAI, as a qualitative measure, is significantly correlated with A. posticalis defoliation disturbance in the context of a managed P. sylvestris stand;

• LAI can be modeled from vegetative index values in a way that is significantly similar to ground-observed LAI values;

• High spatial resolution (≤50 cm) is absolutely required to investigate insect defoliation disturbance using remote sensing methods in the context of small time scales and forestry application.

Author Contributions

Conceptualization, J.S., M.B. and M.S.; methodology, J.S., M.B. and M.S.; software, J.S. and M.B.; validation, J.S., M.B. and M.S.; formal analysis, J.S.; investigation, J.S., M.B. and M.S.; resources, J.S., M.B. and M.S.; data curation, J.S.; writing—original draft preparation, J.S., M.B. and M.S.; writing—review and editing, J.S., M.B. and M.S.; visualization, J.S. and M.B.; supervision, M.B. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data for this research can be shared upon request to the author.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Locator map of the study area. The (A) map presents the sample plot location; the (B) map presents the control plot location. The black dot represents the location of the study area in the context of Poland.

Figure 2. The diagram illustrating the complete workflow of the presented research.

Figure 3. Histogram of log-transformed observed LAI values.

Figure 4. Observed LAI distributions by plot and date (A) and observed frass (g/0.5m2) distributions by plot and date (B).

Figure 5. Example images displaying canopy segmentations at all plots for the date 15 June 2021.

Figure 6. Example images displaying vegetative indices (from left to right) NDVI, GNDVI, EVI, and MSAVI2 at all plots for the date 15 June 2021.

Figure 7. Compound plot displaying VI total pixel DN (digital number) value distribution vs. observed LAI. (Upper) From left to right: (A) NDVI and (B) GNDVI. (Lower) From left to right: (C) EVI and (D) MSAVI2. Rho values, as a result of the Spearman’s rank correlation between VI and observed LAI, are displayed for all plots.

Figure 8. (A) Model 1 LAI predictions vs. observed LAI values. (B) Model 2 LAI predictions vs. observed LAI values. (C) Model 1 LAI predictions vs. observed LAI values. Including: Model fitted, adj.R2, RMSE, and MAE on all graphs.

Figure 9. (A) Pearson’s residual plot displaying model 1 residuals plotted against fitted LAI. (B) Pearson’s residual plot displaying model 2 residuals plotted against fitted LAI 1 (C) Pearson’s residual plot displaying model 3 residuals plotted against fitted LAI.

Figure 10. Spatial Resolution Comparison for Sentinel-2 (A), Planet Labs satellite system (B), and DJI Mavic Pro 2 imagery (C).

Table 1

Dataset of observed values used during the course of study, with reference to plot sampling locations presented in the locator map. Bold text indicates equivalent dates for which DSMs were available.

Data	Observed LAI				Frass (g/0.5 m2)
	Infested Plots			Control Plot	Infested Plots			Control Plot
	103a	104a	104b	236a	103a	104a	104b	236a
24.04.21	0.960	0.847	0.780	0.845	n/a	n/a	n/a	n/a
12.05.21	0.969	1.015	0.764	1.099	0	0.135	0	0.148
24.05.21	0.960	1.042	0.841	0.927	0.826	1.275	0.998	0.442
31.05.21	0.983	0.889	0.856	1.195	4.453	6.542	3.514	0.445
08.06.21	0.977	1.068	0.920	1.197	32.58	51.47	29.38	1.652
15.06.21	0.918	0.793	0.780	1.198	33.27	97.19	46.94	1.833
22.06.21	0.934	0.771	0.701	1.052	5.478	29.057	34.172	2.447
06.07.21	0.973	1.065	0.866	1.385	1.362	7.297	4.139	0.471
13.07.21	0.963	0.980	0.879	1.470	0.732	1.317	1.873	0.735
20.07.21	0.979	1.045	0.941	1.321	n/a	1.35	1.604	0.606

Table 2

Summary statistics of observed frass (g/0.5m²), observed LAI at plots, modeled LAI at plots, and digital number value (DNV) for indexes: NDVI, GNDVI, EVI, and MSAVI2. Bold text presence represents an aid for comparing summary statistics of observed and modeled LAI values at the stand (total) level. N indicates sample size. Mean indicates mean values for said variable. St.Dev indicates standard deviation. Min indicates minimum value, 25%—Q indicates the 25th quartile, and 75%—Q indicates the 75th quartile. Max indicates variable max value.

Variable	Plot	N	Mean	St.Dev	Min	Q1	Q3	Max
Observed Frass (g/0.5 m2)	All Plots	35	12.4	21.6	0	0.7	7.3	97.2
Observed LAI	104a	10	0.951	0.11	0.771	0.847	1.045	1.068
	104b	10	0.833	0.071	0.701	0.78	0.879	0.941
	103a	10	0.962	0.02	0.918	0.96	0.977	0.983
	236a	10	1.169	0.186	0.845	1.052	1.321	1.47
	All Plots	40	0.968	0.165	0.701	0.847	1.045	1.47
M3 Modeled LAI	104a	390	1.043	0.155	0.724	0.916	1.159	1.482
	104b	264	0.843	0.092	0.607	0.81	0.909	1.044
	103a	228	1.276	0.162	0.914	1.184	1.404	1.736
	236a	255	1.248	0.176	0.907	1.133	1.396	1.688
	All Plots	1137	1.005	0.169	0.8	0.884	1.044	1.47
Pixel DNV	NDVI	1137	0.726	0.065	0.589	0.664	0.77	0.809
	GNDVI	1137	2521.71	269.021	1993.06	2348.59	2749.32	2896.86
	EVI	1137	4.505	0.261	3.992	4.284	4.765	4.848
	MAVI2	1137	31.763	55.06	74.032	12.382	44.5	80.128

Table 3

Results of the regression analyses presented in this study. The regression column indicates the variables and models used in the analysis. The estimates column indicates the model’s starting value estimates. The standard error column (SE) indicates the calculated standard error. The tvalue column indicates the t-value score of the given model. The p-value column indicates the p-value score of the model, and asterisks indicate the level of significance of p (p ≤ 0.05 *, p ≤ 0.001 ***).

Regressions	Estimate	SE	t-Value	p-Value
LAI~Frass (g/0.5 m2)	−0.052	0.013	−3.873	0.000303 ***
M1~Observed LAI	0.730	0.029	25.62	<2 × 1016 ***
M2~Observed LAI	0.791	0.029	27.19	<2 × 1016 ***
Predicted LAI (236a)~Observed LAI	1.213	0.008	161.36	<2 × 1016 ***
Predicted LAI (104b)~Observed LAI	1.260	0.098	12.93	0.0128 *

Table 4

Spearman’s rank correlation test results. All VIs were compared against observed LAI values. Correlation results (rho) from bottom to top: NDVI, GNDVI, EVI, and MSAVI2. Asterisks (**) indicate statistically significant rho values.

Correlation	rho
NDVI~Observed LAI	0.612 **
GNDVI~Observed LAI	0.374 **
EVI~Observed LAI	0.536 **
MSAVI2~Observed LAI	0.205 **

Table 5

Results of the nonlinear least squares analyses. The parameter column indicates the specific modeling parameter. The estimates column indicates the model’s starting value estimates. The standard error (SE) column describes calculated standard error values for the model. The t-value column indicates the t-value score for the given model. The p-value column indicates the p-value score of the model, and asterisks indicate the level of significance of p (p ≤ 0.001 ***).

Parameter	Estimate	SE	t-Value	p-Value
a1	−4.247	0.229	−18.54	<2 × 1016 ***
a2	3.747	0.177	21.19	<2 × 1016 ***
b1	4.435	0.185	23.93	<2 × 1016 ***
b2	−3.031	0.143	−21.19	<2 × 1016 ***

Table 6

Goodness of fit results for models 1, 2, and 3. Displaying from left to right: Model number, degrees of freedom (n − 1), Akaike information criterion, adjusted R2, root mean square error (LAI value), and mean absolute error (LAI value).

Models	df	AIC	R2	RMSE	MAE
M1	1135	−1045.528	0.3658	0.152	0.112
M2	1135	−1097.196	0.39397	0.148	0.111
M3	1133	−1942.804	0.74	0.098	0.069

References

1. Seidl, R.; Schelhaas, M.-J.; Lexer, M.J. Unraveling the Drivers of Intensifying Forest Disturbance Regimes in Europe. Glob. Change Biol.; 2011; 17, pp. 2842-2852. [DOI: https://dx.doi.org/10.1111/j.1365-2486.2011.02452.x]

2. Kautz, M.; Meddens, A.J.H.; Hall, R.J.; Arneth, A. Biotic Disturbances in Northern Hemisphere Forests—A Synthesis of Recent Data, Uncertainties and Implications for Forest Monitoring and Modelling. Glob. Ecol. Biogeogr.; 2017; 26, pp. 533-552. [DOI: https://dx.doi.org/10.1111/geb.12558]

3. Voolma, K.; Pilt, E.; Õunap, H. The First Reported Outbreak of the Great Web-Spinning Pine-Sawfly, Acantholyda posticalis (Mats.) (Hymenoptera, Pamphiliidae), in Estonia. For. Stud.; 2009; 50, pp. 115-122. [DOI: https://dx.doi.org/10.2478/v10132-011-0072-z]

4. Viitasaari, M. Sawflies (Hymenoptera, Symphyta) I. Zoosystematics Evol.; 2004; 80, 135. [DOI: https://dx.doi.org/10.1002/mmnz.20040800123]

5. Głowacka, B.; Skrzecz, I.; Bystrowski, C. Reducing the abundance of great pine web−spinning pine sawfly Acantholyda posticalis Mats. in pine stands. Sylwan; 2014; 158, pp. 323-330.

6. Sławski, M.; Mokrzycki, T.; Perliński, S.; Rutkiewicz, A.; Sławska, M. Distribution of wintering pre-imaginal stages of the great web-spinning pine sawfly Acantholyda posticalis Mats. in Scots pine stands being the outbreak centres. Sylwan; 2019; 163, pp. 556-563. [DOI: https://dx.doi.org/10.26202/SYLWAN.2019004]

7. Voolma, K.; Hiiesaar, K.; Williams, I.H.; Ploomi, A.; Jõgar, K. Cold Hardiness in the Pre-imaginal Stages of the Great Web-spinning Pine-sawfly Acantholyda posticalis. Agric. For. Entomol.; 2016; 18, pp. 432-436. [DOI: https://dx.doi.org/10.1111/afe.12172]

8. Raffa, K.F.; Aukema, B.H.; Bentz, B.J.; Carroll, A.L.; Hicke, J.A.; Turner, M.G.; Romme, W.H. Cross-Scale Drivers of Natural Disturbances Prone to Anthropogenic Amplification: The Dynamics of Bark Beetle Eruptions. BioScience; 2008; 58, pp. 501-517. [DOI: https://dx.doi.org/10.1641/B580607]

9. Senf, C.; Seidl, R.; Hostert, P. Remote Sensing of Forest Insect Disturbances: Current State and Future Directions. Int. J. Appl. Earth Obs. Geoinf.; 2017; 60, pp. 49-60. [DOI: https://dx.doi.org/10.1016/j.jag.2017.04.004]

10. Mngadi, M.; Germishuizen, I.; Mutanga, O.; Naicker, R.; Maes, W.H.; Odebiri, O.; Schroder, M. A Systematic Review of the Application of Remote Sensing Technologies in Mapping Forest Insect Pests and Diseases at a Tree-Level. Remote Sens. Appl. Soc. Environ.; 2024; 36, 101341. [DOI: https://dx.doi.org/10.1016/j.rsase.2024.101341]

11. Frost, C.J.; Hunter, M.D. Insect Canopy Herbivory and Frass Deposition Affect Soil Nutrient Dynamics and Export in Oak Mesocosms. Ecology; 2004; 85, pp. 3335-3347. [DOI: https://dx.doi.org/10.1890/04-0003]

12. Zheng, G.; Moskal, L.M. Retrieving Leaf Area Index (LAI) Using Remote Sensing: Theories, Methods and Sensors. Sensors; 2009; 9, pp. 2719-2745. [DOI: https://dx.doi.org/10.3390/s90402719] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/22574042]

13. De Bei, R.; Fuentes, S.; Gilliham, M.; Tyerman, S.; Edwards, E.; Bianchini, N.; Smith, J.; Collins, C. VitiCanopy: A Free Computer App to Estimate Canopy Vigor and Porosity for Grapevine. Sensors; 2016; 16, 585. [DOI: https://dx.doi.org/10.3390/s16040585] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27120600]

14. Marques Ramos, A.P.; Prado Osco, L.; Elis Garcia Furuya, D.; Nunes Gonçalves, W.; Cordeiro Santana, D.; Pereira Ribeiro Teodoro, L.; Antonio Da Silva Junior, C.; Fernando Capristo-Silva, G.; Li, J.; Henrique Rojo Baio, F. et al. A Random Forest Ranking Approach to Predict Yield in Maize with Uav-Based Vegetation Spectral Indices. Comput. Electron. Agric.; 2020; 178, 105791. [DOI: https://dx.doi.org/10.1016/j.compag.2020.105791]

15. Morales-Gallegos, L.M.; Martínez-Trinidad, T.; Hernández-de La Rosa, P.; Gómez-Guerrero, A.; Alvarado-Rosales, D.; Saavedra-Romero, L.D.L. Tree Health Condition in Urban Green Areas Assessed through Crown Indicators and Vegetation Indices. Forests; 2023; 14, 1673. [DOI: https://dx.doi.org/10.3390/f14081673]

16. Congedo, L. Semi-Automatic Classification Plugin Documentation. Release 5.0.0.1; 2016; 4, 29.

17. Daniel, W.W.; Cross, C.L. Biostatistics: A Foundation for Analysis in the Health Sciences; 11th ed. Wiley: Hoboken, NJ, USA, 2018; ISBN 978-1-119-49657-1

18. Akoglu, H. User’s Guide to Correlation Coefficients. Turk. J. Emerg. Med.; 2018; 18, pp. 91-93. [DOI: https://dx.doi.org/10.1016/j.tjem.2018.08.001]

19. Oruba, A.; Ryczywolski, M.; Wajda, S.; Bosy, J.; Graszka, W.; Leończyk, M.; Somla, J. The Implementation of the Precise Satellite Positioning System ASG-EUPOS. Proceedings of the ENC-GNSS 2008 Conference Proceedings, EUGIN; Toulouse, France, 23–25 April 2008.

20. Mamouei, M.; Budidha, K.; Baishya, N.; Qassem, M.; Kyriacou, P.A. An Empirical Investigation of Deviations from the Beer–Lambert Law in Optical Estimation of Lactate. Sci. Rep.; 2021; 11, 13734. [DOI: https://dx.doi.org/10.1038/s41598-021-92850-4]

21. Saito, K.; Ogawa, S.; Aihara, M.; Otowa, K. Estimates of LAI for Forest Management in Okutama. Proceedings of the 22nd ACRS Conference on Remote Sensing; Singapore, 5–9 November 2001; Rissho University: Kumagaya, Japan, 2001; Volume 1, pp. 600-605.

22. Fan, L.; Gao, Y.; Brück, H.; Bernhofer, C. Investigating the Relationship between NDVI and LAI in Semi-Arid Grassland in Inner Mongolia Using in-Situ Measurements. Theor. Appl. Clim.; 2009; 95, pp. 151-156. [DOI: https://dx.doi.org/10.1007/s00704-007-0369-2]

23. Tan, C.-W.; Zhang, P.-P.; Zhou, X.-X.; Wang, Z.-X.; Xu, Z.-Q.; Mao, W.; Li, W.-X.; Huo, Z.-Y.; Guo, W.-S.; Yun, F. Quantitative Monitoring of Leaf Area Index in Wheat of Different Plant Types by Integrating NDVI and Beer-Lambert Law. Sci. Rep.; 2020; 10, 929. [DOI: https://dx.doi.org/10.1038/s41598-020-57750-z]

24. Olsson, P.-O.; Kantola, T.; Lyytikäinen-Saarenmaa, P.; Jönsson, A.; Eklundh, L. Development of a Method for Monitoring of Insect Induced Forest Defoliation—Limitation of MODIS Data in Fennoscandian Forest Landscapes. Silva Fenn.; 2016; 50, 1495. [DOI: https://dx.doi.org/10.14214/sf.1495]

25. Zhao, Y.; Chen, X.; Smallman, T.L.; Flack-Prain, S.; Milodowski, D.T.; Williams, M. Characterizing the Error and Bias of Remotely Sensed LAI Products: An Example for Tropical and Subtropical Evergreen Forests in South China. Remote Sens.; 2020; 12, 3122. [DOI: https://dx.doi.org/10.3390/rs12193122]

26. Rosso, P.; Nendel, C.; Gilardi, N.; Udroiu, C.; Chlebowski, F. Processing of Remote Sensing Information to Retrieve Leaf Area Index in Barley: A Comparison of Methods. Precis. Agric.; 2022; 23, pp. 1449-1472. [DOI: https://dx.doi.org/10.1007/s11119-022-09893-4]

27. Debeurs, K.; Townsend, P. Estimating the Effect of Gypsy Moth Defoliation Using MODIS. Remote Sens. Environ.; 2008; 112, pp. 3983-3990. [DOI: https://dx.doi.org/10.1016/j.rse.2008.07.008]

28. Eriksson, H.M.; Eklundh, L.; Kuusk, A.; Nilson, T. Impact of Understory Vegetation on Forest Canopy Reflectance and Remotely Sensed LAI Estimates. Remote Sens. Environ.; 2006; 103, pp. 408-418. [DOI: https://dx.doi.org/10.1016/j.rse.2006.04.005]

29. Duarte, A.; Borralho, N.; Cabral, P.; Caetano, M. Recent Advances in Forest Insect Pests and Diseases Monitoring Using UAV-Based Data: A Systematic Review. Forests; 2022; 13, 911. [DOI: https://dx.doi.org/10.3390/f13060911]

30. De Wit, C.T. Photosynthesis of Leaf Canopies; PUDOC: Wageningen, Germany, 1965; pp. 1-57.

31. Jia, K.; Liang, S.; Zhang, L.; Wei, X.; Yao, Y.; Xie, X. Forest Cover Classification Using Landsat ETM+ Data and Time Series MODIS NDVI Data. Int. J. Appl. Earth Obs. Geoinf.; 2014; 33, pp. 32-38. [DOI: https://dx.doi.org/10.1016/j.jag.2014.04.015]

32. Li, Y.; Li, M.; Li, C.; Liu, Z. Forest Aboveground Biomass Estimation Using Landsat 8 and Sentinel-1A Data with Machine Learning Algorithms. Sci. Rep.; 2020; 10, 9952. [DOI: https://dx.doi.org/10.1038/s41598-020-67024-3]

33. Lovynska, V.; Sytnyk, S.; Stankevich, S.; Holoborodko, K.; Tkalich, Y.; Nikovska, I.; Bandura, L.; Buchavuy, Y. Aboveground Biomass Estimation in Conifer and Deciduous Forests with the Use of a Combined Approach. Biosys. Divers.; 2024; 32, pp. 210-216. [DOI: https://dx.doi.org/10.15421/012422]

34. Maselli, F. Monitoring Forest Conditions in a Protected Mediterranean Coastal Area by the Analysis of Multiyear NDVI Data. Remote Sens. Environ.; 2004; 89, pp. 423-433. [DOI: https://dx.doi.org/10.1016/j.rse.2003.10.020]

35. Miraki, M.; Sohrabi, H.; Fatehi, P.; Kneubuehler, M. Detection of Mistletoe Infected Trees Using UAV High Spatial Resolution Images. J. Plant Dis. Prot.; 2021; 128, pp. 1679-1689. [DOI: https://dx.doi.org/10.1007/s41348-021-00502-6]

36. Motohka, T.; Nasahara, K.N.; Oguma, H.; Tsuchida, S. Applicability of Green-Red Vegetation Index for Remote Sensing of Vegetation Phenology. Remote Sens.; 2010; 2, pp. 2369-2387. [DOI: https://dx.doi.org/10.3390/rs2102369]

37. Santoro, M.; Cartus, O.; Carvalhais, N.; Rozendaal, D.M.A.; Avitabile, V.; Araza, A.; De Bruin, S.; Herold, M.; Quegan, S.; Rodríguez-Veiga, P. et al. The Global Forest Above-Ground Biomass Pool for 2010 Estimated from High-Resolution Satellite Observations. Earth Syst. Sci. Data; 2021; 13, pp. 3927-3950. [DOI: https://dx.doi.org/10.5194/essd-13-3927-2021]

38. Schrader-Patton, C.; Grulke, N.; Bienz, C. Assessment of Ponderosa Pine Vigor Using Four-Band Aerial Imagery in South Central Oregon: Crown Objects to Landscapes. Forests; 2021; 12, 612. [DOI: https://dx.doi.org/10.3390/f12050612]

39. Usta, A.; Yilmaz, M. Prediction of Pine Mistletoe Infection Using Remote Sensing Imaging: A Comparison of the Artificial Neural Network Model and Logistic Regression Model. For. Pathol.; 2023; 53, e12783. [DOI: https://dx.doi.org/10.1111/efp.12783]

40. Lacasa, J.; Hefley, T.J.; Otegui, M.E.; Ciampitti, I.A. A Practical Guide to Estimating the Light Extinction Coefficient with Nonlinear Models—A Case Study on Maize. Plant Methods; 2021; 17, 60. [DOI: https://dx.doi.org/10.1186/s13007-021-00753-2]