Diseño Completamente al azar

Paquete - Agricolae - Ggplot2

Planteamiento

Evaluar el contenido de nitrógeno (mg/planta) en plantas de maíz.

Exploración visual de los datos

library(readr)
rhizobium <- read_csv("C:/Users/juanc/Downloads/rhizobium.csv")
## Parsed with column specification:
## cols(
##   Tmt = col_character(),
##   peso = col_double()
## )

Con ggplot, vamos a reordenar los datos por peso, en orden ascendente.

library(ggplot2)
rhz<-ggplot(rhizobium, aes(x= reorder(Tmt, peso), y= peso, fill=Tmt)) 

rhz + geom_boxplot()+ xlab("TRATAMIENTO") + ylab("PESO")

Realizamos un ANOVA

dcar<-aov(peso~Tmt,data=rhizobium)
summary(dcar)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## Tmt          5  847.0  169.41   14.37 1.48e-06 ***
## Residuals   24  282.9   11.79                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

La prueba de F (P<0.01) indica que al menos uno de los tratamientos es diferente. Es decir, el contenido de nitrógeno en las cepas de maíz es diferente en al menos una de las cepas (= tratamientos).

Aplicamos prueba de Tukey

La prueba de Tukey nos permite comparar las medias individuales provenientes de un análisis de varianza de varias muestras sometidas a tratamientos distintos.

TukeyHSD(dcar, "Tmt", ordered = TRUE)
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
##     factor levels have been ordered
## 
## Fit: aov(formula = peso ~ Tmt, data = rhizobium)
## 
## $Tmt
##                   diff         lwr       upr     p adj
## 3DOk4-3DOk13      1.38 -5.33416704  8.094167 0.9870716
## composite-3DOk13  5.44 -1.27416704 12.154167 0.1621550
## 3DOk7-3DOk13      6.66 -0.05416704 13.374167 0.0527514
## 3DOk5-3DOk13     10.72  4.00583296 17.434167 0.0006233
## 3DOk1-3DOk13     15.56  8.84583296 22.274167 0.0000029
## composite-3DOk4   4.06 -2.65416704 10.774167 0.4434643
## 3DOk7-3DOk4       5.28 -1.43416704 11.994167 0.1852490
## 3DOk5-3DOk4       9.34  2.62583296 16.054167 0.0029837
## 3DOk1-3DOk4      14.18  7.46583296 20.894167 0.0000128
## 3DOk7-composite   1.22 -5.49416704  7.934167 0.9926132
## 3DOk5-composite   5.28 -1.43416704 11.994167 0.1852490
## 3DOk1-composite  10.12  3.40583296 16.834167 0.0012341
## 3DOk5-3DOk7       4.06 -2.65416704 10.774167 0.4434643
## 3DOk1-3DOk7       8.90  2.18583296 15.614167 0.0048849
## 3DOk1-3DOk5       4.84 -1.87416704 11.554167 0.2617111

R nos ayudó a calcular la diferencia de medias, para que las medias sean diferentes el intervalo no debe contener cero. Asimismo, se indica el valor de p asociado a dicha diferencia.

Nota:

El análisis de varianza aov supone que nuestro diseño es equilibrado, lo que quiere decir que el número de repeticiones es el mismo en todos los tratamientos.

Aplicando prueba de Tukey con la librería Agricolae

Recordar que “Tmt” es la variable que corresponde a los tratamientos.

HSD <- HSD.test(dcar, "Tmt", group = TRUE)
with(rhizobium, HSD)
## $statistics
##    MSerror Df     Mean       CV      MSD
##   11.78867 24 19.88667 17.26515 6.714167
## 
## $parameters
##    test name.t ntr StudentizedRange alpha
##   Tukey    Tmt   6         4.372651  0.05
## 
## $means
##            peso      std r  Min  Max  Q25  Q50  Q75
## 3DOk1     28.82 5.800172 5 19.4 33.0 27.0 32.1 32.6
## 3DOk13    13.26 1.427585 5 11.6 14.4 11.8 14.2 14.3
## 3DOk4     14.64 4.116188 5  9.1 19.4 11.9 15.8 17.0
## 3DOk5     23.98 3.777168 5 17.7 27.9 24.3 24.8 25.2
## 3DOk7     19.92 1.130044 5 18.6 21.0 18.8 20.5 20.7
## composite 18.70 1.601562 5 16.9 20.8 17.3 19.1 19.4
## 
## $comparison
## NULL
## 
## $groups
##            peso groups
## 3DOk1     28.82      a
## 3DOk5     23.98     ab
## 3DOk7     19.92     bc
## composite 18.70     bc
## 3DOk4     14.64      c
## 3DOk13    13.26      c
## 
## attr(,"class")
## [1] "group"

Si las medias son iguales tienen la misma letra o diferentes letras si son diferentes las medias. Por lo cual los tratamientos 3DOK4 y 3DOk13, composite y 3DOK7 tienen la misma media. El coeficiente de variación (CV) nos indica que hay una dispersión aceptable de los datos (Distribución homogénea)

Representación gráfica de los tratamientos

Utilizaremos una gráfica de puntos

dotr <- ggplot(rhizobium, aes(x= reorder(Tmt, -peso), y= peso)) + geom_dotplot(binaxis = "y", stackdir = "center")

dotr + xlab('TRATAMIENTO') + ylab('PESOS (mg)') + theme_bw()
## `stat_bindot()` using `bins = 30`. Pick better value with `binwidth`.

Con esta gráfica se observa claramente como el tratamiento 3DOK1 posee un mayor peso en cuanto contenido de nitrógeno se refiere. Encontramos que los valores de los pesos menos dispersos se encuentra en el tratamiento 3DOK7 con valores cercanos a la media (19.88) mientras que los más dispersos se encuentran en los tratamientos 3DOK4 y 3DOK1.

Análisis de residuales

Es importante verificar los residuales para detectar posibles anomalias. Un análisis visual es mediante la generación de gráficas con la opción plot(resultados_anova), donde resultados_anova es el objeto que contiene los resultados del análisis de varianza. Se presentan aquí solamente se presentan las que se consideran más importantes, aquellas que contienen los residuales estudentizados.

Gráfica de residuales estandarizados vs. valores teóricos

plot(dcar)

## hat values (leverages) are all = 0.2
##  and there are no factor predictors; no plot no. 5