Introducción

En este documento reportaré el procedimiento seguido para crear un algoritmo que prediga la fecha de la cosecha tomando en cuenta la variedad, carga y fecha de aplicación de cianamida, tomando en cuenta datos de clima histórico.

Descripción del algoritmo objetivo

Para poder cumplir con el objetivo de predecir la fecha de cosecha en base a variedad, carga, clima y fecha de aplicación de cianamida tenemos que hacer lo siguiente.

  1. Calcular la cantidad de horas/calor sobre 16 grados (solo los positivos) para la data histórica. Sumar el total por el día. Promediar las fechas iguales (ignorando el año).

  2. Predecir la cantidad de horas/sol necesarias (del año específico que tenemos la data, no del histórico), la carga y la variedad.

  3. En base a ese modelo, podemos predecir la cantidad de horas que son necesarias para la madurez dependiendo de la carga y variedad

  4. Buscar en la data histórica cuantos días son necesarios desde la aplicación de la cianamida para lograr esa cantidad predicha de horas sol

  5. Sumar esa cantidad de días a la aplicación de la cianamida para tener la fecha estimada de cosecha

Procedimiento

Cargar las bases de datos

Comenzamos por cargar los datos de Clima (Clima) y de maduración de uvas (Uvas). En ambos casos estamos seleccionando únicamente las variables que vamos a emplear. En cuanto a la temperatura, estamos usando la variable Temperatura Exterior como indicador, y dejando de lado la máxima y la mínima que también están disponibles en la base de datos.

Revisión de los datos iniciales

La calidad de un modelo predictivo, dependerá de que los insumos de los que se alimenta sean precisos. Vamos a revisar que no hayan cosas extrañas en los datos de la fuente para asegurarnos de la calidad de los datos.

Estos son los datos que tenemos de temperatura cada intervalos de 15 min. Las temperaturas oscuilan entre los 10 y 35 grados, y es notoria la ausencia de registro durante varios periodos antes del 2017.

En la mayor parte de los datos de clima, están dados por intervalos de 15 minutos. Sin embargo, no hay 94 registros en todos los días. En el gráfico vemos en verde los dias con 96 registros y en rojo los días con más o menos registros. Esto es porque hay días con registros faltantes, o días donde se registra la temperatura cada 5 o 10 minutos. Para uniformizar el registro de clima, promediaremos la temperatura a nivel de hora.

Podemos ver que aún hay 15 días con menos registros. Esto es un problema porque esto puede generar vacíos en el registro del clima promedio, por lo que será necesario imputar estos valores con valores probables.

Después de un proceso de limpieza de datos, ahora tenemos unicamente registros de temperatura promedio a nivel de hora. La falta de estandarización a nivel de intervalo de datos nos impide recuperar la mayor granularidad en los datos. Los registros ausentes fueron completados con el algoritmo mice, dado que no es posible realizar una suma del total del día si el día no cuenta con un registro completo de 24 horas.

Definir el indice de calor más apropiado

El primer paso es definir que indice de calor es el más apropiado. El IC se define como la sumatoria de grados por encima de cierto nivel (usualmente se usa 10C) dentro de cada hora, despues de la aplicación de la cianamida. Esto consistiría (aproximadamente en el area pintada de azul en el siguiente gráfico), con la salvedad de que se calcula desde la aplicación de la cianamida hasta la fecha de cosecha ideal.

El IC más apropiado es aquel que tenga la mayor correlacion con la fecha ideal de cosecha.

Debido a como se relacionan los indices de calor con la fecha de cosecha, usaremos el indice de calor en base a 18 grados.

Calcular el indice de calor en el clima histórico

En toda la data historica calcularemos el indice de calor por día. Para esto: 1. Restamos 18 a la temperatura registrada cada 15 minutos 2. Tomamos solo los valores positivos 3. Los sumamos 4. Los promediamos de acuerdo a la fecha, sin tomar en cuenta el año.

El gráfico muestra el nivel promedio en la data histórica de horas calor acumuladas en cada día del primer mes del año.

Predecir la cantidad de horas calor necesarias para la madurez (año 2019)

Ahora, dependiendo de la variedad y la carga, construiremos un modelo que prediga cuanto IC requiere la planta para alcanzar la madurez.

ModHoras = lm(B18 ~ Variedad + Carga,data = Uvas) 
ModHoras %>% summary()
## 
## Call:
## lm(formula = B18 ~ Variedad + Carga, data = Uvas)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2129.5  -796.5   -88.0   842.8  2733.2 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               11576.027   1610.841   7.186 1.03e-08 ***
## VariedadCandy Snaps       -4513.927   1828.751  -2.468   0.0179 *  
## VariedadCotton Candy      -1653.285   1408.325  -1.174   0.2474    
## VariedadJack's Salute       494.486   1341.733   0.369   0.7144    
## VariedadSweet Celebration   -73.319   1399.333  -0.052   0.9585    
## VariedadSweet Globe        -458.214   1283.715  -0.357   0.7230    
## VariedadSweet Jubile       -603.359   1794.777  -0.336   0.7385    
## VariedadSweet Sapphire      286.842   1406.685   0.204   0.8395    
## Carga                        -1.605      0.360  -4.458 6.54e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1257 on 40 degrees of freedom
## Multiple R-squared:  0.3943, Adjusted R-squared:  0.2732 
## F-statistic: 3.255 on 8 and 40 DF,  p-value: 0.006

Podemos ver que la cantidad de IC necesaria no es dependiente de manera muy cercana con la variedad y la carga (Solo se explica el 39.43% de la variabilidad en el IC). Además, podemos ver que la carga se relaciona negativamente con la cantidad de horas necesarias, Por cada caja de carga disminuye el IC por 1.61. Esto no corresponde con la información agricultural previa. Esto podría ser un indicador de error en la calidad de los datos. La ecuación completa es: \[B18 = 11576.0 - 4513.9 * VariedadCandy Snaps - 1653.3 * VariedadCotton Candy + 494.5 * VariedadJack's Salute - \\ 73.3 * VariedadSweet Celebration - 458.2 * VariedadSweet Globe - 603.4 * VariedadSweet Jubile + 286.8 * VariedadSweet Sapphire\\ - 1.6 * Carga\]

Usar el modelo para predecir la cantidad de horas calor necesarias

Una vez que tenemos el modelo de regresión, podemos usar la formula para predecir cuantas horas necesita una planta dependiendo de su variedad y su carga. En el ejemplo podemos ver que podemos esperar que una planta de Candy Snaps con una carga de 3500 cajas necesitará de 1445 de IC para alcanzar el punto optimo de madurez.

predict(ModHoras,newdata = data.frame(Variedad = "Candy Snaps",Carga = 3500))
##        1 
## 1445.397

Calcular cantidad de días después de la aplicación de cianamida se lograrán las horas calor necesarias

Por ahora, podemos calcular el IC necesario para que la planta madure, pero necesitamos integrar esta información con nuestro registro histórico de clima, de manera que podamos crear una función que reciba como input la fecha de aplicación de la cianamida, calcule cual es el IC que la planta va acumulando, y devuelva la fecha de cosecha estimada.

FechaCosecha("Candy Snaps",3500,"2019-06-23")
## [1] "2019-09-07"

Evaluación de la predicción

Vamos a calcular que fecha de cosecha hubiera predicho el modelo para los datos que tenemos y luego compararlos con la fecha ideal de cosecha reportada en la base de datos. Depués podemos calcular la diferencia en días entre la fecha real y la fecha estimada en base al algoritmo.

En la figura podemos ver que el algoritmo sistemáticamente subestima la fecha de cosecha (La fecha real siempre es mayor a la estimada).
Relación entre fecha real y fecha predicha de cosecha

Relación entre fecha real y fecha predicha de cosecha

Si calculamos el tamaños de la desviación, podemos ver que el algorimo falla en promedio por 35.2449 días. El histograma muestra la distribución de los errores del promedio (A), y luego la desviacion de los errores si sumamos los 35 días a la predicción para al menos centrar los errores en 0 (B). En ambos casos los resultados son desalentadores.

Un posible problema es la diferencia entre las horas calor reales acumuladas en el 2019 y las horas calor promedio/históricas

El gráfico explica el problema. La data usada para el modelo es la data del año 2019, pero la data histórica es significativamente diferente de esta. Si fueran iguales, todos los puntos deberían caer sobre la línea. Sin embargo los datos reales demuestran que la temperatura históricos es significativamente más alta que la del 2019, y el ángulo de la linea, ni siquiera es paralelo a la de los puntos lo que indica que el problema tampoco se puede solucionar restando una constante a los datos históricos.

Vamos a ver si la diferencia se encuentra en que las temperaturas son demasiado diferentes, o es una cuestión del calculo del IC. Curiosamente, el gráfico muestra las desviaciones en temperatura por hora entre la primera fecha de aplicación de cianamida hasta la última fecha de cosecha. Podemos ver que la diferencia de temperaturas entre el histórico y el real no muestra un sesgo claro. Los puntos rojos muestran las horas en las que la temperatura estuvo por encima en la data histórica, mientras que los rojos muestran las horas en las que la temperatura mfue mayor en el año 2019.