Modelling within-region spatiotemporal variability in grapevine phenology with high resolution temperature data

Aim: While it is well recognised that climate can vary spatially so that grapevine phenology significantly varies within a winegrowing region, the magnitude of temporal variability in climate determining factors to engender a significant difference in phenology over time can only currently be described in general terms. This study’s aim was to quantify and compare the temporal variability, and its interaction with spatial variability, in grapevine phenology and three derived viticulture climate indexes for three viticulture regions in eastern Australia with varying topographies, latitudes and continentalities. Methods and results : Maximum and minimum temperature data were spatially interpolated to produce fine scale topoclimate maps at 30 m resolution for every day of a 20 year period from 1998 to 2018 for the three viticulture regions. Grapevine phenology modelling using a heat accumulation methodology was then applied to estimate growing season temperature, cool night index and post-harvest growing degree days (GDD) for every 30 m map pixel of the three regions in each of the 20 growing seasons. Conclusions : Summary statistics that quantify the spatial and year-to-year temporal variability and the interaction of spatial and temporal variability (i. e. spatiotemporal variability) demonstrated significant differences in grapevine phenology between each of three regions. The key conclusion was that within-region temporal variability in growing season temperature, cool night index and post-harvest GDD exceeded within-region spatial variability. Significance and impact of the study: This is the first study to quantify temporal variability in modelled grapevine phenology, compare this with the level of spatial variability and also consider the interaction of spatial and temporal variability in grapevine phenology within viticultural regions. This study opens many potential lines of further investigation into the effect of temporal variability on wine production and the cultural terroir factors that develop within a wine region as a result, particularly in the context of forecast future increases in inter-annual temperature variability.


INTRODUCTION
Economically sustainable production of winegrapes depends on suitable terroir, i. e. site characteristics combined with tailored viticultural practices (Wilson, 1998).Climate is arguably one of the most important primary components of terroir (van Leeuwen et al., 2004 ;Winkler et al., 1974), with the temperature component of climate, in particular, having an important role in determining regional winegrowing characteristics.Winegrowing regions are generally contiguous geographic areas that have common terroir elements, which combine to produce distinct wine styles that often typify a region.A cultural element of terroir is selection of grapevine varieties more suited to the specific temperature regime of a vineyard site.Comparisons can be made with established regions with a similar climate to select appropriate varieties for new vineyard plantings based on historical temperature data (Jones, 2006a).
Care needs to be taken when considering published historical climate data of a viticultural region, because within-region spatial variability in temperature can be significant.As a result, the climate data for a particular point within a region (collected at a meteorological station, for example) may not be representative of the climate at another location within the region due to climatic variation.Spatial variability of climate within regions has been characterised for Western United States AVAs (American Viticultural Areas ; Jones et al., 2010), Australian GIs (Geographic Indications ; Hall and Jones, 2010), and wine regions in Portugal (Fraga et al., 2014), France (Bois et al., 2018) and New Zealand (Anderson et al., 2012).Many regions span wide ranges of temperatures and climate categories, such as Winkler's (1974) temperature zones.Within-region spatial variability is often significant, and the level of spatial variability tends to correlate with the degree of spatial variability in topography, primarily elevation.
While spatial variability has been characterised in summaries of viticultural indexes used to compare growing conditions between regions, temporal variability in winegrape conditions has not.Climate is commonly characterised by averaging climate data over a 30 year period.In the above-mentioned studies characterising spatial variability in temperature related viticultural growing conditions, temperature data over the 30 year period is averaged to produce daily data that is the average maximum and minimum temperatures of the same Julian day in each of the 30 years.Grapevine phenology modelling then proceeds with this averaged dataset.Any information about temporal variability in temperature across the period is therefore lost from the dataset.The question posed therefore, is what is the level of temporal variability in grapevine phenology of a region?Instead of using averaged climate data and modelling this once, the current study ascertained temporal variability in production by modelling each year in the period individually, and summarising the key temperature derived viticulture indexes across the years.Three viticultural indexes were selected for comparison : growing season temperature (GST), cool night index (CI) and post-harvest growing degree days (post-harvest GDD), with the aim of describing fine scale spatial variability using topoclimate modelling and inter-annual temporal variability by modelling each of the growing seasons individually.The three selected indexes are useful in terms of providing key viticultural comparisons spatially and temporally, and are also likely to be important for describing associated viticultural changes associated within a warming climate, with a view to making future comparisons with global climate model temperature forecasts.Growing season temperature is the mean temperature (of both maximum and minimum temperatures) across a growing season.It is an easily understood metric that can be used to compare regional differences in viticultural suitability and is strongly correlated with annual totals of heat accumulation metrics (Hall and Jones, 2010) including Winkler 's growing degree days (Winkler et al., 1974) and the Huglin Index (Huglin, 1978).The cool night index (CI) (Tonietto and Carbonneau, 1999) is the mean minimum temperature of the ripening period month and is considered important in influencing grape colour (Kliewer and Torres, 1972 ;Mori et al., 2005).Post-harvest GDD (Hall et al., 2016) is the total growing degree days available to grapevines after harvest and is calculated for the period from harvest to the end of the year (in the Southern Hemisphere, the end of the year is taken to be 30 June).Accumulation of carbohydrates stored in the grapevine trunk and roots during the post-harvest period affects the following growing season's bloom and yield, and Andrew Hall and John Blackman covaries significantly with climate (Holzapfel et al., 2010).
Metrics to describe the interaction of spatial and temporal variability, i. e. the spatiotemporal variability, were produced, compared between regions, and discussed in terms of their use in providing a more holistic description of climate in a viticultural context.Consideration of spatiotemporal climatic variability may add value to future analyses of climate change effects on sustainable viticultural production.

Case study region and cultivar
Australia is a major producer of table wine with national crushes of 1.93 million and 1.79 million tonnes in 2017and 2018, respectively (Wine Australia, 2018).Wine production is conducted across a wide range of climates with major viticultural regions located across all five recognised wine production climate types defined by a commonly applied temperature classification system formulated by Winkler et al. (1974) (Hall and Jones, 2010).Three Geographic Indications (GIs, or wine regions) as defined by Wine Australia (Australian Wine and Brandy Corporation, 2008) were included in this study : Orange, Canberra District and Yarra Valley.All three are located in eastern Australia.Their latitudinal separation is climatically significant, with the centre location of Orange at 33.3°S, Canberra District at 35.03°S and Yarra Valley at 37.7°S.Orange and Canberra District have similar levels of seasonality whereas Yarra Valley is closer to the coast, resulting in a less seasonal, more maritime influenced climate.Elevation maps of the three regions describe the varying topography of each of the three regions (Figure 1).Canberra District's median elevation is 644 m a.s.l., ranging from 264 to 1419 m ; Orange's median elevation is 835 m a.s.l., ranging from 376 to 1390 m ; and Yarra Valley's median elevation is 251 m a.s.l., ranging from 17 to 1338 m.The differences in median growing season temperatures of the three regions are largely a product of their differences in latitude, elevation and continentality, with the spatial FIGURE 1. Topographic maps of the three study regions Digital terrain maps of Canberra District, Orange and Yarra Valley Geographic Indications that include a 10 km buffer beyond each GI's boundary.range of growing season temperatures (mean temperature from 1 October to 30 April) provided in Table 1.With the middle elevation, latitude and continentality, the Canberra District's median GST for the period 1971 to 2000 was 17.7 °C.Orange has a higher mean elevation than Canberra District but its more equatorial latitude produced a similar median GST of 17.6 °C.The Yarra Valley has the lowest mean elevation of the three study regions, but was the coolest as it is the most poleward and most coastal (with consequent cooler summers) ; its median GST was the lowest, at 16.6°C.Syrah (also known as Shiraz) accounts for the greatest proportion of any one variety to Australia's national crush, being 23.9 % in 2018 (Wine Australia, 2018).The Syrah cultivar has midrange heat accumulation requirements in order to initiate budbreak (Moncur et al., 1989) and produce mature fruit (Gladstones, 1992), making it a suitable representative variety on which to base modelling.Syrah is described as being suitable to regions that have average growing season temperatures in the range 16 -19 °C (Jones, 2006b), which overlaps significantly with the temporally averaged climate data of all three regions presented in Table 1.

Phenological modelling
Spatially interpolated daily maximum and minimum temperature data were acquired in raster format (i.e. maps made up of small individual units (pixels), each having a numerical value describing the temperature for the coincident actual location stored as 2D map arrays) for each day from July 1 1998 to June 30 2018 from the Australian Bureau of Meteorology (BOM) (2010) for a total of 15336 temperature datasets.These data cover the whole of Australia with a spatial resolution of 0.05 decimal degrees (approximately 5 km at 35°S).The analysis proceeded using only the pixels within each wine region boundary, which were extracted from each temperature dataset.All data processing was conducted with R version 3.5.1 (R Core Team, 2018), with the function packages, sp (Pebesma and Bivand, 2005), raster (Hijmans and van Etten, 2012) and rgdal (Bivand et al., 2019).Due to the volume of data calculations, parallel processing techniques were used on a 48 core minicomputer with 1 Tb RAM.The project generated approximately 12 Tb of intermediate data.
To produce topoclimate maps, temperatures were estimated for each pixel in a high spatial resolution digital terrain model (DTM) using equations based on locally calculated environmental lapse rates (the change in temperature with respect to change in elevation).A high-quality 1-second resolution (approximately, 30 m at 35° latitude) DTM derived from Shuttle Radar Topography Mission (Gallant et al., 2011) was acquired as the base topographic dataset.The pixels within the wine region boundaries of each region were extracted, and equivalent low resolution versions of the DTMs were produced by aggregating pixels to the same resolution and locations as the lower resolution temperature datasets.Environmental lapse rates are known to vary spatially and temporally and differ for maximum and minimum temperatures.Therefore, for each day and each region, environmental lapse rates were calculated for both minimum and maximum temperature datasets.Equations describing the ELRs for each of the 15336 temperature datasets were applied to the DTMs, thereby producing maximum and minimum temperature maps (or topoclimate maps) for each of the three regions for each day in the 20 year period at the high spatial resolution of the DTMs.
The daily topoclimate data were then used as input to model the phenology of Syrah grapevines.Key phenological dates and modified viticultural indexes based on those dates were calculated for every pixel and for each of the 20 years.Phenological modelling followed the method set out in Hall et al. (2016) (Figure 2).For a comparison of grapevine phenology modelling methods see Fila et al. (2014).To estimate budbreak date, GDD 4.5 was summed one day at a time, starting at 1 July, until a 450 threshold (Moncur et al., 1989) was met, where the heat accumulation on each day was calculated using the threshold, θ = 4.5 °C, as follows. ( where T max is the day's recorded maximum temperature and T min is the day's recorded minimum temperature.To estimate harvest date, starting on the budbreak date, GDD 10 (where θ = 10 °C) was summed one day at a time until the target heat accumulation of 1214 C°days was reached.
The 1214 C°days target differs from the 1407 C°days target used in a previous study modelling Syrah used by Hall and Jones (2010), who employed a heat accumulation target for Syrah suggested by Gladstones (1992).The 1407 C°days target is equivalent to a GST of approximately 17 °C, which is a more appropriate target for modelling Syrah phenology across a wide range of climates.Since the regions used for this study are relatively cool, the lower end of the range of GST temperatures suggested by Jones (2006a), i. e. 16 °C, was selected as the basis for the heat accumulation target used in this study.A GST of 16 °C was converted to a GDD heat accumulation target of 1214 C°days using the conversion equation in Hall and Jones (2010).
OENO One, 2019, vol., x For each pixel, the day on which the target heat accumulation is reached was calculated separately to produce a continuous map of estimated modelled harvest date.(Tonietto and Carbonneau, 2004).

RESULTS
The  shown in the figures conform to the relatively high elevation areas of each of the study regions.
All three indexes exhibited significant betweenregion spatial variation, as well as within-region spatial variation.In addition, temporal variability in all indexes was significant for all three regions.It is clear from Figures 3 -5 that the level of temporal variability exceeded that of the spatial variability.The colour variation within individual maps, indicating spatial variability, is generally less than the colour variation between maps, which indicates higher temporal variability.Temporal variability most greatly exceeded spatial variability for all three indexes in the Yarra Valley (Figure 5).The more topographically varied higher elevation area of the Yarra Valley (in the east of the region) was excluded from the results, because in those areas the heat accumulation target was not met in all years.The remaining Yarra Valley area is relatively flat and therefore climatically homogeneous, which resulted in low spatial variability in derived temperature climate derived viticultural indexes.As the coolest of the three regions, year-to-year variability in the viticultural indexes was relatively high, because greater variability in the timing of the estimated maturity date resulted in more temporal variability in the timing of the periods for which the indexes were calculated, particularly those two indexes that are most sensitive to changes in maturity date timing, i. e. CI and post-harvest GDD.
The interaction between spatial and temporal variability, i. e. the spatiotemporal variability, can be described by comparing the differences in the standard deviations and coefficients of variation of the indexes calculated for the different temporal percentiles (Table 2).Large differences in the coefficient of variation between the 5 th and 95 th temporal percentiles indicate that a region's viticultural index exhibits higher levels of spatiotemporal variability, i. e. the level of temporal variability in the spatial variability of an index.The mean pixel values of the maps describing the 5 th and 95 th temporal percentiles of the growing season temperature differed by 2.2 °C, 2.7 °C and 2.0 °C for Canberra District, Orange and Yarra Valley, respectively.Differences of at least 2 °C in GST are enough to change between categories of GST.In addition, spatial variability in GST varied across the maps produced for the different temporal percentiles in the Orange and Yarra Valley regions.The coefficient of variation of GST decreased steadily for these two regions as the temporal percentile increased ; for Orange, at the 5 th temporal percentile CV = 0.0295, compared to CV = 0.0092 at the 95 th temporal percentile ; and for Yarra Valley, at the 5 th temporal percentile CV = 0.0097 compared to CV = 0.0046 at the 95 th temporal percentile.In general, therefore, in relatively cooler years, spatial variability in GST was greater.The Canberra District, however, did not exhibit a similar change in spatial variability, with the coefficient of variation staying similar at each of the five temporal percentiles.Differences in overall regional topographic characteristics may be an explanatory factor for these differences in spatiotemporal variability between regions, with the Canberra region spanning more heterogeneous terrain with subsequently greater levels of spatial variability in output index maps (Figure 3) than both Orange (Figure 4) and Yarra Valley (Figure 5).The spatial mean value of the temporal median cool night index varied more between regions than GST, being 11.5 °C for Canberra District, 14.5 °C for Orange and 13.5 °C for Yarra Valley.
The range of mean pixel values of the maps describing the 5 th and 95 th temporal percentiles of the cool night index differed by 3.7 °C, 4.6 °C and 4.1 °C for Canberra District, Orange and Yarra Valley, respectively.The cool night index has a greater level of temporal variability than GST, because in cooler years, the timing of the modelled ripening period occurred later.The cool night index therefore exaggerates differences in inter-annual mean temperatures, because, in cooler years, the ripening period was delayed to a cooler part of the year.The cool night index exhibited similar spatiotemporal patterns to that shown for growing season temperature.The coefficient of variation of CI decreased steadily for two regions as the temporal percentile increased ; for Orange, at the 5 th temporal percentile CV = 0.0120 compared to CV = 0.0046 at the 95 th temporal percentile ; and for Yarra Valley, at the 5 th temporal percentile CV = 0.0062 compared to CV = 0.0017 at the 95 th temporal percentile.In general, therefore, in relatively cooler years, spatial variability in GST was greater.The Canberra District exhibited a smaller change in spatial variability across the temporal percentiles, CV = 0.058 at the 5 th percentile and CV = 0.046 at the 95 th percentile, indicating that the level of spatial variability was more invariant in response to temporal variability.
Of the three viticultural indexes considered, the spatial mean value of the temporal median postharvest GDD varied the least in terms of difference between regions, being 526 C°days for Canberra District, 578 C°days for Orange and 431 C°days for Yarra Valley.Both Canberra District and Orange exhibited similar levels of temporal variability with the 5 th temporal percentile being 313 and 327 C°days, respectively, and the 95 th temporal percentile being 751 and 799 C°days, respectively, with Yarra Valley having a similar range but at a lower overall magnitude, 221 C°days at the 5 th temporal percentile and 669 C°days at the 95 th temporal percentile.The standard deviation (spatial variability) of the maps of post-harvest GDD increased only at a small rate with respect to increasing temporal percentile, but since the mean values increased rapidly with increasing temporal percentile, the spatiotemporal variability (described by the difference in the coefficient of variation) was significant.Effectively, spatial variability in temporal variability of post-harvest GDD was much greater in cooler years than in warmer years (CV = 0.528, 0.488 and 0.329 at the 5 th temporal percentile for Canberra District, Orange and Yarra Valley, respectively, compared to CV = 0.235, 0.206 and 0.136 at the 95 th temporal percentile, respectively).Thus, of the three indexes, post-harvest GDD exhibited the greatest level of spatiotemporal variability.

DISCUSSION
The results of this study have shown that significant levels of temporal variability exist for estimated viticultural indexes based on phenological modelling using heat unit accumulation in three case study viticulture regions in eastern Australia.Temporal variability exceeded spatial variability in the selected viticultural indexes in all three case study regions.The level of temporal variability may be related to relative differences in overall average temperatures of each of the regions.The relatively cooler and less continental (and hence lower seasonality) area of the Yarra Valley had the greatest level of temporal variability in growing season temperature, cool night index and post-harvest heat accumulation, with the key phenological event of maturity modelled to take place in a later (cooler) part of the year than the other two regions.Coupled with the relatively lower level of topographic variability of the area (that consistently met the target heat  A shortening of the ripening period to a hotter period of the year compared to the GST of the seven months from 1 October (that includes the cooler period after maturity up until 30 April) has similar explanatory weight to an overall warming climate.Comparing relative regional differences in GST between regions is more informative.
There is a clear difference between the two studies in the relative difference between regions, with the Orange GST being 1.5 °C warmer than Canberra District in the current study, compared to being 0.1°C cooler than the Canberra District for the whole of region, 1-October to 30-April period for 1971-2000 data of Hall and Jones (2010).This difference may be partially explained by the spatial pattern of variability in rates of projected warming as illustrated by Hall and Jones (2009), whereby inland regions of Australia are expected to warm more rapidly than areas closer to the climate modulating influence of the Southern Ocean.Investigation of a greater number of regions across more varied climate zones is required to more fully investigate the effect of spatial variability in projected rates of climate change on spatiotemporal changes to grapevine phenology.
An analysis of observed change over time in temperature due to global warming and its subsequent effect on modelled phenology is of great interest (e.g.Fraga et al., 2016 ;Moriondo et al., 2013 ;Ramos and Jones, 2018).Further analysis of the results produced in this study along with extending the period of modelling to produce comparably modelled outputs would enable trend analysis to inform climate change impacts.In some relatively cool wine regions with established grapevine varieties, a shortening of the period between phenological events may affect the ability of grapes to fully develop desirable characteristics (sugar/acid levels, flavour profiles) necessary to produce high-quality wine (Barnuud et al., 2014 ;Mullins et al., 1992).In addition to extending the study to consistently model viticultural indexes on a year-by-year basis over a greater number of years, incorporating a temporal variability element to assess the effect of warming temperatures on grapevine phenology would provide a more holistic assessment of viticultural implications of climate change.However, it is unlikely that the phenological predictions will consistently match the actual observed phenology of the grapevines in any one year.Firstly, modelling only provides estimates of phenological events for which a degree of error should be expected, particularly during cooler periods when phenological progression can be slow.Modelled budbreak date estimates, for example, have been shown to have a RMSE of 5 -10 days (de Cortázar-Atauri et al., 2009 ;Ferguson et al., 2014).Secondly, cultural factors to even out production lead to vineyard management designed to mitigate the effects of climatic variability (e.g.Zheng et al., 2017), which can include management of inter-seasonal carbohydrate reserves (Holzapfel et al., 2010).In addition, individual weather events that are not captured in this type of analysis, such as latespring frost, drought, heavy rain or hail events, and heat events (Greer and Weedon, 2013) in particular, may cause significant deviation from modelled predictions.Improvements to the phenological modelling procedure may improve the overall accuracy of the results using data describing actual phenological events to validate selected model parameters.
The results of this study suggest that warmer growing conditions may lead to a reduction in temporal variability.With the maturity date occurring earlier, when phenological progression is more consistent, there is less year-to-year variability in the periods over which the indexes are calculated.With less variability in the key phenological date of maturity, ripening period temperature and post-harvest GDD are more likely to occur in similar periods of the year.More variability in the timing of harvest, on the other hand, leads to an increase in year-to-year variability in the period of the year over which CI and post-harvest GDD are calculated with a consequent greater level of year-to-year variability in temperature.
results of this study take the form of a suite of fine scale maps of viticultural climate indexes for each of the 20 years of the study period.Maps describing the level of temporal variability across the 20 years for the key indexes of growing season temperature, cool night index and post-harvest heat accumulation are presented in Figure 3 (Canberra District), Figure 4 (Orange) and Figure 5 (Yarra Valley).Temporal variability is conveyed within these figures by using separate maps for the 5 th

FIGURE 2 .
FIGURE 2. Phenological modellingExample of estimation of phenological events (budbreak and harvest), and heat unit accumulation in the Southern Hemisphere.

FIGURE 4 .
FIGURE 4. Temporal variability in maps of viticultural indexes for the Orange GI Maps of Growing Season Temperature (a -c), Cool Night Index (d -f) and Post-harvest GDD (g -i) illustrating temporal variability via the 5 th temporal percentile (a, d and g), 50 th temporal percentile or median (b, e and h) and 95 th temporal percentile (c, f and i) of the distributions of the indexes calculated from 1998 to 2018.

Figure 5 .
Figure 5. Temporal variability in maps of viticultural indexes for the Yarra Valley Maps of Growing Season Temperature (a -c), Cool Night Index (d -f) and Post-harvest GDD (g -i) illustrating temporal variability via the 5 th temporal percentile (a, d and g), 50 th temporal percentile or median (b, e and h) and 95 th temporal percentile (c, f and i) of the distributions of the indexes calculated from 1998 to 2018.
Distribution percentiles are provided to describe variation in calculated indexes across each of the three districts, with 5 % indicating the 1 in 20 coldest years, 25 % indicating the 1 in 4 coldest years, 50 % indicating the median temperature, 75 % indicating the 1 in 4 warmest years and 95 % representing the 1 in 20 warmest years.The mean is the spatially averaged pixel values across each region, SD is the standard deviation of the pixel values across each region, and CV is the coefficient of variation of the pixel values across each region, provided for each of the temporal percentiles.

TABLE 1 .
Topographic and growing Season Temperature data of the three study regions Data are spatial summaries of each entire region for averaged climate data covering 1971-2000 (after Hall and Jones, 2010).

TABLE 2 .
Temporal distribution of viticultural phenological indexes