ropensci / software-review

rOpenSci Software Peer Review.
286 stars 104 forks source link

agroclimatico: Índices y Estadísticos Climáticos e Hidrológicos #599

Open paocorrales opened 11 months ago

paocorrales commented 11 months ago

Nombre de la Persona Encargada: Paola Corrales Usuario GitHub de la Persona Encargada: !--author1-->@paocorrales<!--end-author1-- Usuario GitHub de las Otras Autoras del Paquete: @eliocamp, @yabellini, @NatiGattinoni Repositorio: https://github.com/AgRoMeteorologiaINTA/agroclimatico Versión Enviada: 1.0.0 Tipo de Envio: Estándar Editora: !--editor-->@Pakillo<!--end-editor-- Revisores: @VeruGHub, @pmnatural

Archivo: TBD Versión Aceptada: TBD Idioma: es

Package: agromet
Title: Índices y Estadísticos Climáticos e Hidrológicos
Version: 1.0.0
Authors@R: 
    c(person(given = "Yanina",
           family = "Bellini Saibene",
           role = c("ctb"),
           email = "yabellini@gmail.com"),
      person(given = "Elio",
           family = "Campitelli",
           role = c("aut"),
           email = "elio.campitelli@cima.fcen.uba.ar",
           comment = c(ORCID = "0000-0002-7742-9230")),
      person(given = "Paola",
           family = "Corrales",
           role = c("cre", "aut"),
           email = "paola.corrales@cima.fcen.uba.ar"),
      person(given = "Natalia",
           family = "Gattinoni",
           role = c("aut"),
           email = "gattinoni.natalia@inta.gob.ar"),           
      person(family = "INTA",
           role = c("cph")),
      person("Ruida", "Zhong", role = "cph")
      )
Description: Conjunto de funciones para calcular índices y estadísticos climáticos 
    hidrológicos a partir de datos tidy. Incluye una función para graficar resultados 
    georeferenciados y e información cartográfica.
License: GPL-3 + file LICENSE
Encoding: UTF-8
Language: es
LazyData: true
Roxygen: list(markdown = TRUE)
RoxygenNote: 7.2.3
URL: https://github.com/AgRoMeteorologiaINTA/agromet
BugReports: https://github.com/AgRoMeteorologiaINTA/agromet/issues
Depends: 
    R (>= 2.10)
Imports: 
    readr,
    data.table,
    scales,
    lmomco,
    sf,
    automap,
    png,
    grid,
    ggnewscale,
    tidyr,
    ggplot2,
    rappdirs,
    magrittr,
    kableExtra,
    Rcpp (>= 0.12.0),
    cli,
    SPEI (>= 1.8.1)
Suggests: 
    testthat,
    dplyr,
    vdiffr,
    lubridate,
    covr,
    knitr,
    rmarkdown
VignetteBuilder: knitr
LinkingTo: Rcpp

Alcance

Existen algunos paquetes cómo frost y meteor que incluyen funciones para utilizar datos meteorológicos y calcular indices que en algunos casos están asociadosa la actividad agropecuaria. Otro paquete a mencionar es agroclim que según la documentación también incluye el cálculo de una serie de índices agrometeoroógicos fue archivado por CRAN en abril de este año pero el código esta disponible en un repositorio público.

agromet se diferencia de estos paquetes principalmente por 2 características. Por un lado su documentación está escrita integramente en español lo que facilita su uso en Sudamérica. Por el otro lado, sus funciones estan pensadas para ser compatibles con el ecosistema de paquetes que utilizan la filosofía "tidy", lo que facilita posterior manipulación, visualización y uso de los resultados obtenidos con agromet.

Revisiones Técnicas

Tilda los siguientes items para confirmar que los has completado:

Este paquete:

Opciones de Publicación

Opciones para MEE - [ ] Este paquete es novedoso y será de interés para la mayoría de lectores de la revista. - [ ] El manuscrito que describe el paquete no tiene más de 3000 palabras y está escrito en Inglés. - [ ] Tienes intenciones de archivar el código del paquete en un repositorio a largo plazo, que cumple los requerimientos de la revista (mira las [Políticas de Publicación de MEE (documento en Inglés)](http://besjournals.onlinelibrary.wiley.com/hub/journal/10.1111/(ISSN)2041-210X/journal-resources/policy-on-publishing-code.html)) - (*Alcance: Considera los [Objetivos y Alcance de MEE (documento en Inglés)](http://besjournals.onlinelibrary.wiley.com/hub/journal/10.1111/(ISSN)2041-210X/aims-and-scope/read-full-aims-and-scope.html) para tu manuscrito. No otorgamos garatías de que tu manuscrito esté en el ámbito de MEE.*) - (*Aunque no es requerido, recomendamos tener un manuscrito completamente preparado y en Inglés, al momento de enviar.*) - (*Por favor, no envíes tu paquete de forma separada a Methods in Ecology and Evolution*)

Código de Conducta

ropensci-review-bot commented 4 months ago

@paocorrales, @eliocamp, @yabellini, @NatiGattinoni: please post your response with @ropensci-review-bot submit response <url to issue comment> if you haven't done so already (this is an automatic reminder).

Here's the author guide for response. https://devguide.ropensci.org/authors-guide.html

Pakillo commented 1 month ago

Estimadas @paocorrales, @eliocamp, @yabellini, @NatiGattinoni:

Espero que vaya todo bien. Me preguntaba cómo iba avanzando la revisión del paquete. Si puedo ayudar en algún aspecto o tienen cualquier duda o comentario no duden en decirlo.

Gracias

paocorrales commented 1 month ago

Hola!

En primer lugar, muchas gracias @pmnatural y @VeruGHub por la revisión y los comentarios, y a @Pakillo por ordenarlos y organizar todo.

Incorporé los comentarios y sugerencias al paquete. El detalle fino de los cambios se puede seguir en los issues del repositorio del paquete (1 issue por función o familia de funciones).

Respondiendo a los comentarios específicos:

dias_promedio

No calcula ningun estadístico meteorológico como dice la descripción si no que extrae dia, mes y dia juliano de un set de datos. Previamente el usuario puede haber extraido los días que representan un evento meteorológico de interés, pero no tiene por qué.

La descripción de la funcionalidad no es clara en la ayuda. Sugiero poner "Calcula el valor promedio anual del primer y ultimo dia para una serie de fechas dada". En los ejemplos habría que cambiar "summarise" por "reframe" como aparece en el Readme. En "Value" se dice que el output tiene "variables extras en el caso de hacer el cálculo para distintos grupos", pero ¿cómo se puede especificar para qué grupos quieres el cálculo?

Erratas en la documentación: "seria de fechas"

dias_promedio(datos$fecha) #Ok
dias_promedio(NH0358$fecha) #Ok

helada <- NH0358 |> 
  filter(t_min <= 0) 

dias_promedio(helada$fecha) #Ok

helada |> summarise(dias_promedio(fecha)) #En un pipe, solo si la funcion se escribe dentro de un summarise o reframe se consigue la funcionalidad descrita

#La función no devuelve un output que encaje bien con la filosofía de tidyverse y es difícil de utilizar en un pipe: 

helada |> dias_promedio(fechas = fecha)
helada |> mutate(diasprom = dias_promedio(fecha))

Mejoré la documentación para explicar mejor el uso de la función junto con los verbos de dplyr. Corregí errores en los ejemplos (ahora usa reframe()). Incorporé un nuevo ejemplo que muestra como usar la función agrupando datos. Para esto el paquete tiene un nuevo set de datos que corresponde a una serie temporal de una estación meteorológica similar al que ya se estaba usando.

completar_serie

No se especifica el formato del vector fecha ni de rango en la documentación.

> #Ejecuto el ejemplo #Ok
ejemplo <- data.frame(fechas = seq(as.Date("1985-01-01"), as.Date("2015-12-01"),
                             by = "1 month"),
                    pp = 1)

set.seed(42)
datos_perdidos <- sample(nrow(ejemplo), nrow(ejemplo)/10)
ejemplo$pp[datos_perdidos] <- NA
datos_incompletos <- na.omit(ejemplo) 

completar_serie(datos_incompletos, fechas, resolucion = "1 mes")

#No funciona, pero sigue las instrucciones de la documentación:
completar_serie(datos = datos, 
                fecha = fecha, 
                resolucion = "día",
                rango = as.Date(c("1951-01-01", "1952-12-31")))

Había un error en el código. Ahora funciona tanto si la resolución es "dia" o "1 dia". Mejoré la documentación agregando los detalles a los argumentos mencionados.

decil y anomalia_porcentual

En la descripción se necesita incluir una definición más técnica de qué valores se devuelven, es decir, que hace la función decil y la anomalia_porcentual y referencias. En estadística, los deciles de una serie de datos son estadísticos descriptivos que representan los valores que limitan la serie de datos en 10 partes. Al aplicar la función sin tener mucho detalle en la descripción yo hubiera esperado obtener eso. Sin embargo, en el ejemplo ejecutado debajo obtengo un vector de igual longitud de la serie de datos. Sugiero mejorar la descripción en este sentido. La anomalía porcentual ¿respecto a qué valor de la serie de referencia se calcula? Debería especificarse. ¿Por qué para calcular la anomalia es necesario especificar que no se consideren los NA, pero para decil no? Parece una cuestión de inconsistencia en cómo están escritas las funciones.

decil(datos$t_max) #Ok

anomalia_porcentual(datos$t_max, na.rm = TRUE) #Ok

Mejore la documentación para explicar mejor la salida de la función. Respecto a la diferencia entre las funciones, efectivamente ambas son internamente distintas y deciles() no requiere el argumento na.rm. Las documentamos en conjunto ya que se suelen usar en conjunto y encontramos muchos ejemplos donde funciones documentadas como "familia" comparten algunos pero no todos los argumentos documentados.

scala_temp/pp_XXX

La documentación está incompleta. No hay ejemplos de uso ni se entiende bien que hacen estas funciones. Tampoco se especifica qué valores retorna.

Agregue ejemplos en la documentación y lo que devuelve. Los ejemplos usan los nuevos datos del paquete y son mas simples.

ith

Funciona correctamente y está bien documentada.

NH0358 <- NH0358 %>% #Ok
  mutate(t_media = (t_max + t_min)/2) %>%
  mutate(ith = ith(t_media, hr))

hist(NH0358$ith)

Para esta y otras funciones revisé el output de los ejemplos y cuando eran muy largos ahora solo muestra las primeras 10 lineas.

metadatos_nh

Para facilitar la comprensión a los usuarios convendría especificar que son "estaciones NH". La descripción es incompleta en el sentido de que sólo se especifica una parte de la información que devuelve la función.

Erratas en la documentación: "provincia o provincia", "especificar" debería ser "específicas".


metadatos_nh(codigo = c("0001", "0011")) #Ok
metadatos_nh(provincia ="Chaco") #Ok
metadatos_nh(organismo = "SMN") #Ok
metadatos_nh(lat = c(-30, -20), lon = c(-65, -55)) #Ok

Arreglé los errores de tipeo y extendí la documentación para incluir información sobre metadatos.

kable_inta

Sustituir "Setear" por "Formatear" o similar. Funciona correctamente.


kable_inta(datos[1:10,1:5]) #Ok
kable_inta(metadatos_nh(provincia ="Chaco")) #Ok

Documetación arreglada.

leer_surfer

En el ejemplo se sobreescribe un objeto con el nombre de otra función implementada en el paquete (scala_temp/pp_XXX). Debería estar explicado en algún sitio cual es el link entre ambas.

Erratas en la documentación: "laeer".


escala <- system.file("extdata", "escala_pp_mensual.lvl", package = "agroclimatico")

escala_prueba <- leer_surfer(escala) #Ok
scales::show_col(escala_prueba$colores)

Documentación arreglada y mejores ejemplos (asociadas a las funciones scale_*())

mapa_XXX

Añadiría en la descripción los posibles valores para el argumento "provincias".


ggplot() +
  geom_sf(data = mapa_argentina()) #Ok

ggplot() +
  geom_sf(data = mapa_provincias()) #Ok

ggplot() +
  geom_sf(data = mapa_provincias(provincias = "Chaco", departamentos =  TRUE)) #Ok

ggplot() +
  geom_sf(data = mapa_provincias(provincias = "Rio Negro", departamentos = FALSE)) #Ok

ggplot() +
  geom_sf(data = mapa_argentina_limitrofes()) #Ok

ggplot() +
  geom_sf(data = mapa_departamentos(provincias = "Chaco")) #Ok

Gracias por la sugerencia, ahora incluye la lista de provincias.

mapear

Erratas en la documentación: "veriable", "a el", "masyores", "1500m".

El argumento "breaks" no está bien definido. Habría que especificar que se refiere a que hay que introducir valores (¿uno o varios?) en el rango del argumento "valor".

No se explican las funciones coord_argentina ni theme_inta_mapa en ningún sitio. ¿Se utilizan sólo internamente? En ese caso no deberían estar exportadas.


prueba <- metadatos_nh()

mapear(prueba$altura, #Ok
       prueba$lon, 
       prueba$lat, 
       breaks = c(0, 100, 400, 1000, 1500), 
       cordillera = 1400,
       variable = "mm")

Arreglé los errores de tipeo e incluí las aclaraciones pedidas. Además ahora la función mapear() acepta el data.frame de datos como primer argumento. Las funciones secundarias coord_argentina y theme_inta_mapa ahora tienen documentación.

olas y umbrales

No se explica en la documentación como hay que especificar los umbrales. Es poco habitual que se pueda poner cualquier palabra como nombre de un argumento y una mala praxis porque se puede sobreescribir algo importante en la función. Sería más adecuado que hubiera un argumento para los umbrales y que éstos se pudieran introducir en forma de lista con nombres (ej. umbral = list(calor = "t_max > 20"))

Corregir algunas erratas en la documentación de olas: "mbral", "drtn"

Corregir algunas erratas en la documentación de umbrales: "⁠N (numérico) ocurrencia del evento *⁠ ``prop' "

olas(datos$fecha, heatwave = datos$t_max > 20) #Ok

umbrales(heatwave = datos$t_max > 20) #Ok

Corregí los errores de tipeo (drtn es un tipo de dato que representa diferencia entre tiempos). No cambié el comportamiento de las funciones olas() y umbrales(), siguen aceptando nombres de columnas como argumentos. Si bien entiendo el comentario de @VeruGHub, cuando desarrollamos estas funciones lo tuvimos en cuenta y decidimos mantener este comportamiento ya que hay muchas otras funciones de R que lo usan (por ejemplo data.frame() o summarise()). Por esto y teniendo en cuenta que las funciones no tienen otros argumentos definidos, consideramos que no hay riesgo de sobrescritura.

pdsi/pdsi_ac

No funciona indicandole los argumentos no definidos por defecto.

Habría que incluir en la documentación las referencias que se citan e incluir las unidades de la precipitacion y la etp. Es redundante incluir un argumento para los coeficientes en forma de lista y luego argumentos para cada coeficiente por separado. Sugiero eliminar una de las dos formas.

Corregir algunas erratas de la ayuda: "list de coeficientes", "El el cálculo".


pdsi(precipitacion = c(500, 600), #No funciona
     etp = c(1000, 1200)) 

pdsi_ac(precipitacion = c(500, 600), #No funciona
     etp = c(1000, 1200)) 

Extendí la documentación con todo lo mencionado. Además incluí la necesidad de tener una serie larga de datos para calcular el indice (no es suficiente con 2 datos). La función pdsi_coeficientes() ahora se documenta aparte para no generar confución.

scale_color/fill_inta pdsi/pdsi_ac

El ejemplo no funciona. Tampoco funciona otro ejemplo que he creado siguiendo las instrucciones de la documentación.

Corregir erratas en la documentación: "escala.)", "tiene prioridad por sobre", "deinidos"


prueba <- metadatos_nh()
prueba$altura <- as.numeric(prueba$altura)

ggplot(prueba) + #No funciona
  geom_point(aes(x = lon, y = lat, color = altura)) +
  scale_color_inta(escala = escala_pp_diaria)

set.seed(496) #El ejemplo no funciona
datos_aleatorios <- subset(metadatos_nh())
datos_aleatorios <- data.frame(datos_aleatorios,
pp = rgamma(nrow(datos_aleatorios), 0.5, scale = 1)*70)

ggplot(datos_aleatorios, aes(lon, lat)) +
geom_contour(aes(z = pp, fill = after_stat(level))) +
scale_fill_inta(escala = escala_pp_diaria)

Mejoré la documentación de esta familia de funciones, en primer lugar indicando que son escalas de colores para variables discretas (no estaba explicito) y mejoré los ejemplos usando un nuevo set de datos con valores de precipitación mensual. Esto simplifica los ejemplos y son más realistas.

spei_indice/spi_indice

La funcionalidad está repetida respecto a las funciones spi y spei del paquete SPEI. A menos que se aclare qué aporta esta función sugiero eliminarla y crear una viñeta que incluya cómo se usa el paquete SPEI con datos en el formato del INTA. Si he entendido bien, las fechas deben introducirse como meses ¿no? Por favor, especificar esto si se mantiene la función. En la documentación no se entiende bien que es el argumento "escalas" y no se cita el paquete SPEI.

prueba <- data.frame(fecha = seq(as.Date("1985-01-01"), as.Date("2015-12-01"), by = "1 month"))
set.seed(42)
prueba$pp <- rgamma(nrow(prueba), shape = 2, scale = 10)

spi_indice(prueba$fecha, prueba$pp, escalas = 1)

#Misma funcionalidad que en paquete SPEI:
spi_agro <- spi_indice(fecha = prueba$fecha, precipitacion = prueba$pp, escalas = 1)

spi_spei <- SPEI::spi(data = prueba$pp, scale = 1)

identical(spi_agro$spi, as.vector(spi_spei$fitted))

Modifiqué la documentaición incluyendo lo solicitado. Una de las ventajas principales de estas funciones y que justifican su existencia es que devuelven el resultado como un data.frame que se puede usar de manera directa para el análisis de datos con dplyr.

Además de lo mencionado anteriormente:

Espero no haberme olvidado de nada, aguardo sus comentarios!

Pakillo commented 2 weeks ago

Muchas gracias @paocorrales !

Estoy revisando todos los cambios. Me parece que la página web del paquete no está actualizada a la última versión, ¿cierto? Veo que el GitHub Actions de pkgdown está desactivado. Por favor actívelo o regenere la web con pkgdown manualmente. Nos ayudaría ver la web actualizada para la revisión

¡Gracias!

paocorrales commented 2 weeks ago

La web ahora está actualizada y la GitHub action activada. Gracias!

Pakillo commented 2 weeks ago

¡Muchas gracias @paocorrales por la revisión! Creo que el paquete y la documentación han mejorado muchísimo. Algunos comentarios a la última versión:

decil

Es posible que no esté entendiendo bien esta función, pero encuentro su valor de salida un poco confuso. Por ejemplo, en este caso, decil devuelve el decil (de 1 a 10) en que se encuentra cada valor concreto de la variable:

x = 1:10
decil(x)
##  [1]  1  2  3  4  5  6  7  8  9 10

Pero en este caso devuelve valores con decimales desde 0 hasta 10:

x <- seq(0, 19.5, 0.5)
decil(x)
##  [1]  0.25  0.50  0.75  1.00  1.25  1.50  1.75  2.00  2.25  2.50  2.75  3.00
## [13]  3.25  3.50  3.75  4.00  4.25  4.50  4.75  5.00  5.25  5.50  5.75  6.00
## [25]  6.25  6.50  6.75  7.00  7.25  7.50  7.75  8.00  8.25  8.50  8.75  9.00
## [37]  9.25  9.50  9.75 10.00

Pienso que estos valores con decimales podrían generar confusión. Por ejemplo, alguien que quiera saber en qué decil cae cada dato podría intentar hacer un round de estos valores, obteniendo esto:

round(decil(x))
##  [1]  0  0  1  1  1  2  2  2  2  2  3  3  3  4  4  4  4  4  5  5  5  6  6  6  6
## [26]  6  7  7  7  8  8  8  8  8  9  9  9 10 10 10

cuando probablemente necesite ejecutar

ceiling(decil(x))
##  [1]  1  1  1  1  2  2  2  2  3  3  3  3  4  4  4  4  5  5  5  5  6  6  6  6  7
## [26]  7  7  7  8  8  8  8  9  9  9  9 10 10 10 10

¿Cuál es la ventaja de devolver números decimales como salida de decil? ¿Pienso que tal vez sería deseable devolver el decil como número entero comprendido entre 1 y 10? Creo que @VeruGHub apuntaba hacia la misma sugerencia en su revisión.

@VeruGHub también planteaba esta cuestión: ¿Por qué para calcular la anomalia es necesario especificar que no se consideren los NA, pero para decil no? Entiendo que la función decil ignora los NA por defecto, pero tal vez eso debería ser una opción que eligiera el usuario (esto es, añadir un argumento na.rm).

Por ejemplo, la función quantile da un error si le pedimos los deciles de un vector con NA:

x <- 1:10
x[1] <- NA
quantile(x, probs = seq(0, 1, 0.1))
## Error in quantile.default(x, probs = seq(0, 1, 0.1)): missing values and NaN's not allowed if 'na.rm' is FALSE

Si especificamos na.rm = TRUE entonces sí funciona:

quantile(x, probs = seq(0, 1, 0.1), na.rm = TRUE)
##   0%  10%  20%  30%  40%  50%  60%  70%  80%  90% 100% 
##  2.0  2.8  3.6  4.4  5.2  6.0  6.8  7.6  8.4  9.2 10.0

Al igual que @VeruGHub, pienso que la función decil podría/debería incorporar este argumento por consistencia con quantile y anomalia_porcentual y para mayor control de los usuarios. Si no, la función decil siempre arroja un resultado que podría ser poco robusto si la variable contiene muchos NA.

dias_promedio

Esta función ha quedado muy bien, y la documentación y los ejemplos son muy claros. ¡Gracias!

umbrales

El último ejemplo requiere cambiar de summarise a reframe.

olas

Creo que sería bueno incluir un argumento na.rm o similar para evitar inconsistencias y que el usuario pueda elegir qué hacer con los datos faltantes.

Por ejemplo, esta ola de calor de 22 días de duración:

data(NH0358)
NH0358 %>%
  reframe(olas(fecha, calor = t_max > 20)) %>%
  slice_head(n = 1)
##     ola     inicio        fin longitud
## 1 calor 1951-01-01 1951-01-22  22 days

se transforma en otra de 19 días de duración si falta el dato de temperatura del 20 de enero:

NH0358 %>%
  mutate(t_max = if_else(fecha == "1951-01-20", NA, t_max)) |>
  reframe(olas(fecha, calor = t_max > 20)) %>%
  slice_head(n = 1)
##     ola     inicio        fin longitud
## 1 calor 1951-01-01 1951-01-19  19 days

pero vuelve a tener 22 días de duración si falta la fila completa del 20 de enero (cuando en realidad ese día la temperatura pudo ser \<20ºC):

NH0358 %>%
  filter(fecha != "1951-01-20") |> 
  reframe(olas(fecha, calor = t_max > 20)) %>%
  slice_head(n = 1)
##     ola     inicio        fin longitud
## 1 calor 1951-01-01 1951-01-22  22 days

Pienso que para evitar posibles inconsistencias, podría añadirse un nuevo argumento na.rm o similar, de manera que los usuarios puedan elegir qué hacer con los datos faltantes, y que el resultado sea consistente. Si el usuario quiere que los NA se tengan en cuenta, pienso que la función debería devolver NA como fecha de inicio, fin, o duración, cuando sea imposible determinar el día concreto de inicio o fin (por haber datos faltantes).

Por otro lado, es una decisión de los autores, pero como sugerencia, ¿tal vez utilizar “duración” en lugar de “longitud” para especificar la duración (en días) de la ola?

ith

El ejemplo devuelve todos los valores de ith como NA, dado que hr es NA en las 10 primeras filas. Tal vez podemos filtrar primero:

NH0358 %>%
  filter(!is.na(hr)) |> 
  mutate(t_media = (t_max + t_min)/2) %>%
  mutate(ith = ith(t_media, hr)) %>%
  slice_head(n = 10)
##    codigo codigo_nh      fecha t_max t_min precip lluvia_datos lluvia llovizna
## 1      50      0358 1960-01-01  24.4  10.3    0.0            0      0        0
## 2      50      0358 1960-01-02  26.2  11.8    0.0            0      0        0
## 3      50      0358 1960-01-03  30.3  10.8    0.0            0      0        0
## 4      50      0358 1960-01-04  33.5  16.4    0.0            0      0        0
## 5      50      0358 1960-01-05  29.1  19.0    7.7            1      1        1
## 6      50      0358 1960-01-06  31.4  15.3    0.0            0      0        0
## 7      50      0358 1960-01-07  31.9  16.1    0.0            0      0        0
## 8      50      0358 1960-01-08  33.8  17.4    0.0            0      0        0
## 9      50      0358 1960-01-09  32.8  18.3    0.0            0      0        0
## 10     50      0358 1960-01-10  33.1  18.9    0.0            0      0        0
##    granizo nieve t_aire_max t_aire_min t_suelo_max t_suelo_min heliofania_efec
## 1        0     0        7.7         NA          NA          NA            11.5
## 2        0     0        8.0         NA          NA          NA             6.4
## 3        0     0       11.4         NA          NA          NA            10.7
## 4        0     0       11.1         NA          NA          NA            10.8
## 5        0     0       15.9         NA          NA          NA             0.2
## 6        0     0       13.2         NA          NA          NA            11.0
## 7        0     0       11.5         NA          NA          NA            11.5
## 8        0     0       13.9         NA          NA          NA            11.1
## 9        0     0       15.4         NA          NA          NA            11.4
## 10       0     0       14.1         NA          NA          NA            11.3
##    heliofania_rel p_vapor hr   td rocio viento_10m viento_2m  rad etp t_media
## 1              80     9.2 56  5.8     0          3       2.4 27.9 4.8   17.35
## 2              45    11.1 56  8.5     0          2       1.6 19.1 3.9   19.00
## 3              75    10.4 45  7.5     0         10       8.0 26.5 5.9   20.55
## 4              75    10.7 41  8.0     0          4       3.2 26.6 5.7   24.95
## 5               1    15.7 74 13.7     0          1       0.8  8.4 2.8   24.05
## 6              77    13.0 58 10.9     0          0       0.0 26.9 5.0   23.35
## 7              80    11.4 45  8.9     1          6       4.8 27.8 6.1   24.00
## 8              78    12.3 42 10.0     0          7       5.6 27.1 6.4   25.60
## 9              80    16.1 59 14.1     0         11       8.8 27.6 6.9   25.55
## 10             79    12.4 48 10.1     0         13      10.4 27.4 7.4   26.00
##         ith
## 1  61.96434
## 2  64.21560
## 3  65.66553
## 4  70.77370
## 5  72.81753
## 6  70.32707
## 7  69.99700
## 8  71.67448
## 9  73.48225
## 10 72.85120

PDSI

El ejemplo usa with y $ para crear variables, pero en coherencia con la filosofía ‘tidy’ del paquete ¿tal vez sea preferible utilizar código más tipico del tidyverse?

datos <- data.frame(fecha = seq(as.Date("1985-01-01"), 
                                as.Date("2015-12-01"), 
                                by = "1 month"))
set.seed(42)
datos |> 
  mutate(pp = rgamma(nrow(datos), shape = 2, scale = 10),
         etp = rgamma(nrow(datos), shape = 1, scale = 3),
         pdsi_ac = pdsi(pp, etp)) |> 
  slice_head(n = 10)
##         fecha        pp       etp    pdsi_ac
## 1  1985-01-01 36.489561 4.4884365  0.6035107
## 2  1985-02-01  8.881098 2.4454927 -0.6361097
## 3  1985-03-01 15.592198 0.6427609 -0.6008359
## 4  1985-04-01  4.521825 2.1359691 -1.1479653
## 5  1985-05-01 13.728401 4.8577625 -1.2702747
## 6  1985-06-01  8.028043 3.8066111 -1.7626929
## 7  1985-07-01 49.905627 0.4327422  1.5345299
## 8  1985-08-01 14.241745 1.6288137 -0.1455719
## 9  1985-09-01  4.646857 1.4363694 -0.5643165
## 10 1985-10-01  2.812335 4.3753527 -1.2851817

(modificar también el ejemplo de pdsi_coeficientes).

SPI

Creo que sería bueno que los ejemplos de esta función siguieran el estilo utilizado en el tutorial, esto es, usando reframe en lugar de with.

He corregido algunas erratas en este pull request.

Espero que mis comentarios/sugerencias sean claros y fáciles de solucionar. Quedo a la espera de estos cambios, ¡y ya casi estaríamos!

Gracias