Área bajo la curva (AUDPC) en R.

El área bajo la curva de progreso de la enfermedad (AUDPC), es un medida cuantitativa utilizada para calcular la intensidad de una enfermedad a lo largo del tiempo (APS, s.f). Su cálculo, se basa en el método trapezoidal, donde se calcula la intensidad media de la enfermedad entre cada par de puntos temporales adyacentes (Madden et al. 2007), y se calcula mediante la ecuación:

\[AUDPC = \sum_{i=1}^{n}\left(\frac{y_i+y_{i+1}}{2}\right) (t_{i+1} - t_i)\]

Donde \(y_{i}\) corresponde al dato medido del progreso de la enfermedad, \(y_{i+1}\), a la siguiente medición; \(t_{i}\) al tiempo de inoculación y \(t_{i+1}\) el tiempo de la medición posterior.

Cálculo de la AUDPC en R

El siguiente ejercicio, (tomado de la American Phytopathological Society (APS)), se calcula la AUDPC, con un set de datos.

Las mediciones del progreso de la enfermedad se realizaron en 4 días distintos:

Mediciones

ds0<-1
ds1<-2
ds2<-7
ds3<-7.5

disease.severity<-c(ds0,ds1,ds2,ds3)

A partir de ellos, se crea una varible, llamada disease.severity, donde mediante un vector se guarda la información. Esta variable, corresponde a las mediciones realizadas cada día.

Tiempo

t0<-0
t1<-2
t2<-5
t3<-6

time.period<-c(t0,t1,t2,t3)

Se crea una variable con el tiempo, que va a ser fijo. Si se tiene más de un tratamiento y se deben hacer varios cálculos, esta varible es fija para cada uno de ellos, dado que corresponde a los días en los que se realizaron las mediciones. Usualmente son los días después de inoculación (ddi).

Cálculo AUDPC

Con estas dos variables disease.severity y time.period, se procede a realizar el cálculo de la AUDPC

library(agricolae)
audpc(disease.severity,time.period, type="absolute")
## evaluation 
##      23.75

Cálculo de AUDPC para data set.

Si se quiere calcular la AUDPC para varios tratamientos, la función que se necesita modificar de acuerdo a los datos a trabajar es disease.severity.

Para ello, se tiene el siguiente ejemplo adaptado de Mendiburu,(s.f).

library(agricolae)
data(disease)
View(disease)

Se llama la liberia agricolae, la cual ya tiene preguardada un set de datos (disease). Este, cuenta con 7 columnas, que para el ejercicio, las de interes corresponden a la columna 4, 5 y 6 (E2, E5 y E7).

time.period<-c(14,21,28)

En la varible fija, se ponen los ddi que se realizaron las mediciones.

disease.severity<-disease[,c(4,5,6)]; disease.severity
##    E2 E5 E7
## 1   4  9 10
## 2   3  6  7
## 3   3  4  6
## 4   0  0  2
## 5   0  0  2
## 6   2  2  4
## 7   7 14 23
## 8   6 12 14
## 9   3  7 14
## 10  2  3  5
## 11  2  4  8
## 12  0  0  2
## 13  1  2  4
## 14  6 16 21
## 15  7  8 12
## 16  8 11 14
## 17  3  5 13
## 18  5 12 15
## 19  2  9 15
## 20  3  6 12
## 21 11 21 35

Se crea una varible, con un set de datos que incluyan unicamente las mediciones de interés.

[,c(4,5,6)] este fragmento de código, llama a todas las filas del data set, y unicamente a las columnas de interés.

[nrow,ncol].

  • nrow: Filas, ncol: Columnas.
  • ncol: c(4,5,6) -> Se crea un vector con las columnas de interés para crear el data set.

Cálculo de AUDPC con todos los datos

library(agricolae)
results<-audpc(disease.severity, time.period, type="absolute")
results 
##  [1] 112.0  77.0  59.5   7.0   7.0  35.0 203.0 154.0 108.5  45.5  63.0   7.0
## [13]  31.5 206.5 122.5 154.0  91.0 154.0 122.5  94.5 308.0

Tabla con los resultados

library(agricolae)
trat<-disease[,c(3)];trat #Columna de Tratamientos
##  [1] T1 T2 T3 T4 T5 T6 T0 T1 T2 T3 T4 T5 T6 T0 T1 T2 T3 T4 T5 T6 T0
## Levels: T0 T1 T2 T3 T4 T5 T6
dataaudpc<-data.frame(trat,results) #Data ser de Tratamientos y AUDPC
View(dataaudpc)

Se crea un data set con los datos de AUDPC obtenidos.

  • trat: Corresponde a la variable con las columna de los tratamientos, del data set original.
  • results: La variable con los resultados del cálculo de AUDPC.
  • data.frame: función que crea el data set, con las columnas de interés.

Código completo

library(agricolae)
data(disease)
View(disease)
time.period<-c(14,21,28)
disease.severity<-disease[,c(4,5,6)]; disease.severity
##    E2 E5 E7
## 1   4  9 10
## 2   3  6  7
## 3   3  4  6
## 4   0  0  2
## 5   0  0  2
## 6   2  2  4
## 7   7 14 23
## 8   6 12 14
## 9   3  7 14
## 10  2  3  5
## 11  2  4  8
## 12  0  0  2
## 13  1  2  4
## 14  6 16 21
## 15  7  8 12
## 16  8 11 14
## 17  3  5 13
## 18  5 12 15
## 19  2  9 15
## 20  3  6 12
## 21 11 21 35
results<-audpc(disease.severity, time.period, type="absolute")
results
##  [1] 112.0  77.0  59.5   7.0   7.0  35.0 203.0 154.0 108.5  45.5  63.0   7.0
## [13]  31.5 206.5 122.5 154.0  91.0 154.0 122.5  94.5 308.0
trat<-disease[,c(3)];trat #Columna de Tratamientos
##  [1] T1 T2 T3 T4 T5 T6 T0 T1 T2 T3 T4 T5 T6 T0 T1 T2 T3 T4 T5 T6 T0
## Levels: T0 T1 T2 T3 T4 T5 T6
dataaudpc<-data.frame(trat,results) #Data ser de Tratamientos y AUDPC
View(dataaudpc)

Gráfico de los datos

Se reestructuran los datos, para la construcción del gráfico, que permita una mejor visualización y comprensión de los datos analizados, donde para cada tratamiento se calcula el promedio por días de medición, para observar mejor su comportamiento.

library(dplyr)
## 
## Adjuntando el paquete: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats   1.0.1     ✔ readr     2.1.6
## ✔ ggplot2   4.0.2     ✔ stringr   1.6.0
## ✔ lubridate 1.9.5     ✔ tibble    3.3.1
## ✔ purrr     1.2.1     ✔ tidyr     1.3.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
data(disease)
time_points <- c("d14", "d21", "d28")
colnames(disease)[4:6] <- as.character(time_points)
View(disease)


disease1 <- disease %>%
  pivot_longer(cols = c("d14", "d21", "d28"), 
               names_to = "Etapa", 
               values_to = "Valor") %>%
  # Agrupamos para sacar el promedio de las repeticiones por tratamiento y etapa
  group_by(trt, Etapa) %>%
  summarise(Media = mean(Valor), .groups = 'drop')
View(disease1)

ggplot(disease1, aes(x = Etapa, y = Media, color = trt, group = trt)) +
  geom_line(size = 1) +      # Dibuja las líneas
  geom_point(size = 2) +     # Añade puntos en cada etapa
  theme_minimal() +          # Un estilo limpio
  labs(title = "Evolución de Tratamientos a lo largo del Tiempo",
       x = "Etapa de Medición",
       y = "Valor Promedio",
       color = "Tratamiento") +
  scale_color_brewer(palette = "Set1") # Colores más fáciles de distinguir
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Para observar las diferencias de la AUDPC cálculada entre tratamientos, se construye un boxplot a partir de los datos de la tabla dataaudpc, donde se compara cada tratamiento, con el resultado obtenido como se muestra a continuación:

# Graficar las curvas de progreso
ggplot(dataaudpc, aes(x = trat, y = results, fill = trat)) +
  geom_boxplot(alpha = 0.7) +
  labs(
    title = "Comparación del AUDPC entre tratamientos",
    x = "Tratamiento",
    y = "AUDPC"
  ) +
  theme_minimal() +
  theme(legend.position = "none")

Análisis de datos

Si se quisiera realizar el análisis estadístico clásico frecuentista de los resultados y datos anteriormente gráficados, se necesita realizar el análisis de varianza (ANOVA), para ver si existen diferencias entre los tratamientos de acuerdo a la severidad presentada, en el tiempo estudiado.

Para ello, se realiza primero la prueba ANOVA, para ver si existen diferencias entre los tratamientos

dataaudpc$trat <- as.factor(dataaudpc$trat)
modelo <- aov(results ~ trat, data = dataaudpc)
summary(modelo)
##             Df Sum Sq Mean Sq F value Pr(>F)   
## trat         6  81926   13654   5.531  0.004 **
## Residuals   14  34561    2469                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Se observa con el resultado de la prueba, que si existen diferencias entre los tratamientos por lo que se procede a realizar la comprobación de supuestos.

Prueba de Normalidad

residuos <- residuals(modelo)
shapiro.test(residuos)
## 
##  Shapiro-Wilk normality test
## 
## data:  residuos
## W = 0.91401, p-value = 0.06593

Aunque el valor es muy cercano a 0.05, de acuerdo a los rangos establecidos, no se rechaza la hipótesis nula, por lo qué, de acuerdo con el resultado, se dice que los datos se comportan de manera normal (Distribución de Gauss).

Prueba de Homocedasticidad

bartlett.test(results ~ trat, data = dataaudpc)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  results by trat
## Bartlett's K-squared = 4.3308, df = 6, p-value = 0.632
library(car)
## Cargando paquete requerido: carData
## 
## Adjuntando el paquete: 'car'
## The following object is masked from 'package:purrr':
## 
##     some
## The following object is masked from 'package:dplyr':
## 
##     recode
leveneTest(results ~ trat, data = dataaudpc)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  6   0.295 0.9293
##       14

Ambas pruebas que buscan determinar si existen homogeneidad de varianzas, muestran que no hay heterogeneidad en las varianzas.

Al cumplir ambos supuestos, se pueden realizar las pruebas de Post-hoc, en este caso la prueba de Tukey

library(agricolae)

posthoc<-HSD.test(modelo, "trat", group = TRUE)
posthoc
## $statistics
##    MSerror Df Mean       CV      MSD
##   2468.667 14  103 48.23852 138.5236
## 
## $parameters
##    test name.t ntr StudentizedRange alpha
##   Tukey   trat   7         4.828954  0.05
## 
## $means
##      results      std r       se   Min   Max    Q25   Q50    Q75
## T0 239.16667 59.63710 3 28.68604 203.0 308.0 204.75 206.5 257.25
## T1 129.50000 21.85749 3 28.68604 112.0 154.0 117.25 122.5 138.25
## T2 113.16667 38.71154 3 28.68604  77.0 154.0  92.75 108.5 131.25
## T3  65.33333 23.30415 3 28.68604  45.5  91.0  52.50  59.5  75.25
## T4  74.66667 74.19119 3 28.68604   7.0 154.0  35.00  63.0 108.50
## T5  45.50000 66.68396 3 28.68604   7.0 122.5   7.00   7.0  64.75
## T6  53.66667 35.40598 3 28.68604  31.5  94.5  33.25  35.0  64.75
## 
## $comparison
## NULL
## 
## $groups
##      results groups
## T0 239.16667      a
## T1 129.50000     ab
## T2 113.16667     ab
## T4  74.66667      b
## T3  65.33333      b
## T6  53.66667      b
## T5  45.50000      b
## 
## attr(,"class")
## [1] "group"

Donde se observan diferencias entre los tratamientos con respecto al control (T0), siendo menor la severidad entre los tratamientos T3, T4, T5, y T6. Los resultados y los grupos anteriores se visualizan en la última gráfica.

library(agricolae)
posthoc<-HSD.test(modelo, "trat", group = TRUE)
posthoc
## $statistics
##    MSerror Df Mean       CV      MSD
##   2468.667 14  103 48.23852 138.5236
## 
## $parameters
##    test name.t ntr StudentizedRange alpha
##   Tukey   trat   7         4.828954  0.05
## 
## $means
##      results      std r       se   Min   Max    Q25   Q50    Q75
## T0 239.16667 59.63710 3 28.68604 203.0 308.0 204.75 206.5 257.25
## T1 129.50000 21.85749 3 28.68604 112.0 154.0 117.25 122.5 138.25
## T2 113.16667 38.71154 3 28.68604  77.0 154.0  92.75 108.5 131.25
## T3  65.33333 23.30415 3 28.68604  45.5  91.0  52.50  59.5  75.25
## T4  74.66667 74.19119 3 28.68604   7.0 154.0  35.00  63.0 108.50
## T5  45.50000 66.68396 3 28.68604   7.0 122.5   7.00   7.0  64.75
## T6  53.66667 35.40598 3 28.68604  31.5  94.5  33.25  35.0  64.75
## 
## $comparison
## NULL
## 
## $groups
##      results groups
## T0 239.16667      a
## T1 129.50000     ab
## T2 113.16667     ab
## T4  74.66667      b
## T3  65.33333      b
## T6  53.66667      b
## T5  45.50000      b
## 
## attr(,"class")
## [1] "group"
###Preparación letras de Tukey
# Extraer letras
letras <- posthoc$groups
# Pasar a data frame y conservar el nombre del tratamiento
letras$trat <- rownames(letras)
# Reordenar por tratamiento (importante)
letras <- letras[, c("trat", "groups", "results")]
letras
##    trat groups   results
## T0   T0      a 239.16667
## T1   T1     ab 129.50000
## T2   T2     ab 113.16667
## T4   T4      b  74.66667
## T3   T3      b  65.33333
## T6   T6      b  53.66667
## T5   T5      b  45.50000
library(dplyr)
pos_letras <- dataaudpc |>
  group_by(trat) |>
  summarise(y = max(results)) |>
  left_join(letras, by = "trat")

Gráfica final

library(ggplot2)

ggplot(dataaudpc, aes(x = trat, y = results, fill = trat)) +
  geom_boxplot(alpha = 0.7) +
  geom_text(
    data = pos_letras,
    aes(x = trat, y = y + 10, label = groups),
    size = 5
  ) +
  labs(
    title = "Comparación del AUDPC entre tratamientos",
    x = "Tratamiento",
    y = "AUDPC"
  ) +
  theme_minimal() +
  theme(legend.position = "none")

Nota: Se realizaron actualizaciones al respecto documento en febrero de 2026.

Bibliografía

American Phytophatology Society, (APS). (2024). Audpc. AUDPC; American Phytophatology Society. https://www.apsnet.org/edcenter/disimpactmngmnt/topc/EcologyAndEpidemiologyInR/DiseaseProgress/Pages/AUDPC.aspx

Mendiburu, F. (s.f.). AUDPC function—Rdocumentation. Rdocumentation. Recuperado 5 de septiembre de 2024, de https://www.rdocumentation.org/packages/agricolae/versions/1.3-7/topics/audpc