Inuence of climatic Parameters on The Twospotted Spider Mite population based on Remote Sensing in Southeast of the Caspian Sea

Tetranychus urticae Koch (Acari: Tetranychidae) is a serious pest in cotton elds worldwide. Monitoring of T. urticae with time-series of vegetation index and climatic factors based satellite data was applied to near real-time assessing. The current study aimed to determine correlations between T. urticae population dynamic and effects of Aerosol index of Sentinel-5, Sentinal-2-NDVI (10m), The Land Surface Temperature (LST), MODIS-Evapotranspiration (ET) and CHIRPS-precipitation. Spider mite out-breaking has coincided with the wheat harvesting and where experienced several dusty days with high aerosol index 0.167. Rainfall had a signicant negative correlation with T. urticae population (R 2 = 0.378) and a threshold precipitation level was estimated at least 2 mm to clean up the canopy. We could not nd a signicant pattern between temperature and T. urticae population until August 2020 and then the signicant positive relationships were observed during August 2020, R² = 0.3519, 0.1283, 0.1675 and 0.178, weekly. Evapotranspiration depicted a statistically synchronous relationship R 2 = 0.637 with T. urticae dynamism. There was a positive correlation between increasing NDVI and T. urticae population until August 2020 and then was shifted to negative pattern R 2 = 0.273 and 0.139. These ndings, aerosol index of sentinel-5 and MODIS-evapotranspiration have potential to forecast spider mite population with high temporal resolution.


Introduction
Tetranychus urticae Koch, is a serious economic pest in cotton elds worldwide because of polyphagous behavior with a wide range of plant hosts about 900 plant species (Brown et al. 2017). The main spider mite fauna, Tetranychus spp., is consisting of Tetranychus paci cus, Tetranychus turkestani and T. urticae, and known as a comment pest of cotton. T. urticae distribution is extended in all territory of Iran (Golestan, Tehran, Azerbaijan, Khorasan, Ardebil, Fars ) on cotton elds (Honarparvar et al. 2012). Yield loss due to dusty spider mite infestation of cotton, Gossypium hirsutum L., was measured by entomological researchers such as during 2010 and 2011 to determine two-spotted spider mite population and injury ratings (Scott et al. 2013). Among abiotic parameters, climatic parameters play a key role in pest population dynamism especially temperature, relative humidity, rainfall, sunshine hours, etc. The climate factors affected not only mites (Ahmed et al. 2012) but also ticks such as a 14-year population study on Scabies in Taiwan (Liu et al. 2016). Five years of systematic sampling program in Bangladesh was studied the effect of weather parameters on Oligonychus coffeae N. (Acarina: Tetranychidae), Red Spider Mite, and showed positive (temperature, relative humidity, and sunshine hours) and negative response (rainfall and cloud coverage) (Ahmed et al. 2012). Demographic parameters of spider mite pests are temperature-dependent such as the intrinsic rate of increase (r m ) of T.
paci cus , Eotetranychus willamettei (McGregor) (Acari: Tetranychidae) (Stavrinides and Mills 2011) as well as development time, sex ratio, and fecundity of T. urticae (Margolies and Wrensch 1996). The lifetables of T. urticae signi cantly was in uenced by different levels of macronutrient N, P and K (Wermelinger et al. 1991). Population dynamics of T. urticae were seasonally investigated under acaricide constraint on eggplant in Bursa Province, Turkey. T. urticae population was positively and negatively responded with mean temperature and mean humidity, respectively (Kumral and Kovanci 2005). According to the targets of the Paris agreement (1.6 °C warming by 2050), until 2050 for tomato, a suitability modeler based on climate change (CC) which was equipped with irrigation facilities (AEI) predicted unsuitable conditions for tomato production and increasing outbreak risk of two-spotted spider mite globally, because of failure in biological control (Litskas et al. 2019). Therefore, climate conditions would also be affected the pest population dynamic, especially in large-scale projections. Using re ectance spectroscopy in two common ways with satellite or unmanned aerial vehicle (UAV) equipped by multispectral imagery. Diffuse re ectance spectroscopy (Visible/Near Infrared Re ectance Spectroscopy) identi ed infestation regions by T. urticae and also quantitatively assessed T. urticae damages in Strawberries (Fraulo et al. 2009). Some researches such as, Reisig and Godfrey (2006), reported that, Aerial and satellite images to distinguish infested cotton by aphid (Aphis gossypii Glover) and spider mite (Tetranychus spp.) from healthy plant. Martin and Latheef, (2017) evaluated a groundbased multispectral optical sensor for detecting spider mite damage in greenhouse condition on cotton production. The supervised classi cation approaches such as Support Vector Machine (SVM) and a transferred Convolutional Neural Network (CNN) was reported for mite-infestation using UAV multispectral imagery (Huang et al. 2018). Species composition are another GIS and remote sensing approaches applied for modeling ecological niche of Tetranychoid mites (Acari: Tetranychoidea) in in different climates of Tehran Province, Iran (Ghasemi Moghadam et al. 2016). The current study aimed to determine the potential effects of ve climate and vegetation characters including air pollution (dust), The Normalized Difference Vegetation Index (NDVI), The Land Surface Temperature (LST), Evapotranspiration (ET) and precipitation on Spider mite population.

2-1-Study area
The study area presents unique climatic conditions of northeast of the Caspian Sea including the Hyrcanian forest (southern), desert (northern), fertile lands and vast paddy elds (central areas) where 14 adjacent districts located in the northeast of Iran (Figure 1), with latitude ranging from 36° 30' to 38° 10' N and longitudes from 53° 50' to 56° 20' E. Climatic diversity of study areas distinguishes the role of climatic effects on hotspot formation of spider mite. The area occupies 21400 km 2 approximately, and the altitude oscillates from -39 to 3780 meters above sea level. Average values of annual temperature and precipitation are 16.88 °C and 454 mm, respectively. In the study area, most of the farmer communities have small landholding and wheat, barley, canola, and broad bean are the most important autumn crops. Soybean, cotton, rice, and sorghum also are the most important summer crops in this province. Golestan province is one of the top three cotton-producing provinces in Iran as it is called the "land of white Gold". Our study has strongly covered all cotton agricultural areas within the Golestan province. The highest area under cotton cultivation belongs to Aqqala (33%), western part of Gonbad-e-Qabous (31%), and Gorgan (10%), while the crop is barely present in Minoodasht and Maravehtapeh at all. To achieve more accurate data, all satellite images were masked by of The Shuttle Radar Topography Mission (SRTM, 30m) (Farr et al. 2007). The mask areas were high-dense forests (Hyrcanian forests) where is not the ecological niche of T. urticae.

2-2-1-Precipitation
The Climate Hazards Group Infrared Precipitation with Stations (CHIRPS) is one of reliable precipitation dataset which provided globally high resolution of precipitation by interpolation approaches (around a 0.05°) for a long period (daily, 1981-present) based on infrared Cold Cloud Duration (CCD), ImageCollection ID on Google Earth Engine (GEE), cloud-based image processing platform, "UCSB-CHG/CHIRPS/DAILY", (Funk et al. 2015).
The Terra-MODIS land surface temperature (LST) product (MOD11A1, 1 km, daily) form 09/06/2020 to 17/09/2020 was manipulated by Google Earth Engine (GEE). For all fteen monitoring windows (Table.1), LST time series provided by mean of LST (day, night) temperature for every window. Net Evapotranspiration was provided by The MOD16A2/V6 product, Evapotranspiration/Latent Heat Flux, is an 8-Day Global 500m. The algorithm was embedding the MOD16 data collection is based on Penman-Monteith equation (Allen 1996), which combination of different input sources including daily meteorological reanalysis data, vegetation property, albedo, and MODIS-land cover.

2-2-4-The Normalized Difference Vegetation Index (NDVI)
Normalized Difference Vegetation Index (NDVI) is one of the most applicable vegetation indexes in remote sensing projection (Lamqadem et al. 2018). NDVI is worked on the red and NIR re ectance of soil.
The NDVI values are in range from −1 (water resources and snow body) to +1 (full vegetation coverage). In the current research, a 7-10-day multi-temporal Sentinal-2 was provided NDVI index during monitoring windows.

2-2-5-Air pollution (aerosol index)
The "NRTI/L3_AER_AI" dataset from Sentinel-5 provides near real-time high-resolution imagery of the UV Aerosol Index (UVAI or AAL). This index is worked according to wavelength-dependent changes in Rayleigh scattering in the UV spectral range for a pair of wavelengths (at the 354 nm and 388 nm). The positive AAI indicates the presence of phenomena which can absorb UV like dust and smoke. The pair of selected wavelengths have been absorbed very low by ozone (Soleimany et al. 2021).
2-3-spider-mite sampling and spatial analysis for T. urticae distribution map (Dependent Factor) 2-3-1-spider mite sampling as ground truth data Ground-truth data of spider-mite population was reported by mean percentage of leaf area which were infected with various types of spider mite symptoms (dust-silk webbing, yellow spots, etc.). 100 cotton leaves were randomly observed by experienced exports and estimated infestation scoring of spider-mite according to equation (1) at each eld.
See formula 1 in the supplementary les.
where n and p were the number of leaves and their infestation percentage, respectively. i, j, l and etc. indicted same infestation percentage. For the unity of procedure, The estimated percentages of infestation were reported based on speci c intervals including 0 for no evidence of spider-mite, 1-5 for 1-10%, 11-20%, 21-30%, 31-50% and 50%≥ infected leaf area, respectively. In this research, collected data from 65 cotton elds were used throughout Golestan Province; 6 (min) * 6 (min) grid cells in the Degrees, Minutes, Seconds (DMS) coordinate system ( Figure 2). These elds were measured the percentages of spider-mite involvement during a xed time (Saturday every monitoring window). After that, distribution maps of spider-mite were spatially analyzed for 15 monitoring windows by Inverse Distance Weighted (IDW) interpolation using ArcMap/ GIS 10.6 software. The deterministic interpolator method by 'Spatial Analyst Tools' in ArcMap software package (version 10.6, ESRI, Redlands, CA) was used to draw spidermite distribution maps from the collected data. This process was performed and tried in every monitoring eld located inside a cell grid.

2-3-2 Correlation analysis:
Various satellite data, such as AAL index of Sentinel-5P, LST-MODIS (Terra), NDVI-Sentinel-2, CHIRPS and ET-MODIS, were assessed as main factors (independent) on spider-mite out breaking. For each window, 200 random points were extracted between an interpolation map of spider-mite (dependent) and ve satellite-based parameters (Independent) by 'ee.FeatureCollection.randomPoints' by Google Earth Engine platform. The relations between ve parameters and the distribution map were statistically examined using ANOVA regression analysis in SPSS (version 23).

2-3-3 Spatial Autocorrelation
Moran's Index has known the more applicable statistic for spatial autocorrelation. Global Moran's I estimate the possibility of spatial correlation at study region. The amplitude of uctuations of the Moran's I values is between 1 and −1. Positive autocorrelation (clustered) and negative autocorrelation (dispersed) in the data translates into positive and negative values of Moran's I, respectively. Random distribution of a variable (no autocorrelation) results in a value close to 0 (Overmars et al. 2003). The relationship between a pixel and its surrounding pixel was estimated by weights matrix. Therefore, a distance-based weight matrix (a threshold distance 5000 m) was applied to consider "neighbors" (just non zero value) for all pixels located at a certain distance. The normal approximation for global Moran's I could be standardized to and (Legendre and Fortin 1989). The signi cance level of was a threshold (1.96) so that the spatial autocorrelation can be consider signi cant if uctuates between 1.96 and -1.96 (Zhang and McGrath 2004). The spatial correlogram shows patterns of spatial autocorrelation when increasing the distance between observations. The spatial correlogram is drawn by two common shapes including Moran's I or (standardized Moran's I) (Legendre and Fortin 1989) plotted in ordinate, against distances (in abscissa). Although, the standardized correlogram represented the spatial correlation distance that is the rst positive peak (Zhang et al. 1998). Local Moran's I is computed to identify the locations of spatial clusters and outliers (Anselin 2010

2-3-4 Geostatistic Method
In the grid system, the classical Inverse Distance Weighted (IDW) interpolation method had the lowest root mean square error (RMSE) in pre-evolution (Gorgan data) than the rest of the methods used to build the abundance maps which is consistent with the results of Al-Kindi et al., (2017). The IDW method estimated values of an unknown pixel by nearby pixels predicted but restricted in the range of maximum and minimum values of true pixels. Nearest Neighbor Statistical (NNS) analysis was spatially detected statistical moth distribution including absence, random, regular, or aggregation possibility in each area (Vinatier et al. 2011).

Results And Discussion
3.1. Spatial Pattern Analysis of Spider-mite using the Spatial Autocorrelation Analysis Generally, the higher Moran's I in absolute value presents the greater the spatial correlation and also, more signi cant spatial autocorrelation shows by the higher standardized form of Moran's I . is able to compare statically spatial patterns of different phenomenon or different calculating parameters of the same phenomenon. At the global autocorrelation, Table 2 depicted the rst three windows (there was not data for spider-mite population at rst window) did not show the signi cance correlation (random distribution) for the standardized form of global Moran's I (≥1.96), and fourth window was at the signi cance level during June 2020. Those periods, the cotton plant did not complete canopy but at fourth window was symmetrical with the wheat harvesting calendar in Golestan province. Wheat harvesting causes huge local dust. After that, a sinusoidal pattern was shown in Global autocorrelation (randomness to aggression, vice versa) ( Table 2). By beginning in August 2020 (nine and eleventh windows), study areas had faced with a spider mite out-breaking and the strongest spatial structure ( Figure.3 h, i). Descending spatial structure at tenth, thirteen and fourteen windows could be related to pesticide application at cotton elds. Figure.S1 represent the standardized spatial correlograms of spider-mite distribution at all monitoring windows (15 windows) and the threshold distance of weight matrix where the Moran's I and the standardized Moran's I, Z(I) reached a maximum value. At rst windows (May 30, 2020 -June 9, 2020), spider-mite was not observed in pilot farms. Therefore, any data was reported by local experts. In generally, the optimal distance was 10 km to reach maximum Moran's I for detecting local spatial pattern.
Positive value of standardized Moran's I values at a distance from 5 km to 15 km indicated spatial clusters of similar spider-mite population at these distance ranges. The interpolation maps of spider-mite population were performed using the Inverse Distance Weighted (IDW) method by the cross-validation of parameters. Evaluation indices from cross-validation of IDW maps for all monitoring windows are given in Table 3 and Figure S2. Cross-Validation of spatial interpolation was estimated model accuracy.
Common parameters which could measure errors are Mean Error (ME), Mean Absolute Error (MAE), Mean Squared Error (MSE) and Root Mean Squared Error (RMSE). RMSE which sensitive to outliers is known as an optimal model evaluator by measuring error size (Willmott 1982;Hernandez-Stefanoni and Ponce-Hernandez 2006). The smaller RMSE showed that semivariogram parameters calculated by tting of experimental values are suitable and geostatistical prediction works more accurate. According to interpolation maps of spider mite distribution, some regions struggled with high population density of T. urticae (Figure 3). The pattern of infestation of T. urticae was the temporal hot spots in central and eastern areas of Golestan province where grains were mechanized harvested with dust pollution (Figure 3 b-f). In late July, by establishing spring and summer crops (such as cotton), the pattern of T. urticae was changed to other regions of study areas. From the eighth window (July 21, 2020-July 29, 2020), Aqqala was showed the main area of T. urticae (Figure 3 g). This nding was con rmed by the Agriculture Administration which named Aqqala is the niche of T. urticae at Golestan province. On fteenth window, T. urticae population has declined sharply because of completion of the growing season of cotton and the arrival of the harvest dates (Figure 3 n). Spatial analysis of T. urticae during growing season of cotton showed ups and downs in T. urticae populations but the climatic aspects which in uencing this pattern was remained hidden. Based on empirical observations, climatic variables have a signi cant impact on T. urticae population, especially Aerosol and dust in air.

Effect of climatic factors on Spider-mite population
As mentioned before, until the fth windows (Figure 3d), there was not enough canopy to establish T. urticae on cotton elds, after that, spider mite out-breaking was coincided with the wheat harvest and high aerosol index 0.167 on 21 June 2020. In the rst three windows, study area where experienced several dusty days, faced with rst T. urticae population peak (Figure 3d and Figure 4). Statically tight relationship between dusty days and spider mite was repeated in transition mode from the ninth and thirteenth to their next windows (w9 to w10 and w13 to w14) (Figure 4). The effect of dust or different environmental conditions on spider mites was reported by several studies (Hodek 1987;Thomas 2001;Guerena and Sullivan 2003). Results of previous studies conducted by Flint (1998) and Guerena and Sullivan (2003) supported current results. For mentioned studies also represented that the dusty conditions almost cause to increase T. urticae population on farms. According to Demirel and Cabuk (2008) nding, spider mite densities were 1.72, 1.75, 4.39 and 2.65 times higher on cotton farms in the vicinity of dirt roads than asphalt roads. Therefore removal of dust population is considering one of applicable approach in organic cotton production (Guerena and Sullivan 2003) to control T. urticae including insecticidal soaps and water washing by complete coverage. Precipitation effect which has the ability to remove dust on cotton canopies on T. urticae population presented in (Figure 5). In rst windows, it could be observed a hidden aspect of rainfall on low population of T. urticae, because of not only low density canopy but also about 2.5 mm/day precipitation at this period. This negative relationship was repeated in the ninth and twelfth windows (approximately 2 mm). Based on current result, even if rain did not fall about 2 mm on July 6 and 16, 2020 you would experience severe outbreaking during July 2020. The threshold precipitation level was estimated at least 2 mm to clean up canopy. Rao et al., (2018)  The absence of precipitation when coupled with suitable temperature was introduced a main contributing environmental factor for the rise in T. urticae population by Rao et al., (2018). Evaluating the relationship between spider-mite population and the means of temperature was another part of the climate studies. The main source of temperature measurement was MODIS-LST imagery twice a day ( Figure 6). We could not nd a signi cant pattern between temperature and T. urticae population until August 2020. The signi cantly tight relationships were observed in the ninth (y = 8.4748x + 38.298; R² = 0.3519; P-value= 0.000), tenth (y = 7.5261x + 38.43; R² = 0.1283; p-value= 0.008), eleventh (y = 3.7942x + 33.272; R² = 0.0859; P-value= 0.041), twelfth (y = 6.6459x + 35.306; R² = 0.1675; P-value= 0.004), thirteenth (y = 6.322x + 36.522; R² = 0.178, P-value= 0.002) windows. In many studies, negative correlated with temperature and spider mite population were reported (Majeed et al. 2016;Fahim and El-Saiedy 2021). The current result is compatible by Fahim and El-Saiedy (2021) who reported a non-signi cant relationship with respect to mean temperature and T. urticae on the beginning of the season. However, there many ndings from previous literatures supporting the positive relationship between the T. urticae population and temperature (Meena et al. 2013;Chauhan and Shukla 2016). Seasonal abundance of T. urticae is in uenced by biotic and abiotic factors. Parasitoid and predators are known as biological factor which are suppressed by temperature and drought (Romo and Tylianakis 2013). Another climate parameter, evapotranspiration was predicted that it affects the abundance, distribution and activity of pests. Increasing evapotranspiration has potential to simulate drought conditions (Mullan et al. 2005). Correlation between spider mite population and evapotranspiration depicted a strong relationship statistically (Mean Square= 0.349; F= 21.038, R 2 = 0.637; p-value= 0.001) in Figure 7. The current result had similarity with nding of Litskas et al., (2019) who reported relationship between evapotranspiration and the T. urticae and its natural enemy, P. persimilis, with R 2 0.46 and 0.60 , respectively. Since that, all climate factors induced their effects to the plant, NDVI is could support a correlation with spider mite population. During monitoring windows, we observed a negative and positive relationship between NDVI and spider mite scoring. There was positive correlation between increasing NDVI and T. urticae population until the tenth window (August 2020) (Figure 8 a-g). The fth and sixth windows (middle of July 2020) showed signi cant relation R 2 =0.107 (p-value= 0.016) and R 2 =0.110 (p-value= 0.015), respectively. Beginning in August, 2020, the type of relationship was shifted to negative especially ninth and thirteenth windows with R 2 =0.273 (p-value= 0.000) and R 2 =0.139 (p-value= 0.006), respectively. This phenomenon was interpreted that by enhancing of cotton canopy, T. urticae had opportunity to establish their communities. After that, negative correlations were due to harmful effect of T. urticae population on severity NDVI decreasing (Figure 8). The NDVI decreasing created by T. urticae was con rmed by numerous studies which assess T. urticae population or their damages (Lan et al. 2013;Martin et al. 2015;Martin and Latheef 2017). During multi-temporal NDVI series was observed severe reduction on the eleventh window (3-11 August, 2020) because of high density of aerosol index and low rainfall (Figure 8 j).

Conclusion
We found that T. urticae had diverse responses to the climatic factors. In a pre-judgment at the beginning of study, among of climatic factors, aerosol index or dusty days was predicted severity affect to increase the spider mite density on monitoring elds. But evapotranspiration was exactly synched with T. urticae dynamic population. Indeed, our ndings aerosol index sentinel-5 and MODIS-evapotranspiration have suitable potential to predict spider mite population with high temporal resolution. Studying these drivers offers a realistic view of what exports design accurate model under a regime of climate change.

Declarations
Funding: The author received no speci c funding for this work Distribution of spider-mite throughout Golestan Province; 6 (min)* 6 (min) grid cells in the DMS coordinate system (yellow points are monitoring eld) Figure 3 Distribution maps of spider mite based on IDW model during monitoring windows, a-n are the sequence windows form June 9, 2020 to September 17, 2020. (First window, May 30 to June 9 was not spider mite population data)

Figure 4
The relationship between UV Aerosol Index extracted from Sentinel-5 imagery and spider mite population (mean score of each window) form June 9, 2020 to September 17, 2020. (First window, May 30 to June 9 was not spider mite distribution data)

Figure 5
The relationship between daily CHIRPS-precipitation and spider mite population (mean score of each window) form June 9, 2020 to September 17, 2020. (First window, May 30 to June 9 was not spider mite distribution data) Page 20/22 Figure 6 The relationship between daily Land Surface Temperature and spider mite population (mean score of each window) form June 9, 2020 to September 17, 2020. (First window, May 30 to June 9 was not spider mite distribution data), … and --are Day and night LST, respectively.