Modeling Phenols, Anthocyanins and Color Intensity of Wine Using Pre-Harvest Sentinel-2 Images

The inclusion of technological innovation and the development of remote sensing tools in wine production are an efficient and productive factor that supports the production and improves the quality of the wine produced. In this study we explored models based on Sentinel-2 image bands and spectral indices to estimate key wine quality variables, such as phenols (TP), anthocyanins (TA) and color intensity (CI), providing different sensory characteristics of wine. Two Cabernet Sauvignon wine harvest seasons were studied, 2017 and 2018, and models with coefficients of determination (R2) higher than 60% were obtained for color intensity and total anthocyanins during the first season, both in a period very close to harvest during the first days of April, so the high periodicity of Sentinel 2 becomes strategic. In addition, homogeneous sectors can be identified in the plots for selective harvesting and thus the winery space can be programmed appropriately. These results suggest further work on the number of samples in order to transform it into a useful tool with the potential to define a differentiated harvest and estimate the accumulation of phenolic compounds and the intensity of wine color, key elements in the final quality of the wine.


Introduction
In a context of increasing competition in international markets, it has become of utmost importance to achieve higher quality standards in the vineyard. This has led to a renewal of viticulture, a revision of agricultural techniques [1] and the inclusion of technological innovation [2] in the production of wine and in the control of the factors that improve its quality. In red wine, tannins and anthocyanins are the most important phenolic classes. Tannins contribute to the mouthfeel of wines and provide the pigments necessary to give red wine its long-term color stability [3,4]. Anthocyanins are directly responsible for the bluish-red color of red grape skins and naturally for the color of red wine. Phenols are responsible for the color, astringency and bitterness of red wine and contribute to the olfactory profile [5]. Since anthocyanins are located in the skin tissue of most grape cultivars, fermentation and maceration (processes in which the skin is used) have an important effect on the concentration of anthocyanin present in the final wine [6]. In red wines, color is one of the main qualitative parameters. On the one hand, it represents the first organoleptic factor perceived by the taster, and on the other hand, high positive correlations have been determined between color and overall wine quality [7,8]. The main sources of red color in wines come from anthocyanins or their additional derivatives that are extracted or formed during the winemaking process [4,9,10]. Wine color not only provides information about possible defects, the type or state of evolution of the wine, but also has an important influence on acceptability [8]. Even the price of wine is assigned not In the process carried out at Viña Montes the harvest is mixed, which means this is done manually in lots that are historically identified as of better quality, and mechanized In the process carried out at Viña Montes the harvest is mixed, which means this is done manually in lots that are historically identified as of better quality, and mechanized in those of lower quality; however, even this differentiation is not strict and may have Remote Sens. 2021, 13, 4951 4 of 15 some modifications. Therefore, this practice is eventually decided during the course of the harvest process. The lots, or predefined zonings of a lot, are evaluated weekly after veraison (the development period in which grapes change color and begin to ripen), in order to measure the sugar and acid content, two of the most common parameters on which winemakers base their decisions on when to harvest [32]. Such evaluations help to control the grape harvest and winemaking processes and thus meet the demands of winemakers [33]. Defining when and which plots will be harvested is based solely on the organoleptic evaluation that is carried out in the field daily, a task that takes between 10 and 30 min per sampling site. Therefore, only one sample is collected per visit, which is expected to be representative of the entire plot.

Vine Cycle Climate of the Study Area
Most of the world's wine production takes place in Mediterranean climates, characterized by high temperatures during the summer, high light intensity and humidity levels that decrease sharply throughout the day [34]. Soil and climate conditions are key for vines, so that the same grape variety, at a similar degree of ripeness but grown in two different areas, may result in two different wines [35][36][37][38]. This is the reason why terroir is important, which we can define as the extension of land whose natural characteristics are a unique set of factors (soil, terrain and climate); therefore, it is the terroir that imprints a distinctive seal onto the wine [35,36]. The study area is dominated by a warm Mediterranean-type temperate climate [38] with a dry season of six months and a rainy winter, with a total annual precipitation between 400 and 600 mm. Additionally, the soils are of alluvial origin, with silt loam textures [39].
The vegetative cycle begins when the buds start growing in early spring (September or October) and generally ends from April to May. The cycle may have some variations that are influenced by differences in temperature and precipitation conditions [40]. An increase in average temperatures and a reduction in thermal oscillation influence varietal aroma and wine color.
The 2017 phenological cycle is marked by the 2016 rainfall, scarce compared to the cumulative rainfall of 2017 (2018 grapevine phenological cycle), which doubled compared to 2016 (Figure 2a). In the 2018 season, the air temperature was relatively lower than in 2017 during the grapevine ripening cycle. The mean temperature difference reached 2 • C, which is equivalent to 11% decrease in that period. Figure 2b shows the comparison of these differences in degrees Celsius with respect to the DOY (day of year).
in those of lower quality; however, even this differentiation is not strict and may have some modifications. Therefore, this practice is eventually decided during the course of the harvest process. The lots, or predefined zonings of a lot, are evaluated weekly after veraison (the development period in which grapes change color and begin to ripen), in order to measure the sugar and acid content, two of the most common parameters on which winemakers base their decisions on when to harvest [32]. Such evaluations help to control the grape harvest and winemaking processes and thus meet the demands of winemakers [33]. Defining when and which plots will be harvested is based solely on the organoleptic evaluation that is carried out in the field daily, a task that takes between 10 and 30 min per sampling site. Therefore, only one sample is collected per visit, which is expected to be representative of the entire plot.

Vine Cycle Climate of the Study Area
Most of the world's wine production takes place in Mediterranean climates, characterized by high temperatures during the summer, high light intensity and humidity levels that decrease sharply throughout the day [34]. Soil and climate conditions are key for vines, so that the same grape variety, at a similar degree of ripeness but grown in two different areas, may result in two different wines [35][36][37][38]. This is the reason why terroir is important, which we can define as the extension of land whose natural characteristics are a unique set of factors (soil, terrain and climate); therefore, it is the terroir that imprints a distinctive seal onto the wine [35,36]. The study area is dominated by a warm Mediterranean-type temperate climate [38] with a dry season of six months and a rainy winter, with a total annual precipitation between 400 and 600 mm. Additionally, the soils are of alluvial origin, with silt loam textures [39].
The vegetative cycle begins when the buds start growing in early spring (September or October) and generally ends from April to May. The cycle may have some variations that are influenced by differences in temperature and precipitation conditions [40]. An increase in average temperatures and a reduction in thermal oscillation influence varietal aroma and wine color.
The 2017 phenological cycle is marked by the 2016 rainfall, scarce compared to the cumulative rainfall of 2017 (2018 grapevine phenological cycle), which doubled compared to 2016 (Figure 2a). In the 2018 season, the air temperature was relatively lower than in 2017 during the grapevine ripening cycle. The mean temperature difference reached 2 °C, which is equivalent to 11% decrease in that period. Figure 2b shows the comparison of these differences in degrees Celsius with respect to the DOY (day of year).

Data Collection
The images used were acquired from the Sentinel 2A and Sentinel 2B platforms and downloaded at processing Level-1C, then atmospherically corrected to obtain surface
The Sentinel-2 constellation is formed by 2 satellites in heliosynchronous orbit with a temporal resolution of 10 days for each satellite, or 5 altogether. Both satellites carry a multispectral image sensor (MSI) able to acquire images in 13 spectral bands, from 433 nm to 2280 nm. The red (665 nm) and near infrared (842 nm) bands are of particular interest for agricultural applications, since they allow for the calculation of several vegetation indices at a 10 m spatial resolution [41]. The productive cycle of the grapevine has several phases: pruning, weeping, budding, flowering, fruit set, veraison, maturity and leaf fall. The budding phase, in which new stems appear and the first leaves appear, occurs in the first fortnight of October; therefore, the date of the beginning of the season was established as 15 October. Figure 3 shows the acquisition date of a set of 8 images used in this work for the two years of study, referred to the day of the season (DOS).
Remote Sens. 2021, 13, x FOR PEER REVIEW 5 of 16 reflectivities, using the Sen2cor application in the SNAP environment (Sen2Cor v2.9 developed by ESA, last access on: http://step.esa.int/main/snap-supportedplugins/sen2cor/sen2cor-v2-9/ (accessed on 20 November 2021)). The Sentinel-2 images were downloaded directly from the Copernicus Open Access Hub platform (https://scihub.copernicus.eu/dhus/#/home (accessed on 30 August 2021)). The Sentinel-2 constellation is formed by 2 satellites in heliosynchronous orbit with a temporal resolution of 10 days for each satellite, or 5 altogether. Both satellites carry a multispectral image sensor (MSI) able to acquire images in 13 spectral bands, from 433 nm to 2280 nm. The red (665 nm) and near infrared (842 nm) bands are of particular interest for agricultural applications, since they allow for the calculation of several vegetation indices at a 10 m spatial resolution [41]. The productive cycle of the grapevine has several phases: pruning, weeping, budding, flowering, fruit set, veraison, maturity and leaf fall. The budding phase, in which new stems appear and the first leaves appear, occurs in the first fortnight of October; therefore, the date of the beginning of the season was established as 15 October. Figure 3 shows the acquisition date of a set of 8 images used in this work for the two years of study, referred to the day of the season (DOS).  To establish the ripening of the plots, a total of 14 grape samples were collected in the field after veraison. The location of the sampling points was established by identifying strategic and homogeneous response areas. Zones that share common characteristics of vegetative expression and topography were selected, and the samples were homogeneously distributed in the plots. A sample was positioned in each plot except for plot 7.08, which presents two sampling positions, which were averaged ( Figure 1). These positions were considered to carry out the micro-vinification in which phenols, anthocyanins and color intensity were measured.

Satellite Data Preparation
Sentinel images were co-registered using ground control points and applying the nearest neighbor resampling method preserving the original pixel size. Eight of the thirteen Sentinel-2 bands were used, corresponding to B3 (Green band-centered on 0.  To establish the ripening of the plots, a total of 14 grape samples were collected in the field after veraison. The location of the sampling points was established by identifying strategic and homogeneous response areas. Zones that share common characteristics of vegetative expression and topography were selected, and the samples were homogeneously distributed in the plots. A sample was positioned in each plot except for plot 7.08, which presents two sampling positions, which were averaged ( Figure 1). These positions were considered to carry out the micro-vinification in which phenols, anthocyanins and color intensity were measured.

Satellite Data Preparation
Sentinel images were co-registered using ground control points and applying the nearest neighbor resampling method preserving the original pixel size. Eight of the thirteen Sentinel-2 bands were used, corresponding to B3 (Green band-centered on 0. The term vegetation index refers to any transformation of two or more spectral bands for assessing vegetation properties either at the leaf or canopy levels, developed for different types of vegetation [42]. NDVI ((NIR − Red)/(NIR + Red)) is a simple indicator closely related to crop vegetative and productive features [30]. GNDVI ((NIR − Green)/(NIR + Green)) is a variant of NDVI more sensitive to water and nitrogen uptake into the plant canopy [18]. In the case of the Chlorophyll Index ((Red − Blue)/Green), the effect of photosynthetic pigments is dominating the reflectance in the visible spectral region [42]. Carotenoids contribute to strong absorption at 445 nm, anthocyanin at 529 nm and chlorophyll at 650 nm [43]. The NDMI ((NIR − SWIR1)/(NIR + SWIR1)) is sensitive to the moisture content in the canopy [44]. The SWIR1 band corresponds to the wavelength comparatively more similar to the Landsat shortwave band, the base of this indicator [45][46][47][48].
The visualization of the behavior and evolution of the vineyard allows to overview the state of the plots and the development in both seasons. For this purpose, the NDVI was used, since this is considered an indicator related to the quantity, quality and state of vegetation [48]. The best vegetative expression, or the maximum potential of NDVI, occurs in summer shortly before the veraison, since just after this, the foliage begins to change, and near the harvest the vines begin to lose their leaves. The development cycle of the vines begins with pruning and new plantings in July/August, and then starts the flowering in August/September, budding in September/October, flowering in October/November, fruit set in November/December, veraison in January/February and finally the harvest in March/April. Therefore, to characterize the entire cycle of each season, the set of NDVI images used was from September to April. Figure 4 shows the NDVI maps displayed on a 10-rank scale. The plots are grouped by year of planting, 2007 and 2010. The 2007 plantation, located north, is made up of 12 plots and the 2010 plantation, located south, by 11 plots. In order to avoid mixed pixels in the boundary of the plots, a buffer of 15 m was applied. The term vegetation index refers to any transformation of two or more spectral bands for assessing vegetation properties either at the leaf or canopy levels, developed for different types of vegetation [42].
NDVI ((NIR − Red)/(NIR + Red)) is a simple indicator closely related to crop vegetative and productive features [30]. GNDVI ((NIR − Green)/(NIR + Green)) is a variant of NDVI more sensitive to water and nitrogen uptake into the plant canopy [18]. In the case of the Chlorophyll Index ((Red − Blue)/Green), the effect of photosynthetic pigments is dominating the reflectance in the visible spectral region [42]. Carotenoids contribute to strong absorption at 445 nm, anthocyanin at 529 nm and chlorophyll at 650 nm [43]. The NDMI ((NIR − SWIR1)/(NIR + SWIR1)) is sensitive to the moisture content in the canopy [44]. The SWIR1 band corresponds to the wavelength comparatively more similar to the Landsat shortwave band, the base of this indicator [45][46][47][48].
The visualization of the behavior and evolution of the vineyard allows to overview the state of the plots and the development in both seasons. For this purpose, the NDVI was used, since this is considered an indicator related to the quantity, quality and state of vegetation [48]. The best vegetative expression, or the maximum potential of NDVI, occurs in summer shortly before the veraison, since just after this, the foliage begins to change, and near the harvest the vines begin to lose their leaves. The development cycle of the vines begins with pruning and new plantings in July/August, and then starts the flowering in August/September, budding in September/October, flowering in October/November, fruit set in November/December, veraison in January/February and finally the harvest in March/April. Therefore, to characterize the entire cycle of each season, the set of NDVI images used was from September to April. Figure 4 shows the NDVI maps displayed on a 10-rank scale. The plots are grouped by year of planting, 2007 and 2010. The 2007 plantation, located north, is made up of 12 plots and the 2010 plantation, located south, by 11 plots. In order to avoid mixed pixels in the boundary of the plots, a buffer of 15 m was applied.

Samples and Micro-Vinification
A Unicam Helios Gamma 9423 1000E spectrophotometer was used to measure both color intensity and phenolic compounds in general. Due to the wide chemical diversity of phenolic compounds, total phenols in must and wines are generally presented in arbitrary units of a phenolic standard, such as the amount of gallic acid needed to produce the same analytical response or gallic acid equivalents [5], as was performed in this case. When

Samples and Micro-Vinification
A Unicam Helios Gamma 9423 1000E spectrophotometer was used to measure both color intensity and phenolic compounds in general. Due to the wide chemical diversity of phenolic compounds, total phenols in must and wines are generally presented in arbitrary units of a phenolic standard, such as the amount of gallic acid needed to produce the same analytical response or gallic acid equivalents [5], as was performed in this case. When measuring the phenolic compounds with spectrophotometers, no specific compounds but Remote Sens. 2021, 13, 4951 7 of 15 rather total compounds are measured, which are expressed on the basis of the specific compounds that recur most frequently among phenols and anthocyanins.
The method for measuring total phenols was measuring absorbance at 280 nm (a gallic acid calibration curve was used). The measurement of the total anthocyanins was carried out by decolorization and measured at 520 nm. The method to measure color intensity was the Glories method, which corresponds to the sum of the absorbances at 420 nm, 520 nm and 620 nm [49].
Micro-vinification was made from the abovementioned 14 field samples. Approximately 1000 kg of grapes per micro-vinification were crushed and blended to enter the alcoholic fermentation process and this ended when 2 g of sugar per liter was reached. The wine was then transferred to barrels. Then came malolactic fermentation, which was carried out in a heated place to increase the speed of fermentation.

Oenological Evaluation
Since phenolic compounds are important in the overall quality of wine, they are the subject of intense research worldwide [6]. They provide sensory characteristics such as flavor, aroma, color, and astringency, among others [31].
Viña Montes evaluates wines in three categories: Limited (LIM), Premiun (PRE) and Superior (SPR). Each of these has three variants with a lower (−) and a higher (+) subcategory, indicating that a taster can evaluate a wine as Premium, but below or above the average of this category.
The factors used to evaluate the wines were as follows: (i) Visual factors: mainly associated with color; therefore, anthocyanins and color intensity will predominate them. Cleanliness was also observed (cloudy, clear, clean, bright, etc.), which is associated to physicochemical stability. Color intensity varies according to the grape variety but seeks to identify colorimetric intensities such as red, intense red, ruby red, violet, intense violet, etc. (ii) Olfactory factors: associated with volatile properties and it is described through aromas (wood, vanilla, fruits, etc.). The greater the number of descriptors and complexity, the higher the classification of a wine. (iii) Flavor factors: tannins, compounds found in the skin and seeds of berries, which produce the sensation of astringency. Not only a high concentration of tannins is sought, but they must be of good sensory quality. Together with the alcohol, they confer the sensation of structure and body. A wine with more body and better structure is classified as superior quality. (iv) Finish of a wine: indicates a higher or lower category; the finish is the persistence of a wine in the mouth, so a wine with a persistent, long, medium or short finish is also a factor associated with categories.
Attending to these factors, each wine sample was evaluated by four winemakers who assigned an individual evaluation. Then, they were averaged to give a final category to each sample. Table 1 shows the final quality category assigned to each plot and the color intensity value, total phenolic compounds and anthocyanins measured with the spectrophotometer. Figure 5 shows the comparison of the quality factors in the two seasons.

Modeling Phenolic Compounds (PC) and Color Intensity (CI)
The modeling of the variables requires the values measured in the micro-vinification and the response of the foliage that will be reflected in the spectral bands and indices. Since the variables are obtained post-harvest, it will be necessary to know the relationship between the samples and the series of images to represent the vegetative expression from veraison to harvest. The average value of pixels considered as representative of the images was taken to the set of pixels within a 15 m radius around the location corresponding to the sampling and micro-vinification (Figure 1).
For each sampling point a table was prepared, including the name of the plot, the values obtained in the microvinification corresponding to total phenols (TP), the total anthocyanins (AT) and the color intensity (CI), which were used as the dependent variables. Similarly, the independent variables correspond to the average value of pixels around the sample in the original Sentinel 2 bands and in the four spectral indices, on the four dates of Remote Sens. 2021, 13, 4951 8 of 15 each season, from January to April. A multiple ordinary linear regression model was used to determine the relationship between the image-derived data and the variables (phenolic compounds and color intensity). Table 1. Post-harvest evaluation of the wine. Final average tasting: color intensity (absorbance); total phenols (mg/L gallic acid equivalent) and total anthocyanins (mg/L malvidin equivalent).

Modeling Phenolic Compounds (PC) and Color Intensity (CI)
The modeling of the variables requires the values measured in the micro-vinification and the response of the foliage that will be reflected in the spectral bands and indices. Since the variables are obtained post-harvest, it will be necessary to know the relationship between the samples and the series of images to represent the vegetative expression from veraison to harvest. The average value of pixels considered as representative of the images was taken to the set of pixels within a 15 m radius around the location corresponding to the sampling and micro-vinification (Figure 1).
For each sampling point a table was prepared, including the name of the plot, the values obtained in the microvinification corresponding to total phenols (TP), the total anthocyanins (AT) and the color intensity (CI), which were used as the dependent variables. Similarly, the independent variables correspond to the average value of pixels around the sample in the original Sentinel 2 bands and in the four spectral indices, on the four dates of each season, from January to April. A multiple ordinary linear regression model was used to determine the relationship between the image-derived data and the variables (phenolic compounds and color intensity).
To select the best model, a descriptive analysis of the variables was carried out, for which the correlation matrix of the bands on the different dates was reviewed. In order to identify the strength of the linear relationship between the variables, as well as the statistical significance of the estimated correlations, a confidence level of 95% was considered; therefore, p-values less than 0.05 will be sought. To select the best model, a descriptive analysis of the variables was carried out, for which the correlation matrix of the bands on the different dates was reviewed. In order to identify the strength of the linear relationship between the variables, as well as the statistical significance of the estimated correlations, a confidence level of 95% was considered; therefore, p-values less than 0.05 will be sought.
To avoid severe multicollinearity between the variables, variance inflation factors (VIF) were considered, which allow measuring the correlation between the predictor variables of the model.
Ordinary multiple linear regression models were used to find the best coefficient of determination (R 2 ) between January and April. To do this, the best combination of variables or the variable that was significant both for the model and individually was selected, followed by the collinearity check of the variables using the variance inflation factor (VIF).

Results
The selected models included only 1 or 2 variables. Some models composed of variables with a high degree of significance and a slight degree of collinearity were retained (e.g., CI February 2017), as well as one of the models with low collinearity that included a variable with a p-value greater than 0.05 (TP April 2017).
During the 2017 season, color intensity presented high R 2 values for the four dates, April being the best of them, as well as the results for total anthocyanins. In both variables, Remote Sens. 2021, 13, 4951 9 of 15 during the 2018 season, the models reached very low R 2 values. The mean absolute error (MAE) corresponding to the average value of the residuals was higher in the 2018 season. The total phenols variable presented more similar errors in both seasons and reached a better R 2 in the 2018 season, February being the best of them. Table 2 shows the summary statistics for each dependent variable for the four dates under study-color intensity (CI), total phenols (TP), total anthocyanins (TA)-including the coefficient of determination (R 2 ), the variables included in each model, their mean absolute error (MAE), the variance inflation factor (VIF) and the significance (p-value) for the model and the individual variables, which constitute all the information necessary to select the month with the best indicators. Considering the small number of samples studied and in order to avoid using variables with high collinearity, those models with VIF values less than 10 were chosen. Figure 6 shows the coefficient of determination obtained at each date and the number of variables included in the model.
Attending to the results obtained, April is the best month to predict the intensity of the coloration and anthocyanins, unlike phenols, whose best R 2 is presented in February.
A summary of the complete statistics of the selected models as well as a scatter plot of the individual models, including the equations of the models in both seasons for the three dependent variables, color intensity (CI), total phenols (TP) and total anthocyanins (TA), are shown in the Figure 7.
Considering the small number of samples studied and in order to avoid using variables with high collinearity, those models with VIF values less than 10 were chosen. Figure 6 shows the coefficient of determination obtained at each date and the number of variables included in the model.
Attending to the results obtained, April is the best month to predict the intensity of the coloration and anthocyanins, unlike phenols, whose best R 2 is presented in February. A summary of the complete statistics of the selected models as well as a scatter plot of the individual models, including the equations of the models in both seasons for the three dependent variables, color intensity (CI), total phenols (TP) and total anthocyanins (TA), are shown in the Figure 7.
As shown in Figure 7, variables are not repeated in the models with respect to the comparison by variable between the two years. The 2017 season presents two models with one variable (CI and TA) and TP with two variables. In the case of the 2018 season, only CI is selected with one variable and TP and TA with two variables. Regarding the coefficients of determination, only the 2017 season shows reliable models above 60%.
Since these models are different depending on the date, their precision is limited in time, which implies that they should be obtained monthly. Although the study brings us closer to a date in which the greatest relationship occurs, the moment could vary slightly in each campaign, due to external factors such as the meteorological conditions. This temporal limitation does not restrict the benefits of this tool in the objective of this study, since these models would serve to identify the zoning where the next sampling or field visit will take place. This would contribute to carry out the harvest differentiation within a plot, allowing to improve the determination of the maturity and optimize the harvest moment in each plot. By applying phenolic compounds and color intensity estimation models, it is possible to obtain representative maps of the quality status at a given date, showing the variability of each plot in terms of the ripening potential represented by the canopy. This would allow, through the use of a limited number of field samples and Sentinel-2 images, to estimate these wine quality variables in the entire plot and, subsequently, to identify ripeness variability and zonation to determine a selective harvesting. Figure 8 shows comparative maps between the NDVI and estimated total anthocyanins. In agreement with Cortell et al. (2007) [50], who determined that vines with low vigor show higher levels of anthocyanin and flavanol in the grape skin, an inverse relationship between the estimated anthocyanins and NDVI can be noticed. As shown in Figure 7, variables are not repeated in the models with respect to the comparison by variable between the two years. The 2017 season presents two models with one variable (CI and TA) and TP with two variables. In the case of the 2018 season, only CI is selected with one variable and TP and TA with two variables. Regarding the coefficients of determination, only the 2017 season shows reliable models above 60%.
Since these models are different depending on the date, their precision is limited in time, which implies that they should be obtained monthly. Although the study brings us closer to a date in which the greatest relationship occurs, the moment could vary slightly in each campaign, due to external factors such as the meteorological conditions. This temporal limitation does not restrict the benefits of this tool in the objective of this study, since these models would serve to identify the zoning where the next sampling or field visit will take place. This would contribute to carry out the harvest differentiation within a plot, allowing to improve the determination of the maturity and optimize the harvest moment in each plot.
By applying phenolic compounds and color intensity estimation models, it is possible to obtain representative maps of the quality status at a given date, showing the variability of each plot in terms of the ripening potential represented by the canopy. This would allow, through the use of a limited number of field samples and Sentinel-2 images, to estimate these wine quality variables in the entire plot and, subsequently, to identify ripeness variability and zonation to determine a selective harvesting. Figure 8 shows comparative maps between the NDVI and estimated total anthocyanins. In agreement with Cortell et al. (2007) [50], who determined that vines with low vigor show higher levels of anthocyanin and flavanol in the grape skin, an inverse relationship between the estimated anthocyanins and NDVI can be noticed. By applying phenolic compounds and color intensity estimation models, it is possible to obtain representative maps of the quality status at a given date, showing the variability of each plot in terms of the ripening potential represented by the canopy. This would allow, through the use of a limited number of field samples and Sentinel-2 images, to estimate these wine quality variables in the entire plot and, subsequently, to identify ripeness variability and zonation to determine a selective harvesting. Figure 8 shows comparative maps between the NDVI and estimated total anthocyanins. In agreement with Cortell et al. (2007) [50], who determined that vines with low vigor show higher levels of anthocyanin and flavanol in the grape skin, an inverse relationship between the estimated anthocyanins and NDVI can be noticed.

Discussion
Defining the quality of wine is relative; however, there are some standard variables to determine the technological maturity (state of greater accumulation of sugar and low acidity), phenolic maturity responsible for the color of the berries [51], sensorial characteristics, such as color, astringency and bitterness, and aging capacity, which are closely related to the perception of wine quality [52]. For this purpose, it is key to consider the variability of a parcel, since the mixture of grapes from different quality indicators do not result in an intermediate quality wine, since some properties found in low-quality grapes (immature phenolic compounds, herbaceous aromas, etc.) prevail in the wine even when a small proportion of grapes with those characteristics are present [53]. All these considerations are key in the quality of wine, allowing a better use of resources from an environmental and food safety point of view [54].
Our results show that during the 2017 season, the R 2 values of the studied quality parameters are better starting in April. The results obtained in both seasons do not allow to establish a common pattern of behavior in both seasons. This difference is appreciated in the comparison of the values obtained in the microvinification; in the case of anthocyanins, in some lots the content practically doubles, not a very strange condition if we consider that the typical concentrations of free anthocyanins in full-bodied young red wines are around 500 mg/L, but in some cases can be higher than 2000 mg/L [55][56][57][58].
Several studies have established that anthocyanins begin to be produced more intensely from veraison onwards [58], but the process is complex as it is affected by multiple factors that could affect berry composition and anthocyanin accumulation, such as environmental conditions [59] and agricultural management [60]. These factors can significantly influence the anthocyanin content of grapes [61]. For example, water deficit has been shown to have a positive effect on anthocyanin accumulation during ripening [62,63] and others have shown that water limitation has negative effects [59,64].
The fact that the highest modelling potential of phenolic compounds occur later in 2017 may be justified by the study of Mori et al. (2007) who showed that high tempera-tures reduced the total anthocyanin content to less than half of what was expected [64]. A number of other studies agree that temperature is a critical variable in the final amount of anthocyanins present in wine [6,65,66]. The comparison of the analysis of the microvinifications for both seasons shows higher absorbance values for both color intensity and phenolic compounds for the 2018 season, in agreement with other studies that obtained better results with cooler temperatures [66]. This seems to be one of the key factors to consider [67], since vegetative growth and phenolic compound synthesis are affected by temperature [68]. Global warming produces an increase in the gap between technological maturity, corresponding to the stage of higher sugar accumulation and low acidity (which comes earlier) and phenolic and aromatic maturity (which comes later) [69,70] which is ratified by the results obtained in this study, where the best correlations occur in April the first season (warmer). These works suggest the convenience to include temperature variables in the models, which could eventually improve the results obtained.
The study of the variables directly involved in wine quality, such as phenolic compounds and color, contributes to the interest of researchers to delve into studies of nondestructive and affordable techniques based on the use of satellite images that allow understanding the factors that determine the concentration of these quality variables. These models might be useful for the differentiated management of harvests, providing spatial information for selective harvesting, which could be used to categorize different types of wine, and might be complemented with the very extensive information that is progressively available in a vineyard, related to agronomic management, such as nutrition, thinning, leaf removal and grafting, among others.

Conclusions
The purpose of this research was to determine models from Sentinel-2 image series that allow estimating, prior to harvest, wine quality variables such as phenols (TP), anthocyanins (TA) and color intensity (CI). These variables provide sensory characteristics that are evaluated in the wine. Two harvest seasons were studied and models with R 2 coefficients of determination above 60% were obtained for CI in 2017 (75%) and TA in 2017, at 67%. These results transform the estimated TP, TA and CI maps into a useful tool, with the potential to define selective harvests.
Sentinel 2 images play a key role, as their high periodicity allows modeling of preharvest vine conditions. The final stage of ripening and the harvest time is a very short but strategic stretch that would allow identifying homogeneous sectors for harvesting and thus scheduling the winery space.
The results also show the convenience to deepen this study on two specific aspects in the future: (i) to increase the number of samples to improve the dispersion and variety of the data; and (ii) to include a climatic variable, such as cumulative temperature, which could help in understanding the behavior of phenolic compounds and the changes in color in order to obtain more robust models.
The optimization of all processes in the wine production chain aims to raise the quality level of winemaking. Therefore, the production of maps that allows spatially differentiating the variability and quality potential of wine produced in terms of color intensity and phenolic compounds accumulation presents a practical and attractive potential to identify aging capacity and the style of wine expected from different plots as well as their productive aptitude. This work suggests that remote sensing and selective harvesting can be used to better manage the style and quality of wine produced.