Aeration to Manage Insects in Wheat Stored in the Balkan Peninsula: Computer Simulations Using Historical Weather Data

Wheat is one of the major crops throughout the Balkan peninsula of Europe. Specific harvest and binning dates can vary depending on the specific geographic region. Grain aeration, wherein ambient air is used at low-volume airflow rates to cool a grain mass to levels that will suppress insect population development, is an under-utilized component of pest management plans for stored wheat. The successful use of aeration can potentially reduce fumigation of stored wheat, which will contribute to the amelioration of increasingly prevalent phosphine resistance. Historical weather data were used from 19 sites in the Balkan region to predict how quickly grains could be cooled through the use of aeration, using a web-based aeration model, and three different starting dates, including 1, 15, and 30 July. The model was used to predict population growth and development of Sitophilus oryzae, the rice weevil, with and without the use of aeration. Results show that, in the northern regions of the Balkans, aeration implemented at the start of binning reduced insect populations far below pest levels in unaerated wheat, and may potentially eliminate the need for fumigations. In more southerly regions, additional chemical inputs, such as fumigation or grain protectants, may be necessary in conjunction with aeration. Results provide guidelines for the increased potential of using aeration for the management of wheat produced and stored in the Balkan peninsula.


Introduction
Aeration is used as part of many integrated pest management plans for stored grains in the temperate regions of North and South America, continental Europe and Asia, and Australia [1][2][3][4][5]. The practice of using aeration is generally defined as utilizing ambient air to cool a grain mass to temperatures that will not support insect population growth and development [6]. Most stored product insects do not develop at temperatures below 15 • C [7,8], and this temperature is often a target for aeration management [9,10]. Typical aeration airflow rates in English and metric units are about 0.1 to 0.3 ft 3 /minute (min)/bushel, or 0.12 to 0.36 m 3 /min/metric ton (MT), respectively [9]. These airflow rates are independent of the size of storage bin or silo, and assume the aeration equipment is properly sized Agronomy 2020, 10, 1927 2 of 16 to the bin so that the desired airflow rate can be delivered throughout the grain mass. Aeration should not be confused with grain drying, which is done to remove excess moisture from grains after harvest, and involves airflow rates several orders of magnitude greater than those used for aeration [11].
In most temperate climatic regions, hard red winter wheat is planted in autumn, then harvested and stored in late spring and summer, depending on the specific region and climate. The ambient temperature at the time wheat is harvested and stored can be conducive to rapid insect pest population growth. Hence, pest management practices should be followed to prevent economic loss [12]. The most optimal method of utilizing aeration is through the use of a controlling device that is wired into the electrical system of a storage bin, and is set to operate when temperatures fall below a specified threshold [13,14]. Theoretical engineering calculations state that, at an airflow rate of 0.12 m 3 /min/MT, 120 h below a specified temperature are required to cool a storage bin or silo to that temperature [9]. A more complicated system utilizes a controller that operates whenever the outside ambient temperature is below the temperature at a sensing point inside the bin, and often involves the concept of "wet bulb" temperatures [15,16].
When assessing the potential of using aeration in a given geographic region, temperatures after harvest and binning must be cool enough to provide the necessary hours to cool a storage bin. Historical weather data can be used to predict how quickly wheat could be cooled using aeration, at different airflow rates and aeration start dates [17][18][19]. Aeration can also be coupled with insect growth models to predict insect population growth with and without the use of aeration [20][21][22]. Scientists at the USDA-ARS-Center for Grain and Animal Health Research (CGAHR), in Manhattan, KS, USA, developed an Expert System called Stored Grain Advisor that would predict how quickly wheat could be cooled with aeration, and was coupled with an insect population model that would predict growth of either the lesser grain borer, Rhyzopertha dominica (F.), or the rice weevil, Sitophilus oryzae (L.) [23,24]. The last operating system that this program would run on was Windows XP, so this Expert System is now obsolete. Scientists at Texas A&M AgriLife Research Center in Beaumont, TX, USA, developed an easy to use web-based system for aeration management of stored rough rice [25]. This system has been utilized to predict how quickly stored rice could be cooled in the southern US, and also to predict population growth of R. dominica and S. oryzae in stored rough rice, with and without the use of aeration [25,26]. This is the only web-based program currently available in the scientific literature for use in predictive aeration modeling studies.
Wheat is grown throughout the Balkan peninsula of Europe, which can be roughly defined as the countries of Greece, Albania, the European part of Turkey, Moldova, Bulgaria, Romania, and the countries of the former Yugoslav Republic (Slovenia, Croatia, Serbia, Bosnia, Kosovo, Montenegro, and North Macedonia). This is a large geographic region with varied climatic patterns from south to north. The objectives of this study were to: (1) utilize historical weather data to predict how quickly stored wheat could be cooled throughout the Balkan region at different binning/aeration start dates, and (2) to predict the population development of S. oryzae, with and without the use of aeration.

Collecting Weather Data
Daily historical weather data for selected sites throughout the world are available from the website of the USA National Oceanic and Atmospheric Administration (NOAA) National Climatic Data Center (https://www.ncdc.noaa.gov). The Texas A&M Beaumont website (https://beaumont. tamu.edu/climaticdata/WorldMap.aspx) also has all NOAA data processed and updated twice weekly. Previous methodologies for utilizing daily weather data to predict cooling hours below specified temperature thresholds were based on a written Q-basic computer program that also required daily sunrise and sunset times to estimate hourly temperatures for a given day [17,27]. A newer simpler method for examining weather data is based on a modified sine wave model in Microsoft Excel that estimates hourly weather data using just the temperature, hereafter termed Baskerville-Emin [28]).
Agronomy 2020, 10, 1927 3 of 16 This Baskerville-Emin modeling approach was previously used in degree-day models for field crop pests [29][30][31], and was recently adapted for use in a modeling study to classify wheat stored in Slovenia into different temperature zones for aeration management [32]. A more detailed description of the Baskerville-Emin model is found in Arthur et al. [32].
The USA-NOAA website mentioned above was used to examine data for weather stations in each country in the Balkan region, using data from 2010 to 2019. There were several criteria used to select weather stations: (1) weather station data in the center of large cites were eliminated to avoid urban heat-island effects, (2) stations with extensive missing data were eliminated because of difficulties in interpolating between missing data points or extrapolating beyond dates with missing data, and (3) stations were selected to provide geographic representation throughout the region, thus if two stations were close together in a particular country only one was used. In some countries only one station could be used because there was too much missing data for the 10-year period, or there were limited available sites regardless of data completeness. For example, there were no stations in Albania that could be used in the study, and only one station in Greece and Moldova met the data set requirements for completion of data for the 10-year period. The final total of stations used in the study was 20, depicted in the map in Figure 1. However, even with this filtering process, it was sometimes necessary to adjust the data to cover missing data points within a given year. For example, if the average temperature for a given date was missing, but either the maximum or the minimum temperature was missing, then the missing data point was calculated by subtracting the given maximum or minimum temperature from the average. If the average temperature was given but both the maximum and minimum temperature were missing, a 12 • C (20 • F) difference was assumed between the maximum and minimum temperature. If there was a missing value in either the maximum or minimum temperature, but the temperatures on the surrounding dates were given, missing values were interpolated. The final result was complete yearly data sets for the dates of 1 July to 30 November for each of the listed 20 weather stations.
Agronomy 2020, 10, x FOR PEER REVIEW 3 of 17 in Microsoft Excel that estimates hourly weather data using just the temperature, hereafter termed Baskerville-Emin [28]). This Baskerville-Emin modeling approach was previously used in degreeday models for field crop pests [29][30][31], and was recently adapted for use in a modeling study to classify wheat stored in Slovenia into different temperature zones for aeration management [32]. A more detailed description of the Baskerville-Emin model is found in Arthur et al. [32]. The USA-NOAA website mentioned above was used to examine data for weather stations in each country in the Balkan region, using data from 2010 to 2019. There were several criteria used to select weather stations: (1) weather station data in the center of large cites were eliminated to avoid urban heat-island effects, (2) stations with extensive missing data were eliminated because of difficulties in interpolating between missing data points or extrapolating beyond dates with missing data, and (3) stations were selected to provide geographic representation throughout the region, thus if two stations were close together in a particular country only one was used. In some countries only one station could be used because there was too much missing data for the 10-year period, or there were limited available sites regardless of data completeness. For example, there were no stations in Albania that could be used in the study, and only one station in Greece and Moldova met the data set requirements for completion of data for the 10-year period. The final total of stations used in the study was 20, depicted in the map in Figure 1. However, even with this filtering process, it was sometimes necessary to adjust the data to cover missing data points within a given year. For example, if the average temperature for a given date was missing, but either the maximum or the minimum temperature was missing, then the missing data point was calculated by subtracting the given maximum or minimum temperature from the average. If the average temperature was given but both the maximum and minimum temperature were missing, a 12 °C (20 °F) difference was assumed between the maximum and minimum temperature. If there was a missing value in either the maximum or minimum temperature, but the temperatures on the surrounding dates were given, missing values were interpolated. The final result was complete yearly data sets for the dates of 1 July to 30 November for each of the listed 20 weather stations.

Weather Station Data Summation Analysis
Yearly data for each of the 20 stations for the dates of 1 July to 30 November were averaged over the prior 10 years to obtain a final data set for analysis. The Baskerville-Emin model was then used to predict hours below 15 °C on each date from 1 July to 30 November, and to first obtain the number of hours on each date that were below this temperature threshold in August, September, and October. An initial tally was made based on the number of hours below 15 °C in August and September, and the calendar date and Julian Date by which the threshold of 120 h was accumulated for each station. The hours below 15 °C in August and September, and the Julian Date of 120 h accumulation for the 20 stations, were then plotted on a map of the Balkan region, obtained from the NUTS 2013 dataset for political boundaries on the Geoportal of the European Commission (EuroSTAT), using procedures described in detail in Arthur et al. [32]. Briefly, this procedure involved using QGIS v.3.12.3 to plot the temperature data from the stations onto the map by interpolating between the data for the stations.

Weather Station Data Summation Analysis
Yearly data for each of the 20 stations for the dates of 1 July to 30 November were averaged over the prior 10 years to obtain a final data set for analysis. The Baskerville-Emin model was then used to predict hours below 15 • C on each date from 1 July to 30 November, and to first obtain the number of hours on each date that were below this temperature threshold in August, September, and October. An initial tally was made based on the number of hours below 15 • C in August and September, and the calendar date and Julian Date by which the threshold of 120 h was accumulated for each station. The hours below 15 • C in August and September, and the Julian Date of 120 h accumulation for the 20 stations, were then plotted on a map of the Balkan region, obtained from the NUTS 2013 dataset for political boundaries on the Geoportal of the European Commission (EuroSTAT), using procedures described in detail in Arthur et al. [32]. Briefly, this procedure involved using QGIS v.3.12.3 to plot the temperature data from the stations onto the map by interpolating between the data for the stations.

Modeling Analysis and Statistical Analysis
The average data sets for each station were then sent to cooperators at Texas A&M University, Beaumont, TX, USA for input into the web-based application for estimation of S. oryzae population development with and without the use of aeration. The values for the simulations were as follows: bin diameter, 3.7 m; bin height, 7.4 m; bin headspace, 0.6 m; grain depth, 6.1 m; starting grain temperature, 28.3 • C. Simulations were run using three different aeration start dates for each station, 1, 15, and 30 July. The starting value for S. oryzae was 1.5 mixed sex-adults/MT entering the storage bin on each start date. The data from the temperature threshold study described above were then used to classify the weather stations into different groupings or zones. One station was selected from each zone, from warmest to coolest: aeration cooling patterns on each start date for these selected stations were analyzed to compare the temperature patterns in unaerated wheat. Data from the web-based aeration model output were plotted using Sigma-Plot software (Version 11, Systat, San Jose, CA, USA). A correlation analysis was also done using Correlation Procedure in the Statistical Analysis System (SAS, Version 9.2, Cary, NC, USA) to plot predicted numbers of S. oryzae on each start date with the latitude of individual stations.

Grouping of Weather Station Data
The mapping procedure showed there were few hours below 15 • C in August throughout much of the southern and eastern Balkan region, but cooling could be initiated in the northernmost regions, including the country of Slovenia ( Figure 1A). The threshold of 120 h below 15 • C could be achieved throughout most of the Balkan region in September, except for southern Greece and far western Turkey ( Figure 1B). The Julian dates by which 120 h below 15 • C are shown in Figure 1C, with later dates predicted for the southern and eastern Balkans.
The 20 stations were grouped into four zones based on hours below 15 • C in August and September, and the calendar date and Julian date by which 120 h below 15 • C were accumulated ( Table 1). The sites of Larissa Greece, and Edirne, Turkey, were the two warmest sites (the first zone), with 120-h accumulation dates reached on 8 and 6 October, respectively. The second zone was comprised of three sites in Bulgaria and Romania, with varied hours below 15 • C in August and September but with 120-h accumulation dates reached on 29 and 26 September. The third zone consisted of seven stations with hours below 15 • C ranging from 3 to 26 in August and 175 to 221 in September, and 120-h accumulation dates reached from 20 to 24 September. The fourth and final zone consisted of eight stations with 207 to 359 h below 15 • C in September, and 120-h accumulation dates reached from 30 August to 18 September.

Prediction of S. oryzae Populations
In general, the predicted numbers at each station declined in both unaerated and aerated wheat with advancing start date. However, the predicted impacts of aeration were most apparent in the stations from the warmest temperature zones. The obtained results show a huge difference in predicted populations between unaerated and aerated wheat, at each start date, especially in those stations in zones 1 and 2 ( Table 2). For example, the predicted populations were about 300, 114, and 90 times greater in unaerated versus aerated wheat in Larissa, Greece, with start dates of 1, 15, and 30 July, respectively. Using Vidin, Bulgaria as an example from Zone 2, the predicted populations were about 172, 236, and 133 times greater in unaerated versus aerated wheat with the three successive July start dates.   1 The longitude and latitude of each station is given in Table 1.
For most of the stations from zones 2-4, the impact of aeration on population suppression seemed to be greater with the 15 and 30 July start dates compared to the 1 July start date, primarily because aeration suppressed population development to the level that there was no difference between the predicted starting and ending adult populations. With the exception of Vidin and Chisinau, all stations in zones 1 and 2 had more than 1.5 adults/MT predicted with the latest aeration start date of 30 July. There was no increase in predicted adult populations with the aeration start date of 30 July for all stations in zone 3, as the starting and ending predicted adult population was 1.5 adults/MT. In zone 4, there was no predicted population increase for any station with aeration start dates of 15 or 30 July, Agronomy 2020, 10, 1927 8 of 16 and no increase in some stations with an aeration start date of 1 July. Aeration clearly had an impact on predicted populations in all zones, though the predicted decrease in populations with aeration compared to no aeration was more dramatic in zones 1 and 2.

Correlation Analysis of Predicted Adult Populations
There was little to no correlation between longitude or elevation with respect to adult population levels in unaerated or aerated wheat for any start date (Proc Corr, SAS). The, r and p values for adults in unaerated wheat with respect to longitude ranged from 0.21-0.23, and 0.33-0.36, respectively, while r and p values for adults with respect to elevation ranged from −0.24 to −0.40 and 0.08 to 0.33, respectively. The r and p values for adults in aerated wheat with respect to longitude ranged from 0.16 to 0.38 and 0.09 to 0.50, respectively. The r and p values for adults in aerated wheat with respect to elevation ranged from −0.28 to −0.37 and 0.11 to 0.24, respectively. However, there was a reasonably strong negative correlation of adult populations with latitude in unaerated wheat, with lower populations from south to north with each start date of 1, 15, and 30 July (p = 0.0010, r = −0.678; p = 0.0009, r = −0.68; p = 0.0013, r = −0.684, respectively). There was also a strong negative correlation of adult populations with latitude in aerated wheat with each start date of 1, 15, and 30 July (p = 0.0055, r = −0.597; p = 0.0045, r = −0.607, and p = 0.0029, r = −0.630, respectively).

Aeration Cooling Patterns-Example from Each Zone
Four stations were selected from each zone to plot aeration cooling patterns and predicted populations of eggs, larvae, and adult S. oryzae: Larissa Greece, Vidin Bulgaria, Osijek Croatia, and Murska Sobota Slovenia. Data for Larissa are plotted in Figure 2. The predicted temperature in unaerated wheat never dropped below 15 • C during the entire time between July and 30 November regardless of the start date, while aeration after binning immediately began cooling the wheat with the temperature dropping below 15 • C around late September-early October (Figure 2A-C). The high temperatures in unaerated wheat led to the correspondingly high populations of the different life stages of S. oryzae ( Figure 2D-F). Estimated numbers of eggs, larvae, and adults in unaerated wheat on the earliest start date of 1 July were about 80,000, 60,000, and 18,000, respectively, which were multiple orders of magnitude greater than the predicted numbers of those life stages in aerated wheat ( Figure 2D). The predicted numbers of each life stage in unaerated wheat were much lower with the later start dates ( Figure 2E,F). However, even with aeration, S. oryzae populations were predicted to increase beyond the starting value of 1.5 adults/MT regardless of start date. Regardless, the benefits of immediate aeration were apparent as the predicted population levels of all life stages were greatly reduced compared to the predicted levels in unaerated wheat.
The station selected to represent Zone 2, Vidin Bulgaria, is show in Figure 3. In unaerated wheat, the predicted temperature approached, but did not drop below 15 • C regardless of start date, while the predicted temperature in aerated wheat fell below this threshold in mid-early September, depending on aeration start date ( Figure 3A-C). The maximum numbers of eggs, larvae, and adults in unaerated wheat occurred at the earliest start date of 1 July, then declined with each successive later start date ( Figure 3D-F). Predicted numbers in unaerated wheat were much lower than those for the Larissa Greece station, but there was still extensive population growth at each start date. With aeration, there was predicted development of eggs and larvae at each start date, but at the 30-July start date the numbers of adults were not predicted to increase beyond the initial starting value of 1.5 adults/MT ( Figure 3D-F). Even in this zone, aeration limited the population growth and development of S. oryzae.  The station selected to represent Zone 2, Vidin Bulgaria, is show in Figure 3. In unaerated wheat, the predicted temperature approached, but did not drop below 15 °C regardless of start date, while the predicted temperature in aerated wheat fell below this threshold in mid-early September, for the Larissa Greece station, but there was still extensive population growth at each start date. With aeration, there was predicted development of eggs and larvae at each start date, but at the 30-July start date the numbers of adults were not predicted to increase beyond the initial starting value of 1.5 adults/MT ( Figure 3D-F). Even in this zone, aeration limited the population growth and development of S. oryzae.   The station representing Zone 3, Osijek, Croatia, is shown in Figure 4. In unaerated wheat, the predicted temperature again did not drop below 15 • C regardless of start date, similar to the stations in zones 1 and 2 ( Figure 4A-C). In contrast, the predicted temperature in aerated wheat fell below this threshold in mid-early September, similar to patterns for the Viden, Bulgaria, depending on aeration start date ( Figure 4A-C). In unaerated wheat, there was predicted population growth of all life stages, but numbers declined with each successive start date ( Figure 4D-F). With aeration, there was predicted development of all life stages at each start date, but regardless of start date, adults were not predicted to increase much beyond the initial starting value of 1.5 adults/MT if at all The station representing the coolest zone 4, Murska Sobota, Slovenia, is shown in Figure 5. The predicted temperature in unaerated wheat fell below 15 • C in early November, while in aerated wheat this temperature was achieved in mid-early September with the successive start dates ( Figure 5A-C.). Predicted populations of each life stage in unaerated wheat were much lower than those predicted for the other three stations/zones, but adult populations increased beyond the starting value at each start date ( Figure 5D-F). In aerated wheat, adult populations did not expand beyond the starting value at any start date ( Figure 5D-F

Discussion
The threshold level of 120 h below 15 • C utilized in this study for aeration cooling was chosen to be consistent with previous studies conducted for geographic regions in the US [18][19][20][21]. The results of this simulation study for the Balkan region of Europe clearly demonstrate the benefits of using aeration to cool stored wheat. Analysis of historical weather data showed that in the warmer regions the dates by which 120 h below 15 • C hours were accumulated for utilizing aeration occurred much later in the year compared to the sites in the northern zones. However, coupling the historical weather data with the bin cooling and S. oryzae population model (Yang et al., 2017) provided a method to compare predicted bin cooling and S. oryzae population levels for a geographic region, similar to studies conducted for maize stored in the southern and northern US [20,21]. In the absence of aeration, natural cooling occurred far too slowly to limit S. oryzae population growth, particularly in the warmer zones of the Balkan region.
Our aeration modeling methodology was recently used to classify the country of Slovenia into different climatic regions for wheat storage, with a start date of 30 July for storage and aeration [32]. In this expanded study for the Balkan region, earlier start dates of 1 July and 15 July were used, which indicated that wheat stored and binned early in July would be much more vulnerable to S. oryzae infestation compared to later in July or in early August. This is also consistent with previous studies in the US examining early storage for crops, including rough rice [22].
Our results also showed that in the warmer areas of the Balkan region, aeration alone may not be sufficient to limit S. oryzae population growth, particularly with the early July start dates. The fumigant phosphine is utilized world-wide for management of stored grains, including wheat, and can be used along with aeration in a pest management program. However, multiple studies have shown increasing resistance to phosphine in stored product insects throughout much of the world [33][34][35][36][37], and now including sites within the Balkan region and in continental Europe [38]. The more efficient usage of phosphine, through the development of new management methodologies and extensive monitoring of gas concentrations during a fumigation, is being advocated for mitigating phosphine resistance [39][40][41][42][43]. The fumigant sulfuryl fluoride is registered under the trade name of Profume ® in several countries, and has been successfully evaluated as an alternative to phosphine for fumigation of stored grains [44,45]. Aeration is also compatible with other integrated management techniques, including the use of grain protectants at the time grains are loaded into storage, sanitation and cleaning prior to grain storage, and monitoring to determine the extent of pest populations [46].
While we used historical weather data to come to conclusions about the feasibility of aeration, weather patterns in Europe are quickly changing under climate change [47]. In particular, under a projected intermediate scenario of climate change (e.g., RCP 4.5) for southern Europe, the mean annual temperature is projected to rise by 1.9-2.7 • C by 2100, while there are projected to be 27 more growing days per season. Under the worst-case scenario (e.g., RCP 8.5) in the same region, the mean annual temperature is projected to rise by 3.9-5.4 • C by 2100, while there will be a median of 49 more growing days per season. Though it should be confirmed in future studies, we expect aeration to remain possible in the Balkans under an intermediate scenario of climate change, but it may become much more difficult under the worst-case scenario. Further, extreme weather events are predicted to increase, and to disproportionately disadvantage southern Europe compared to northern Europe [48]. Extreme weather events are more likely to make grain aeration less reliable in controlling stored product insects, especially if they happen in the autumn and tend towards warmer temperatures.
The modeling approach used in this study can easily be adapted to assess the potential for wheat stored in geographic locations world-wide, depending on the availability of historical weather data. However, the weather data on the NOAA website are not complete for many countries, though an individual country may have more data available than are presented on the NOAA website. The lack of a current web-based Expert System for coupling stored product insect growth with bin cooling models specifically developed for wheat is problematic, which can be only be overcome at the moment by using the web-based bin cooling model developed for stored rice. This model is acceptable for S. oryzae, but the development of R. dominica is limited on rough rice [49]. Updating the obsolete Stored Grain Advisor into a web-based Expert System will require computation skills and expertise to adapt it to modern operating systems.

Conclusions
Aeration is an under-utilized component for management plans developed to control insect pests in stored wheat. Historical weather data were used to determine the potential for using aeration to cool wheat in the Balkan region of Europe. The results of modeling simulations indicated that, in more northerly regions, aeration alone may limit the population growth of S. oryzae, but in the southern Balkans, the temperatures are too warm to rely on aeration alone. Additional management using fumigants or grain protectants may be necessary. Regardless of location, aeration dramatically reduced the population development of S. oryzae compared to unaerated wheat, and there is excellent potential for the inclusion of aeration in pest management programs for stored wheat. However, as climate change progresses, aeration may become less reliable during extreme hot weather events, especially if these occur in autumn, but this should be specifically modeled in future work.