---
title: "Guía metodológica para el uso del paquete R devRate"
author: "Francois Rebaudo y Camila Benavides "
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{main_spanish}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

<style>
body {text-align: justify}
</style>

## Recordatorios sobre biología de insectos

Antes de comenzar, es necesario hacer algunos recordatorios sobre la biología 
de los insectos. Los insectos son organismos __ectotérmos__, lo que significa 
que la temperatura de su cuerpo es poco o no está controlada por un proceso 
biológico activo: en otras palabras, producen poco o ningún calor y son 
totalmente dependientes de la temperatura de su entorno. Se oponen a 
los organismos __endotérmos__, como los mamíferos, que producen 
en general su propia fuente de calor para regular su temperatura.

En la mayoría de los ambientes, la temperatura varía en el tiempo -por 
ejemplo, es más fría por la noche que durante el día- y en el espacio, por 
ejemplo, es más fría a la sombra que al sol. Los organismos ectotermos 
evolucionan en ambientes donde los cambios de temperatura varían, a veces 
considerablemente. Hablamos entonces de organismos __poikilotérmos__, entre 
los que encontramos la mayoría de los insectos, que contrastan con los 
organismos __homotérmos__ que mantienen una temperatura relativamente 
constante independientemente de las condiciones externas.

La importancia de la temperatura ambiente en los procesos biológicos de los 
insectos es clave, aunque no es el único factor involucrado. 
Por ejemplo, otros factores pueden ser la humedad relativa, o el fotoperiodo (la 
duración del día y la noche). Aqui nos vamos a enfocar unicamente en la 
temperatura.

Los insectos tienen una preferencia térmica, es decir un valor o conjunto de 
valores de temperatura, dentro de los cuales pueden alcanzar un 
desarrollo óptimo. Esta preferencia es variable dependiendo de la especie, 
con especies que se desarrollarán en un amplio rango de temperaturas. A estas se 
denomina especies generalistas, a diferencia de las especies que pueden 
desarrollarse solo en un rango muy estrecho de temperaturas, que son las 
llamadas especialistas (los insectos se encuentran a lo largo de 
este gradiente entre especialistas y generalistas). Las variaciones de preferencia 
también se pueden observar dentro de una especie, especialmente entre las diferentes 
etapas de desarrollo. Para entender estas diferentes estrategias y 
caracterizarlas, es necesario estudiar la relación entre la temperatura y 
el desarrollo de los insectos. 

Nos interesamos en esta occasión en los insectos poikilotérmos, que, recordamos, 
corresponden a los insectos cuya temperatura cambiará como consecuencia de los 
cambios de la temperatura exterior. Estudiaremos particularmente la tasa 
de desarrollo de insectos, generalmente expresada como la inversa de la 
cantidad de días que pasan de una etapa de desarrollo a otra. Para este 
proposito, usaremos datos disponibles en la literatura.

En la mayoría de los ecosistemas, las condiciones ambientales son variables. 
Este es el caso de la temperatura. Para los insectos poikilotérmos, esto 
dará como resultado un cambio en la temperatura de su cuerpo. Mas alla del 
rango en el cual pueden desarrollarse, sufrirán un estrés que puede llevarles 
a la muerte. Por lo 
tanto, existen umbrales críticos, un umbral mínimo y un umbral máximo, más 
allá de los cuales los insectos no pueden continuar su desarrollo. Esto límites 
se encuentran en la literatura bajo los acrónimos __CTmin__ y __CTmax__ para 
_Critical Thermal minima_ _Critical Thermal maxima_, respectivamente. Los 
mismos varían ampliamente de una especie a otra, con algunas especies como el 
escarabajo ártico capaz de soportar hasta -70 grados, mientras que otras 
admiten hasta mas de 45 grados. Como se mencionó anteriormente, la 
sensibilidad a la temperatura, que clasifica los insectos de acuerdo con un 
gradiente de generalistas a especialistas varía mucho de una 
especie a otra e incluso dentro de una especie entre sus etapas fenológicas.
Más allá de estos dos valores umbral, la respuesta de los insectos variará 
según el tiempo frecuencia y tasa de exposición, la velocidad a la 
que la temperatura se lleva más allá de los valores umbral, 
la historia de vida de los insectos y, por supuesto, las especies 
consideradas. En esta guía nos enfocaremos en lo que sucede entre estos 
umbrales.

Para estudiar qué ocurre entre estos dos valores umbral, mediremos el tiempo de
desarrollo del insecto, expresado como la tasa de desarrollo. La tasa de 
desarrollo de los insectos comienza en el primer valor umbral, luego aumenta 
progresivamente hasta alcanzar un valor óptimo, luego desciende rápidamente al 
segundo valor umbral. La respuesta de los insectos a la temperatura es, por lo 
tanto, no lineal.

La temperatura a la cual el desarrollo es máximo se llama temperatura óptima 
de desarrollo u óptimo térmico (__Topt__). Es importante enfatizar que esta es 
la temperatura a la cual el desarrollo es óptimo, no la temperatura óptima para 
el insecto. De hecho, la temperatura óptima depende de otros criterios 
distintos del desarrollo. Un buen ejemplo es la supervivencia: a la 
temperatura óptima para el desarrollo, la supervivencia observada es diferente 
de la temperatura óptima para la supervivencia.

Los estudios sobre la respuesta de los insectos a la temperatura no son 
nuevos: el primer científico que estudio y publico sobre el efecto de 
la temperatura en los insectos es René-Antoine Ferchault de 
Réaumur^[https://en.wikipedia.org/wiki/Ren%C3%A9_Antoine_Ferchault_de_R%C3%A9aumur]: 
sus observaciones sobre las mariposas estan registradas en sus "Mémoires pour 
servir à l'histoire des insectes" entre 1734 y 1742. Despues de Réaumur, muchos 
científicos buscaron y siguen buscando cuantificar la relación entre la 
temperatura y la tasa de desarrollo.

Para hacer esto, el método más común consiste en un primer paso en colocar 
insectos en jaulas a diferentes temperaturas constantes, para medir el tiempo 
de desarrollo de cada etapa fenológica. Mientras haciendo este experimento, 
podemos observar algunas características de la respuesta de los insectos a la 
temperatura:

1. La primera que resulta muy difícil, si no imposible, de obtener datos de 
desarrollo cercanos a los mínimos y máximos, ya que a estas temperaturas la 
gran mayoría de los insectos de la misma especie y población no sobreviven
2. La segunda característica es que se observa variabilidad dentro de la misma 
población de insectos en cuanto a su respuesta a la temperatura
3) La tercera característica es que para el rango de temperaturas que se 
encuentra en el hábitat del insecto, la relación entre la temperatura y la 
tasa de desarrollo es casi lineal, pero fuera de este rango la relación es 
no lineal

El segundo paso es el ajuste de un modelo 
matemático a las mediciones realizadas en el laboratorio. El objetivo del 
ajuste es encontrar estimadores de parámetros de manera que el modelo esté lo 
más cerca posible de los puntos experimentales y al mismo tiempo retenga un 
significado biológico.

El objetivo de este documento es de proveer una guía para facilitar y analizar 
el ajuste de un modelo matemático a las mediciones realizadas en el laboratorio.

## Datos del laboratorio

Los datos de laboratorio que vamos a usar son los de _Bactrocera dorsalis_
^[https://en.wikipedia.org/wiki/Bactrocera_dorsalis], un insecto del Orden de 
los Dipteros, que hemos obtenidos a partir del articulo de Messenger and 
Flitters (1958)^[Messenger, P.S., and Flitters, N.E. (1958). Effect of constant 
temperature environments on the egg stage of three species of Hawaiian fruit 
flies. Annals of the Entomological Society of America, 51(2), 109-119.]. Los 
datos coresponden al tiempo de desarrollo de los huevos expresado en horas a 
distintas temperaturas expresadas en Fahrenheit.

Para seguir la guía, es nesesario contar con el programa R instalado en su 
computadora. Para usar el paquete devRate, hay que instalarlo mediante:

```
install.packages("devRate")
```

Una vez el paquete instalado, no es necesario repetir a futuro la linea de 
instalacion. Sin embargo instalar el paquete no es suficiente para poder 
usarlo dentro de R. Hay que especificar que vamos a usarlo cada vez que usamos 
R. Para esto se puede usar: 

```{r}
library("devRate")
```

A partir de ahora vamos a trabajar desde un archivo que podriamos llamar 
`B_dorsalis_devRate.R`. Este archivo empeza con un comentario indicando de que 
se trata para acordarse de su contenido. Todo lo que viene despues del 
simbolo `#` son comentarios.

```{r}
### script para analizar los datos de tasa de desarrollo de B. dorsalis en 
### función de la temperatura.
require("devRate") # para cargar el paquete devRate
```


Para que los datos sean reconocidos en R, tenemos que usar un formato especial. 
Los datos de temperatura son almacenados en un formato vector que se escribe 
con `c()` en R y que corresponde a una colección ordenada de elementos del 
mismo tipo:

```{r}
c(55.0, 56.0, 57.0, 58.0, 60.0, 62.5, 65.0, 67.5, 70.0, 75.0, 80.0, 85.0, 87.5, 
  90.0, 92.5, 95.0, 96.0, 97.0, 97.5)
```

Los numeros usan el punto como separator de decimales y la coma como separator 
de valores. Aqui las temperaturas estan en Fahrenheit, tenemos que convertirlas 
en Celsius con la formula `T(°C) = (T(°F) - 32) / 1.8`.

```{r}
(c(55.0, 56.0, 57.0, 58.0, 60.0, 62.5, 65.0, 67.5, 70.0, 75.0, 80.0, 85.0, 87.5, 
  90.0, 92.5, 95.0, 96.0, 97.0, 97.5) - 32) / 1.8
```

Para guardar las temperaturas, podemos almacenarlas en un objeto. Para almacenar 
los valores de temperaturas en un objeto, se usa los simbolos `<-`. Si queremos 
guardar los valores de temperatura en un objeto llamado `temp`:

```{r}
temp <- (c(55.0, 56.0, 57.0, 58.0, 60.0, 62.5, 65.0, 67.5, 70.0, 75.0, 80.0, 
           85.0, 87.5, 90.0, 92.5, 95.0, 96.0, 97.0, 97.5) - 32) / 1.8
```

De la misma manera, vamos a guardar los valores de tasa de desarrollo en un 
objeto llamado `devRate`. La tasa de desarrollo corresponde al opuesto del 
numero de días. Aquí los valores estan expresados en horas, entonces dividimos 
el numero de horas por 24 y despues calculamos el opuesto.

```{r}
devRate <- 1/(c(263.0, 232.0, 170.5, 148.0, 121.3, 95.5, 74.0, 62.5, 51.5, 
                38.0, 30.5, 27.0, 25.0, 24.0, 23.5, 25.0, 26.5, 29.3, 34.3)/24)
```

Ahora podemos almacenar los datos de temperatura y de tasa de desarrollo en una 
tabla `datosLab` donde temperatura y tasa de desarrollo son las columnas. La 
función `data.frame` permite hacer la tabla. 

```{r}
datosLab <- data.frame(temp, devRate)
```

Si queremos hacer la tabla y llenarla, hay que especificar el nombre que 
queremos dar a las columnas con el simbolo `=` :

```{r}
datosLab <- data.frame(
  temp = (c(55.0, 56.0, 57.0, 58.0, 60.0, 62.5, 65.0, 67.5, 70.0, 75.0, 80.0, 
            85.0, 87.5, 90.0, 92.5, 95.0, 96.0, 97.0, 97.5) - 32)/1.8, 
  devRate = 1/(c(263.0, 232.0, 170.5, 148.0, 121.3, 95.5, 74.0, 62.5, 51.5, 
            38.0, 30.5, 27.0, 25.0, 24.0, 23.5, 25.0, 26.5, 29.3, 34.3)/24))
```

Podemos verificar el contenido del objeto `datosLab` con la función `print` :

```{r}
print(datosLab)
```

Nuestros datos estan listos para empezar el análisis con el paquete devRate.

## Seleccionar un modelo y ajustarlo

### El modelo lineal

Aunque tenga limitaciones, el método más común y simple es usar los 
resultados obtenidos en la parte casi lineal de la relación entre la 
temperatura y la tasa de desarrollo. El modelo matemático utilizado es 
un modelo lineal, la línea de ecuación `y = bb * x + aa` que todos conocen de 
sus clases de estadistica. El uso del modelo lineal fue popularizado por el 
articulo cientifico de Campbell et al. (1974)^[Campbell, A., Frazer, B.D., 
Gilbert, N.G.A.P., Gutierrez, A.P., Mackauer, M. (1974). Temperature 
requirements of some aphids and their parasites. Journal of applied ecology, 
431-438.]. 
Este modelo se basa en la suposición de que la 
relación entre temperatura y desarrollo en insectos es lineal en un cierto 
rango de temperaturas y sirve de base para el concepto de "grado-dias".
La extensión de la línea de regresión lineal fuera de la zona de linealidad 
intersecta el eje x. El valor de la temperatura en esta intersección se 
denomina temperatura base, que suele ser diferente del CTmin. La temperatura 
base es también conocido como el "umbral de desarrollo" y puede calcularse 
simplemente dividiendo el opuesto de la intersección en y por el coeficiente 
"a". Pero primero hay que ajustar el modelo a los datos de laboratorio.

Para seleccionar los datos que vamos a usar, puede ser util 
visualizar los datos de laboratorio en un grafico. Para esto vamos a usar la 
función `plot`.

```{r}
plot(x = datosLab$temp, y = datosLab$devRate, 
     xlab = "Temperatura", ylab = "Tasa de desarrollo", 
     xlim = c(0, 40), ylim = c(0, 1.2))
```

Podemos observar que los ultimos cinco puntos salen de la zona de linearidad.
Para no usar estos ultimos puntos, vamos a seleccinar una parte de los 
datos de laboratorio. Quitar los ultimos cinco elementos es igual a seleccionar 
unicamente las 14 primeras lineas de nuestra tabla. Se usan corchetes para 
seleccionar una parte de los datos, y una coma para especificar que se trata 
de lineas.

```{r}
datosLab14 <- datosLab[1:14,]
plot(x = datosLab14$temp, y = datosLab14$devRate, 
     xlab = "Temperatura", ylab = "Tasa de desarrollo", 
     xlim = c(0, 40), ylim = c(0, 1.2))
```

Aun que se podria ajustar directamente el modelo lineal con la función `lm`, 
vamos a usar la función generica del paquete devRate.

```{r}
modLin <- devRateModel(eq = campbell_74, dfData = datosLab14)
plot(x = datosLab$temp, y = datosLab$devRate, 
     xlab = "Temperatura", ylab = "Tasa de desarrollo", 
     xlim = c(0, 40), ylim = c(0, 1.2)) # datos del laboratorio
abline(modLin) # añadir modelo lineal en el grafico
print(modLin) # imprimir resultados del ajuste del modelo
```

El concepto de grados-días se basa en la hipótesis de que los 
insectos necesitan un mínimo de acumulación de temperatura 
para completar su desarrollo. El ciclo de vida de los insectos esta divido en 
fases. Cada fase del desarrrollo va a necesitar su mínimo de acumulación de 
temperatura para que el inscto pueda llegar a la fase siguiente y completar su 
desarrollo. La acumulación de temperatura no empeza desde cero grados, sino de 
una temperatura base que depende del insecto y de la fase de desarrollo.

Esta guía no discutirá los numerosos métodos 
disponibles para calcular los grados-días. Sin embargo, el método más simple es 
agregar el valor de temperatura mínima a la temperatura máxima de un día y 
luego dividir el resultado de la suma entre dos para luego restar a este 
resultado el valor de la temperatura base. Es decir, para un día determinado 
cuando la temperatura máxima es de 25 grados, y la temperatura mínima de 15 
grados, y con una temperatura base de una especie de insecto dada de 3 grados, 
el cálculo de los grados los días acumulados en este día son 25 más 15 es 
igual a 40; dividido entre 2 es igual a 20; menos 3 es igual a 17 grados-días.
El individuo completa su desarrollo cuando la suma 
de los grados-día alcanza el valor de su constante térmica, determinada a 
partir de la línea de regresión. La constante térmica denotada por "K" es igual 
a la inversa del coeficiente director de la línea de regresión, es decir, K 
es igual a 1 por bb con `y = bbx + aa`, donde `y` representa el desarrollo y 
`x` la temperatura.

Aqui `bb = 0.04881` y `aa = -0.55176`, entonces: 

* la temperatura base es `Tbase = -aa/bb = 11.30424`
* la constante térmica es `K = 1/bb = 20.4876`.

Podemos buscar en la base de datos del paquete devRate si hay otros modelos 
que fueron usados para la especie _Bactrocera dorsalis_.

```{r}
devRateFind(species = "Bactrocera dorsalis")
campbell_74$startVal[campbell_74$startVal$genSp == "Bactrocera dorsalis",]
```

Al momento de escribir esta guía, en la base de datos existen dos estudios que 
han trabajado con _Bactrocera dorsalis_ usando el modelo lineal, pero ninguno 
ha trabajado unicamente con los huevos, así que no podemos comparar los 
resultados obtenidos.

Los modelos de grados-días tienen una aplicación directa para predecir, a 
veces con precisión, la aparición de un insecto basándose en las proyecciones 
de temperatura futuras, y así servir como un sistema experto para muchas 
aplicaciones agrícolas o epidemiológicas.
El concepto de grados-días tiene la ventaja de ser simple y, por lo tanto, 
fácil de implementar a partir de un conjunto de datos restringidos, pero 
tiene sus límites.
En primer lugar, tiende a subestimar el desarrollo de insectos a  
temperaturas bajas, en el área entre la temperatura base y la temperatura 
umbral de desarrollo CTmin. Además, no permite determinar la temperatura 
umbral de desarrollo más baja (CTmin), ni superior (CTmax). Luego, sobreestima el 
desarrollo de los insectos tan pronto como uno se acerca a la temperatura 
óptima de desarrollo. Por lo tanto, solo puede usarse en un rango limitado 
de temperaturas y brinda poca información sobre la relación entre la 
temperatura y la tasa de desarrollo.
La alternativa al concepto de grados-días es el uso de modelos no lineales.

### Los modelos no lineales

Una de las principales motivaciones para el desarrollo de ecuaciones no 
lineales es el cálculo de CTmin, CTmax, Topt. Existen varios modelos no 
lineales. La lista de los que estan en el paquete devRate se puede encontrar 
en la documentacion mediante `?devRate` o a traves del objeto `devRateEqList`.

```{r}
names(devRateEqList)
```

Dentro de todos estos modelos, la seleccion de uno en particular no tiene consenso en 
la comunidad cientifica. La mayoría de los modelos no lineales se basan en la 
descripción de la curva de respuesta a la temperatura de los insectos, con un 
enfoque matemático sin base biológica real, por lo que los diferentes modelos 
no pueden separarse en sus supuestos subyacentes, y los criterios estadísticos 
para decidir entre ellos no permiten que un modelo se afirme a sí mismo como 
mejor que los demás.
Otros modelos no lineales, llamados modelos biofísicos, se basan en la 
activación enzimática, como el modelo Sharpe y DeMichèle de 1977, modificado 
en 1981 por Schoolfield et al. Lamentablemente, estos modelos son más 
complejos y requieren conjuntos de datos vastos para ser utilizados. 
Además, algunos estudios cuestionan su base teórica. Al final, se usan poco.

Para cuantificar la respuesta de los insectos a la temperatura, se trata de 
ajustar varios modelos al conjunto de datos experimentales, y seleccionar 
el que mejor responda a la problemática del estudio, y evaluarlo con los 
criterios de comparación estadísticos.
Una posibilidad es de buscar en la base de datos del 
paquete cuales son los mas usados con la función `devRateFind` pero como hemos 
visto anteriormente no existen modelos no lineales que hayan sido usados al momento 
de escribir esta guía. Una alternativa es usar los modelos mas comunes como 
los de Briere, Logan, y Lactin. Para esto, hay que definir valores 
inciales para los parametros de 
los modelos para que el algoritmo que va a hacer el ajuste encuentre una 
solución. En caso de que no se conocen valores pertinentes para iniciar el 
algoritmo, se puede usar los valores promedios que se encuentran en la tabla 
`startVal` de cada modelo (por ejemplo `briere1_99$startVal` para el modelo 
de Briere-1). 

```{r}
devRateFind(species = "Bactrocera dorsalis")
modNoLin_01 <- devRateModel(
  eq = briere1_99, # nombre del modelo
  dfData = datosLab, # nombre de los datos de laboratorio
  startValues = list(aa = 0.01, Tmin = 10, Tmax = 40)) # valores iniciales
print(modNoLin_01)
modNoLin_02 <- devRateModel(
  eq = briere2_99, # nombre del modelo
  dfData = datosLab, # nombre de los datos de laboratorio
  startValues = list(
    aa = 0.01, Tmin = 10, Tmax = 40, bb = 2)) # valores iniciales
print(modNoLin_02)
modNoLin_03 <- devRateModel(
  eq = lactin2_95, # nombre del modelo
  dfData = datosLab, # nombre de los datos de laboratorio
  startValues = list(
    aa = 0.03, Tmax = 30, deltaT = 5.0, bb = -1.5)) # valores iniciales
print(modNoLin_03)
```

Los resultados se pueden visualizar graficamente con la función `devRatePlot`.

```{r, fig.width = 7, fig.height = 3}
par(mfrow = c(1, 3)) # para hacer tres graficos en la misma pagina
devRatePlot(eq = briere1_99, 
  nlsDR = modNoLin_01, 
  xlim = c(10, 40), ylim = c(0, 1.2))
devRatePlot(eq = briere2_99, 
  nlsDR = modNoLin_02, 
  xlim = c(10, 40), ylim = c(0, 1.2))
devRatePlot(eq = lactin2_95, 
  nlsDR = modNoLin_03, 
  xlim = c(10, 40), ylim = c(0, 1.2))
```

Existen varios indices estadisticos para comparar los modelos, uno de esos es 
AIC. El modelo con menor AIC es el que se ajusta mejor a los datos. No 
significa que el modelo es valido, sino que tiene el mejor ajuste dentro de 
los modelos seleccionados.

```{r}
c(AIC(modNoLin_01), AIC(modNoLin_02), AIC(modNoLin_03))
```

En este caso, el modelo de Briere-2 se ajusta mejor a los datos que los otros 
modelos. Ahora podemos calcular los valores de CTmin, CTmax, y Topt. Primero 
vamos a simular datos de temperatura desde 0 hasta 45 grados. Despues vamos a 
hacer una predicción de la tasa de desarrollo para esas temperaturas en base al 
modelo de Briere-2. Despues especificamos que los valores de tasa de 
desarrollo negativos son de cero (porque no pueden ser negativos), y que los 
valores que faltan tambien son de cero (en caso de que haya valores ausentes). 
El Topt es el valor de temperatura a la cual la tasa de desarrollo esta máxima, 
el CTmin es el valor de temperatura mínima a la cual observamos desarrollo y el 
CTmax el valor de temperatura máxima a la cual observamos desarrollo.

```{r}
tempS <- seq(from = 0, to = 45, by = 0.1) # temperaturas simuladas
devRateS <- predict(modNoLin_02, newdata = list(T = tempS)) # predicciones
devRateS[devRateS < 0] <- 0
devRateS[is.na(devRateS)] <- 0
c(AIC(modNoLin_01), AIC(modNoLin_02), AIC(modNoLin_03))
Topt <- tempS[devRateS == max(devRateS)]
CTmin <- tempS[devRateS == min(devRateS[devRateS > 0 & 
  tempS < Topt])]
CTmax <- tempS[devRateS == min(devRateS[devRateS > 0 & 
  tempS > Topt])]
cat(paste0("Topt: ", Topt, "\nCTmin: ", CTmin, "\nTmax: ", CTmax))
```

## Construir modelos fenológicos

Mediante el uso de la tasa de desarrollo, es posible modelar el área de 
distribución teórica de una especie, es decir, el área en la que el desarrollo 
sería posible considerando la temperatura, así como las 
diferentes etapas y el número de generaciones potenciales por año de una 
especie; el voltinismo.
El primer paso es tener datos de temperatura. El segundo paso 
consiste en utilizar la relación entre la temperatura y la tasa de desarrollo 
cuantificada previamente.

De hecho, para cada dato de temperatura, se puede calcular una tasa de 
desarrollo que, al acumularse, conducirá al ciclo completo de desarrollo de una 
especie determinada.
Así obtenemos una distribución potencial basada en un modelo de fenológia, 
que puede modificarse mediante la introducción de escenarios de cambio 
climático como un aumento general de las temperaturas o la introducción de 
eventos extremos más frecuentes.

Para esto vamos a imaginar un entorno teórico. Este entorno tendra 
una temperatura sacada de una ley Normal de parámetro  mu igual a 15 grados 
para el promedio, de parametro sigma igual a 1 para la desviación estandar, 
durante un período de 100 días. 

Nuestro modelo de _Bactrocera dorsalis_ solo es valido para los huevos, mientras que 
necesitamos un modelo para todos los estadios fenologicos. Vamos a usar entonces otros 
datos sacados de la literatura. Los datos seran de _T. solanivora_, 
una polilla de la papa (Lepidoptera:Gelechiidae) del el articulo cientifico 
de Crespo-Perez et al. 2011^[Crespo-Pérez, V., Rebaudo, F., Silvain, J.-F. 
Dangles, O. (2011) Modeling invasive species spread in complex landscapes: the 
case of potato moth in Ecuador. Landscape Ecology, 26, 1447–1461.], usando Web
Plot Digitizer^[Rohatgi, A. (2015) WebPlotDigitalizer: HTML5 Based Online Tool 
to Extract Numerical Data from Plot Images.]. 

```{r}
## cargar datos del laboratorio
datosLabTS_egg <- data.frame(
  temp = c(10.0, 10.0, 13.0, 15.0, 15.0, 15.5, 16.0, 16.0, 17.0, 20.0, 20.0, 
    25.0, 25.0, 30.0, 30.0, 35.0), 
  devRate = c(0.031, 0.039, 0.072, 0.047, 0.059, 0.066, 0.083, 0.1, 0.1, 0.1, 0.143, 
    0.171, 0.2, 0.2, 0.18, 0.001))
datosLabTS_larva <- data.frame(
  temp = c(10.0, 10.0, 10.0, 13.0, 15.0, 15.5, 15.5, 15.5, 17.0, 20.0, 25.0, 
    25.0, 30.0, 35.0), 
  devRate = c(0.01, 0.014, 0.019, 0.034, 0.024, 0.029, 0.034, 0.039, 0.067, 0.05, 
    0.076, 0.056, 0.0003, 0.0002))
datosLabTS_pupa <- data.frame(
  temp = c(10.0, 10.0, 10.0, 13.0, 15.0, 15.0, 15.5, 15.5, 16.0, 16.0, 17.0, 
    20.0, 20.0, 25.0, 25.0, 30.0, 35.0), 
  devRate = c(0.001, 0.008, 0.012, 0.044, 0.017, 0.044, 0.039, 0.037, 0.034, 0.051, 
    0.051, 0.08, 0.092, 0.102, 0.073, 0.005, 0.0002))
## ajustar un modelo a los datos (Lactin-1)
## ver vignette quickUserGuide para mayor informacion
modTs01_egg <- devRateModel(
  eq = lactin1_95, 
  dfData = datosLabTS_egg, 
  startValues = list(aa = 0.177, Tmax = 36.586, deltaT = 5.631))
modTs01_larva <- devRateModel(
  eq = lactin1_95, 
  dfData = datosLabTS_larva, 
  startValues = list(aa = 0.177, Tmax = 36.586, deltaT = 5.631))
modTs01_pupa <- devRateModel(
  eq = lactin1_95, 
  dfData = datosLabTS_pupa, 
  startValues = list(aa = 0.177, Tmax = 36.586, deltaT = 5.631))
```

A partir de los modelos ajustados a _T. solanivora_ y una población de 50 
individuos simulados, podemos representar la distribución de las  
generaciones con el tiempo. Para agregar realismo a la simulación, añadimos 
variabilidad en la respuesta de los insectos a la temperatura que seguirá 
una distribución Normal de parametro mu igual al valor de la tasa de 
desarrollo y de parametro sigma = 0.025. 

```{r, fig.width = 7, fig.height = 3}
simul01 <- devRateIBM(
  tempTS = rnorm(n = 100, mean = 15, sd = 1),
  timeStepTS = 1,
  models = list(modTs01_egg, modTs01_larva, modTs01_pupa),
  numInd = 50,
  stocha = 0.025,
  timeLayEggs = 1)
print(simul01)
par(mar = c(4, 4, 0, 0))
devRateIBMPlot(ibm = simul01)
```

En los resultados obtenemos la fenología de cada uno de los 50 individuos, el 
modelo utilizado, y la serie temporal de temperatura. En la visualización 
de los resultados, cada color representa una generación diferente y cada tipo 
de línea una etapa fenológica diferente. 

Dado que estamos interesados en la respuesta de los insectos al cambio climático, 
es posible modificar el escenario inicial, es decir, modificar la serie temporal 
de temperaturas introduciendo un aumento en las temperaturas, por ejemplo pasando 
de un promedio de 15 a 17 grados.

```{r, fig.width = 7, fig.height = 3}
simul02 <- devRateIBM(
  tempTS = rnorm(n = 100, mean = 17, sd = 1),
  timeStepTS = 1,
  models = list(modTs01_egg, modTs01_larva, modTs01_pupa),
  numInd = 50,
  stocha = 0.025,
  timeLayEggs = 1)
par(mar = c(4, 4, 0, 0))
devRateIBMPlot(ibm = simul02)
```

Tambien podemos aumentar la variabilidad de las temperaturas al cambiar la 
desviación estándar de la ley normal que representa las temperaturas, y pasar 
sigma de 1 a 2. 

```{r, fig.width = 7, fig.height = 3}
simul02 <- devRateIBM(
  tempTS = rnorm(n = 100, mean = 17, sd = 2),
  timeStepTS = 1,
  models = list(modTs01_egg, modTs01_larva, modTs01_pupa),
  numInd = 50,
  stocha = 0.025,
  timeLayEggs = 1)
par(mar = c(4, 4, 0, 0))
devRateIBMPlot(ibm = simul02)
```

Variando los parámetros incluso ligeramente, podemos ver que las consecuencias 
pueden ser importantes para la fenología de los insectos, a pesar de la 
simplicidad del modelo aquí representado. Podemos modificar el modelo para 
cambiar el tamaño de las poblaciones, aumentar el número de días de la simulación 
o modificar la variabilidad intrapoblacional en la respuesta a la temperatura para 
profundizar estos elementos y observar la complejidad de la respuesta de los 
insectos a los cambios globales.

También podemos usar un modelo diferente al utilisado y resaltar la 
importancia de la elección del mismo, y concluir que las incertidumbres son 
múltiples:

* Incertidumbres relacionadas al cambio climático: ¿cuáles serán 
las temperaturas que vendrán y cuánto variarán estas temperaturas?
* Incertidumbres sobre los insectos: ¿cuál es la variabilidad intrapoblacional 
de la respuesta de los insectos a la temperatura? ¿Cuál es el efecto de las 
temperaturas fluctuantes? ¿Cómo transferir estudios de laboratorio a condiciones 
externas?
* Incertidumbres relacionadas con la caracterización de la respuesta 
de los insectos a la temperatura y en particular las incertidumbres relacionadas 
con la elección del modelo matemático.

Aquí el modelo presentado es bastante simple, podría enriquecerse teniendo en 
cuenta la variabilidad espacial en las temperaturas. En el ejemplo que sigue 
vamos a crear un entorno teórico para representar un mapa donde a cada punto 
corresponde un dato de temperatura que vamos a usar para calcular el numero 
de generaciones potenciales simulando 100 días a esta temperatura. 

```{r, fig.width = 5, fig.height = 4}
tempEspacio <- matrix(rnorm(100, mean = 15, sd = 1), ncol = 10) # mapa teorico
myDevRate <- 1/devRateMap(nlsDR = modTs01_egg, tempMap = tempEspacio) +
  1/devRateMap(nlsDR = modTs01_larva, tempMap = tempEspacio) +
  1/devRateMap(nlsDR = modTs01_pupa, tempMap = tempEspacio)
filled.contour(100 / myDevRate, # numero de generaciones en 100 dias
	col = rev(heat.colors(30)), 
	main = "#generaciones", plot.axes = FALSE)
```

Tambien se podria incorporar la respuesta de los insectos a otros factores 
como la humedad relativa, o la disponibilidad y calidad del recurso disponible 
para estos. 

## Límites

Aunque las predicciones de estos modelos son generalmente correctas a escala 
regional o continental, las predicciones a escalas más locales se enfrenten con 
las limitaciones de estos modelos:

1) solo se tiene en cuenta la temperatura, mientras que 
otros factores son importantes para el desarrollo de insectos, como la humedad 
relativa, el fotoperiodo o la fuente y la disponibilidad de comida
2) las curvas de desarrollo dependientes de la temperatura generalmente se 
establecen a temperaturas constantes de laboratorio, pero en condiciones 
externas las temperaturas son fluctuantes, lo que tiene un efecto en 
el desarrollo, incluso cuando las temperaturas no exceden los valores umbral 
de desarrollo. Este punto es aún más importante cuando sabemos que el cambio 
climático no solo influirá en el aumento de la temperatura promedia, sino que 
también exacerbará los extremos, con temperaturas que fluctuarán aún más
3) las tasas de desarrollo calculadas a temperaturas cercanas a Ctmin y 
Ctmax pueden variar significativamente de un modelo a otro, con la 
incertidumbre sobre el efecto de las temperaturas cercanas a estos valores 
umbral
4) dentro de la misma especie y en la misma etapa fenologica, hay una 
variabilidad en la respuesta de los insectos a la temperatura que no es 
tomada en cuenta por los modelos que cuantifican la relación entre la 
temperatura y la tasa de desarrollo. Los modelos mecanísticos de la 
fenología de insectos tenderán a sobreestimar o subestimar el desarrollo 
de una parte de la población.

Estos son los desafíos a los que se enfrentan los estudios de este tipo hoy 
en día y que hay que conocer antes de interpretar los resultados de los modelos 
presentados en esta guía. 

## Conclusión

Los modelos que relacionan temperatura y tasa de desarrollo son la base de 
casi todos los modelos fenologicos predictivos que se pueden encontrar en 
la literatura hoy en día. Se ha escrito la presente guía con el proposito de ayudar a 
analizar sus datos, interpretar los resultados, y construir sus modelos 
fenologicos conociendo las fortalezas y debilidades de este enfoque. 
Cualquier comentario o pregunta son bienvenida(o)s, no dude en ponerse en 
contacto conmigo.
