Application of Remote Sensing Data for Assessment of Bark Beetle Attacks in Pine Plantations in Kirkovo Region, the Eastern Rhodopes

: Intensive forest afforestation with native pine species was developed in the 1960s on degraded and deforested lands in the region of the Eastern Rhodopes (south-eastern Bulgaria). Severe damage by wet snow was registered in the coniferous forests of the Rhodopes in March 2015. In the following years, bark beetle attacks were registered on the broken and felled fresh wood. As a result, bark beetle infestation spots appeared in the pine plantations. In the period 2019–2021, damage caused by bark beetles was assessed in the region of State Forestry Kirkovo (the Eastern Rhodopes, south-eastern Bulgaria). An integrated approach using the data of the information system of the Executive Forest Agency (ISEFA), remote sensing data obtained by an “eBee SQ” unmanned aerial vehicle (UAV) equipped with a “Parrot Sequoia” multispectral camera, and subsequent terrestrial observations, was applied. ISEFA data showed that there was no serious damage caused by abiotic and biotic factors in the pine forests of SF Kirkovo until 2014. Snow damage in 2015 affected 513 ha of pine plantations, and bark beetle infestations reached up to 1316 ha in 2016. In 2019, a total of 226.87 ha of pine plantations were captured in three localities—Fotinovo, Kirkovo, and Kremen. The relative share of damage caused by bark beetles was greater in P. sylvestris plantations (15.3–23.0%), compared to damage in P. nigra (2.3%). Four different categories of normalised difference vegetation index (NDVI) were separated in bark beetle infestation spots—living trees, dead trees, grass and shrub vegetation, stones and rocks. The NDVI values in locations with living trees varied between 0.500 (spaces between tree crowns) and 0.700 (central part of the crown projection) (an average of 0.617). In the locations with dead trees, the average values of NDVI of lying trees was 0.273, and in standing trees, NDVI varied between 0.275 (central part of crown projections) and 0.424 (spaces between tree crowns). In the locations with grass and shrub vegetation, stones and rocks, the average NDVI was 0.436 and 0.329, respectively. In the ﬁeld study, average defoliation of 31.2–32.3% was registered in P. sylvestris plantations, and 47.4% in P. nigra plantations. Defoliations mainly were caused by pine processionary moth ( Thaumetopoea pityocampa ) and fungal pathogens ( Dothistroma septosporum and Lecanosticta acicola ). The damage was caused by Ips acuminatus (in P. sylvestris only), and I. sexdentatus , Tomicus piniperda and T. minor (in P. sylvestris and P. nigra ). Infestations by other xylophages, such as Phaenops cyanea , Rhagium inquisitor , and Pissodes spp., were also found on pine stems.


Introduction
In recent years, adverse phenomena of climate change such as prolonged droughts, floods, and storms have occurred. Under such circumstances, forest stands are subject to more frequent and intense stress and deterioration due to the combination of unfavourable climatic conditions and the impact of damaging biotic factors, mostly bark beetle insect pests [1]. The monitoring of forest health status provides baseline data on the distribution, occurrence and harmfulness of biotic agents or abiotic factors.
In the 1960s, intensive forest planting was developed in the region of the Eastern Rhodopes (south-eastern Bulgaria) under the process of afforestation on degraded and deforested lands. Formerly, severely disturbed habitats had a limiting effect on the restoration of native forest stands in this region. This fact prompted the beginning of large-scale afforestation with species suitable for the deteriorating habitat conditions in the region of Regional Directorate of Forestry (RDF) Kardzhali. From the total forest plantation area (16,820 ha), about 34% (5820 ha) was afforested with anti-erosion functions [2]. Pinus nigra Arn. and Pinus sylvestris L. were the main species used in forest plantations to relieve the eroded terrains in the Eastern Rhodope Mountains. Their rapid growth and tolerance of soil conditions make them suitable for plantations in these areas. However, the increasing age of plantations in these habitats led to an increase in moisture deficits, deterioration of the physiological condition of plants and a reduction in their resistance to abiotic and biotic effects [3]. These factors had an adverse effect on forests, especially on the pine monocultures in the lower forest belt [4].
In recent decades, forest mortality has increased among many tree species in the countries of the Mediterranean region due to dry and warm conditions [5]. In the region of south-eastern Bulgaria, the trend of decreasing precipitation during the vegetation period has an adverse effect on forests, reducing their resistance to biotic and abiotic factors. This creates conditions for the increased population density of biotic agents such as wood-boring insect pests or fungal pathogens. This applies to the region of the RFD of Kardzhali to a great extent, where most of the pine plantations were developed outside of their natural range [2].
Establishing the health status of forests is an important requirement for the application of measures for protection of biological diversity and control of insect pests and fungal pathogens. Reliable, long-term, global-scale forest health monitoring, likely combining remote-sensing and ground-based measurements in a methodologically coordinated and consistent manner, is needed to accurately determine the status and trends of forest stress and mortality [5]. The monitoring of forest health status provides baseline data on the distribution, occurrence and harmfulness of biotic agents or abiotic factors. Assessment and monitoring of forest health represents a key point of environmentally policy and management of forest resources [6].
Unmanned aerial vehicles (UAVs) have shown promising possibilities for the collection of remote sensing data. Furthermore, they are flexible in regard to height and flight time, and produce high-resolution results that are appropriate for plant observation [7]. In Bulgaria, UAV images were used for urban green ecosystem monitoring in Karlovo [8], species classification and assessment of the vital status of forest stands [9], and mapping of Ips typographus infestations in two Norway spruce mountain reserves in the Western Balkan Range [10,11].
Natural disturbances (windthrows, wildfires, and snow and ice damage) have negative effects on the structure and function of forest ecosystems. In mountain areas in Bulgaria, windthrow and Ips typographus outbreaks play important roles in shaping old Norway spruce-dominated forests [12]. In 2001, 60 ha of natural stands of Picea abies in the Bistrishko Branishte reserve on Vitosha Mountain were affected by windthrow. Between 2003 and 2008, I. typographus attacked 200 ha of old Norway spruce forests in the periphery of the windthrow [13,14]. The outbreaks of spruce bark beetle on large, steep and inaccessible terrains in the mountains can be easily observed and mapped using remote sensing techniques and terrain verification [11].
In the present study, an integrated approach using registered data for biotic or abiotic damage, remote sensing data, and subsequent field study, was applied for assessment of bark beetle attacks in pine plantations in the region of State Forestry (SF) Kirkovo.

Study Area
Damage caused by bark beetles in pine plantations was assessed in the period 2019-2021 in the region of SF Kirkovo, located in the southernmost part of the Eastern Rhodopes, Bulgaria ( Figure 1). The study area was selected because of the severe damage caused by wet snow in the spring of 2015 and the subsequent strong attacks by xylophagous insects.
The relief on the territory is mainly hilly foothill, as in the lowest parts it turns into plain hilly, and mountainous to the ridge of the state border with Greece. The terrain of SF Kirkovo is located at an altitude of 260 m to 1463 m. The areas is relatively close to the Aegean Sea. The winter months in the region of Kirkovo are relatively mild, with average temperatures in January above 0 • C. The summer is sunny, hot and dry. The average temperatures vary between 11 • C and 13.2 • C, and the average annual amount of precipitation is 521 mm [15].
The studied pine plantations are located in the lands of Fotinovo, Kirkovo and Kremen villages (Table 1). Two plantations (Kirkovo and Kremen) were afforested mainly with P. sylvestris and single trees of P. nigra, and the third one, Fotinovo, mainly with P. nigra and single trees of P. sylvestris. The density index of the plantations was quite high (0.8-0.9). The age of the plantations varied between 42 yeas (Kremen) and 67 years (Fotinovo).

Data from Information System of Executive Forest Agency in Bulgaria
The registered damage (abiotic windthrows, snowfalls, icebreaks, etc.) and bark beetle attacks for the period 2011-2020 in the territory of SF Kirkovo were analysed on the basis of data obtained from the information system of the Executive Forest Agency (ISEFA) in Bulgaria (https://system.iag.bg/, accessed on 28 February 2022).

UAV Images
A "Parrot Sequoia" multispectral camera, integrated with a specialised professional eBee "Flying Wing" UAV system, was used to take high-resolution images in four channels of the electromagnetic spectrum: green (530-570 nm spectral band), red (640-680 nm), red edge (730-740 nm) and near infrared (770-810 nm), as well as in a standard RGB channel. The camera was equipped with a solar radiation sensor that served for calibration of the obtained reflex images.
The flights were carried out on 3 December 2019. The ground sampling distance (GSD) of the RGB images was 7 cm, while the GSD of the multispectral images was 20 cm, which resulted in a flight altitude of 127.4 m above the ground surface.
To assess the health status of pine stands, the normalised difference vegetation index (NDVI) was used, obtained by digital mixing of imagery captured in the red (R) and near-infrared range (NIR) [16]. It was calculated according to the following standard Formula (1): NDVI is a spectral index used to measure the relationship between the near-infrared range and red bands. It varies between −1.0 and +1.0 and indicates the vegetation condition due to chlorophyll absorption of specific wavelengths of light within the red spectral range and high reflectance within the near-infrared range.

Processing of Remote Sensing Data
The UAV-derived data were processed by a Pix4D platform that yielded an orthophoto mosaic and an NDVI raster map of the captured areas ( Figure 2).

Detecting of Bark Beetle Infestation Spots
We hypothesized that the bark beetle attacks occurred after the wet snow damage in 2015. The presence of bark beetle spots (trees possibly affected by bark beetles) were studied in GIS on the basis of orthophotos taken in 2011 and maps generated from photographic surveys with UAV in 2019.
The data were loaded into QGIS 3.20.0 (a free and open-source program), further refined and classified, clearly distinguishing between healthy pine plantations and infested or dead trees. Subsequently, the classified data from NDVI in combination with the RGB image were used to accurately locate and determine the presence and boundaries of bark beetle spots ( Figure 3).

Detecting of Bark Beetle Infestation Spots
We hypothesized that the bark beetle attacks occurred after the wet snow damage in 2015. The presence of bark beetle spots (trees possibly affected by bark beetles) were studied in GIS on the basis of orthophotos taken in 2011 and maps generated from photographic surveys with UAV in 2019.
The data were loaded into QGIS 3.20.0 (a free and open-source program), further refined and classified, clearly distinguishing between healthy pine plantations and infested or dead trees. Subsequently, the classified data from NDVI in combination with the RGB image were used to accurately locate and determine the presence and boundaries of bark beetle spots ( Figure 3).

Detecting of Bark Beetle Infestation Spots
We hypothesized that the bark beetle attacks occurred after the wet snow damage in 2015. The presence of bark beetle spots (trees possibly affected by bark beetles) were studied in GIS on the basis of orthophotos taken in 2011 and maps generated from photographic surveys with UAV in 2019.
The data were loaded into QGIS 3.20.0 (a free and open-source program), further refined and classified, clearly distinguishing between healthy pine plantations and infested or dead trees. Subsequently, the classified data from NDVI in combination with the RGB image were used to accurately locate and determine the presence and boundaries of bark beetle spots ( Figure 3).

Mapping of Bark Beetle Spots
Shape files and raster masks of the studied pine plantations were created on generated NDVI maps, and the transformation (raster-vector) of the layers was performed. After checking the typology and correcting the erroneous input-output data, vectorization and automatic mapping of areas with relatively homogeneous vegetative characteristics were performed. Shaded areas, openings in the tree canopy (small meadows), as well as single dying or dead trees were excluded by filtering mapped areas of up to 16.5 m 2 (average crown area of 3 trees). Finally, four zones were identified in bark beetle spots:

Mapping of Bark Beetle Spots
Shape files and raster masks of the studied pine plantations were created on generated NDVI maps, and the transformation (raster-vector) of the layers was performed. After checking the typology and correcting the erroneous input-output data, vectorization and automatic mapping of areas with relatively homogeneous vegetative characteristics were performed. Shaded areas, openings in the tree canopy (small meadows), as well as single dying or dead trees were excluded by filtering mapped areas of up to 16.5 m 2 (average crown area of 3 trees). Finally, four zones were identified in bark beetle spots: healthy (living) trees; dead trees; grass and shrub vegetation; stones and rocks. GIS databases and interactive composite maps were created for each of the studied zones affected by bark beetle spots within the areas.

Mapping of Bark Beetle Spots
Shape files and raster masks of the studied pine plantations were created on generated NDVI maps, and the transformation (raster-vector) of the layers was performed. After checking the typology and correcting the erroneous input-output data, vectorization and automatic mapping of areas with relatively homogeneous vegetative characteristics were performed. Shaded areas, openings in the tree canopy (small meadows), as well as single dying or dead trees were excluded by filtering mapped areas of up to 16.5 m 2 (average crown area of 3 trees). Finally, four zones were identified in bark beetle spots: healthy (living) trees; dead trees; grass and shrub vegetation; stones and rocks. GIS databases and interactive composite maps were created for each of the studied zones affected by bark beetle spots within the areas.
Ten separate categories in the bark beetle spots were studied in detail: (1) living trees (central part of the crown projection); (2) living trees (peripheral part of the crown projection); (3) living trees (spaces between tree crowns); (4) dead lying trees; (5) dead standing trees (central part of the crown projection); (6) dead standing trees (peripheral part of the crown projection); (7) dead standing trees (spaces between crowns); (8) grass and shrub vegetation; (9) stones and rocks without mosses; (10) stones and rocks with mosses ( Figure 4). The NDVI values for different categories were generated on the basis of the boundary polygon GIS layers. Data on areas captured, studied, and affected by bark beetles, as well as the number of bark beetle spots, were obtained from the attribute tables of the created shapefiles.

Field Verification
An expedition on 8-9 June 2020 was carried out in the captured areas. A total of six sample plots were set up (four in Fotinovo, one in Kirkovo, and one in Kremen). In each sample plot, the defoliation of 40 trees was assessed according to the methodology of ICP "Forest" [17].
A complex of bark beetles was identified by the types of damage on dead lying and standing trees in bark beetle spots [18,19]. At least 30 trees were studied in each locality. Samples of stems, branches and needles were collected to study the insect pests and fungal pathogens in laboratory conditions at the Forest Research Institute in Sofia.
The collected needles showed typical symptoms of Dothistroma needle blight and Brown spot needle blight diseases [20]. The samples were examined using a Leica stereomicroscope at magnifications up to ×40, and a Zeiss light microscope at ×400 and ×1000 magnification. Needles were immersed in 96% ethanol for 10 s, then rinsed in sterile water. Conidiomata formed on the needles were aseptically excised, rolled onto 2% malt extract agar medium (15 g agar BD Difco™ per litre of distilled water with 100 mg/L streptomycin), in order to release conidia from the conidiomata as described by Barnes et al. [21]. The pure cultures were examined after 2-4 weeks of cultivation at 24 • C. The disease-causing fungi were identified using morphological characters and conidiogenesis of their anamorphic stage by direct examination of infected needles and pure cultures isolated on malt extract agar [20,21].

Statistical Analysis
Statistical evaluation of NDVI values was performed using descriptive statistics of Statistica 12.0 for Windows (StatSoft). The results are cited together with their standard error. To compare the means, the t-test for independent samples was applied.

Damage Caused by Abiotic Factors and Bark Beetles in the Region of SF Kirkovo
The data from the information system of EFA showed that for the period 2011-2014, the damage from abiotic factors in the area of SF Kirkovo varied from 0 to 18.7 ha, and no significant attacks caused by bark beetles were registered ( Figure 5). The data obtained from the information system showed that the most damage was caused by Tomicus piniperda (80.1%), followed by Ips acuminatus (13.6%), Ips sexdentatus (5.9%) and other unidentified bark beetle species (0.4%) ( Figure 6).

Remote Data Analysis and Mapping of Bark Beetle Spots
The analysis of orthophotos taken in 2011 showed no significant bark beetle damage within the studied areas ( Figure 7A,C). In 2019, many dead trees in Fotinovo and bark beetle spots in Kirkovo were observed on the UAV photos ( Figure 7B,D). A total of 226.87 ha were captured, of which 122.22 ha represented forest areas that were surveyed ( Table 2). The number of bark beetle spots varied from 115 (Kremen) to 148 (Fotinovo), and the affected areas ranged from 1.68 ha (Fotinovo) to 6.57 ha (Kirkovo). The average area of bark beetle spots was the largest in the region of Kirkovo, followed by Kremen and Fotinovo. The damage caused by bark beetles was greater in P. sylvestris plantations (Kirkovo and Kremen), compared to that in P. nigra (Fotinovo) (Figure 8). The relative share of bark beetle spots in the region of Kirkovo and Kremen was 23% and 15.3%, respectively, while the bark beetle spots in the region of Fotinovo occupied only 2.3% ( Table 2).

NDVI Values of Studied Categories in the Bark Beetle Spots
The NDVI values in locations with living trees varied between 0.500 (spaces between crowns) and 0.788 (central part of the crown projection). The mean values in the central and peripheral parts of the crown projection of living trees were 0.748 ± 0.002 and 0.645 ± 0.002, respectively, and in the spaces between tree crowns, 0.554 ± 0.001 (Table 3).
The differences in NDVI of living trees were statistically proven for the central and peripheral parts of crown projections and the spaces between crowns (p < 0.001).
In the case of dead trees, there were no significant differences (p = 0.70) in the values of the NDVI of lying trees (0.273 ± 0.002) and the central part of the crown projection of standing trees (0.275 ± 0.002), in contrast to the peripheral part of the crown projection (0.344 ± 0.001) and spaces between crowns (0.424 ± 0.001), where the differences were statistically significant (p < 0.001).
Despite the close mean values of NDVI of spaces between crowns of dead trees (0.424 ± 0.001) and NDVI of grass and shrub vegetation (0.436 ± 0.004), the differences were statistically proven (p < 0.001). Statistically significant differences (p < 0.001) were also found between the locations with grass and shrub vegetation and stones and rocks, both without (0.193 ± 0.001) and with mosses (0.345 ± 0.002). In general, the differences in the mean NDVI values from all subcategories for living trees (0.617 ± 0.003) and NDVI of dead trees (0.343 ± 0.002), grass and shrub vegetation (0.436 ± 0.004) and stones and rocks (0.329 ± 0.004) were statistically proven (p < 0.001). The difference in NDVI of dead trees and stones and rocks was also statistically proven (p = 0.05).

Field Study
The mean defoliation in the four sample plots of P. nigra in the region of Fotinovo varied between 33.33% and 76%. In P. sylvestris plantations, the average defoliation was 32.3% in the sample plot in Kirkovo, and 31.2% in the sample plot in Kremen. The defoliation was caused by larval feeding by the pine processionary moth, Thaumetopoea pityocampa (Denis and Schiffermüller, 1775) (Lepidoptera: Notodontidae) and infestation by fungal pathogens Dothistroma septosporum and Lecanosticta acicola.

Discussion
In recent years, forests have been strongly impacted by the development of diseases and pest attacks, fires and adverse climatic events. In Bulgaria, outbreaks of pests and infections with invasive pathogens have the potential to destroy a large number of trees, in some cases almost all individuals of a local tree species [22]. Among the abiotic factors, damage caused by wind, wet snow and ice has the most negative impact. Snow-and windthrow stem breakage often provoke attacks by bark beetle species from Ips and Tomicus genera [23][24][25][26][27]. Although Tomicus species are less aggressive, they can cause significant damage not only to stems but also to shoots [27,28]. As a result, physiological weakening of trees was observed; the impact caused by invasive pests and pathogens can become devastating, covering vast areas of forests and posing a threat to economically important tree species. The assessment of the harmful impact and spread of the most important insect pests is essential for making decisions about their control. It is necessary to take into account the prevalence of pests, as well as their specific impacts on the affected forests.
Insect pest outbreaks can have severe ecological, socioeconomic, and forest carbon balance impacts by causing mortality and reducing the growth of trees over extensive areas [29]. Information on the extent and severity of tree mortality within plantations is required for a variety of forest management activities, including whether to undertake a non-commercial thinning operation or to salvage log an area, and quantifying the impact on resource volume predictions [30].
The main reason for the deterioration of the health status of pine plantations in the region of Kirkovo is the damage caused by heavy snow in 2015. The presence of a huge amount of fresh wood contributed to the proliferation of bark beetles and provoked severe outbreaks in the Rhodope Mountains and other regions in southern Bulgaria. At the end of 2016, over 1.5 million damaged coniferous trees were reported in the forests of the Western Rhodopes (unpublished). It is well known that the elimination of breeding places of bark beetles (felling and removal of infested trees) is the most important preventive measure to reduce the attacks [24,31,32]. In cases of severe damage, the use of trap trees could lead to rapid suppression of bark beetle outbreaks in managed forests [3,33,34].
The data of the Information system of EFA identify Tomicus piniperda as a major pest in the Kirkovo region (causing about 80% of all damage in pine plantations). However, other studies in Bulgaria have pointed Ips acuminatus as the major destructive pest, especially in Pinus sylvestris stands [3,19].
The higher degree of defoliation in P. nigra plantations compared to P. sylvestris was due to the damage caused by Thaumetopoea pityocampa, which is the most dangerous defoliating pest in pine forests in Bulgaria [35,36].
In assessing the health status of forest stands, NDVI thresholds for healthy trees in the range of 0.55-0.70 are usually accepted [11,28,37]. In this study, the average NDVI values for living trees (0.617), dead trees (0.343), grass and shrub vegetation (0.436) and stones and rocks (0.329) were obtained in December, when no grass vegetation usually grows. However, during the growing season, higher NDVI values should be expected, mainly in the categories of grasses and shrubs and dead standing trees.
The present study included several ways in which bark beetle outbreak spots had been determined, examined and compared after a natural disturbance. This combined approach of using remote sensing data and field verification is extremely suitable for studying the attacks of a wide range of bark beetles, including Ips typographus, outbreaks of which have increased in recent decades in spruce-dominated mountain forests across Europe [38][39][40][41][42].

Conclusions
The NDVI values for different categories of land cover (living trees, lying and standing dead trees, grass and shrub vegetation, stones and rocks) established in this study may be used for accurate mapping of bark beetle spots not only in pine plantations but also in other forest stands.
The use of a complex methodological approach through the application of silvicultural and preventive forest pathological activities in the watersheds, with predominant participation of degraded lands, will lead to the improvement of the health condition of pine plantations. In this regard, the control of Ips and Tomicus species will contribute to the overall improvement of the health of the plantations and reduce the risk of attacks by other species of xylophages.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon request.