Parcial Diseño de Experimentos
Jesús Ortiz
5/11/2020
Nombres: Cristian Danilo Romero Orjuela (cromeroo@unal.edu.co) Jesús David Ortiz Melo (jeortizm@unal.edu.co)
Primer Punto
- En un estudio conducido en ambiente controlado se tuvieron 72 macetas, cada una con una planta a la que a cierta edad se le midió el contenido de clorofila (índice de clorofila) con un sensor (SPAD). El total de macetas se correspondió con 9 tratamientos asociados a estrés hídrico. Se sabe que la varianza de las 72 observaciones es 833.
ANOVA
| Suma de Cuadrados | Grados de Libertad | Cuadrado Medio | Fo | |
|---|---|---|---|---|
| Entre los tratamientos | 6000 | 8 | 750 | 0.8891 |
| Error( Dentro de los tratamientos) | 53143 | 63 | 843.53 | |
| Total | 59143 | 71 |
Datos que poseemos :
\[Varianza_{total}=\frac{SC_{total}}{n-1}\]
Al poseer la varianza procedemos a multiplicar con con la cantidad de datos menos 1 arrojandonos el valor de SC=59143.
Teniendo dicho valor procedemos por calculos aritmeticos simples encontrar el Error( Dentro de los tratamientos) para luego hallar el cuadrado medio.
Una vez que poseemos todos los datos encontramos el Fcalculado para compararlo contra el tabulado
\[ F_{0}=\frac{Ms_{trat}}{Ms_{e}} \]
Si el F tabulado es 2.8. ¿qué puede decirse acerca de la Hipótesis nula de igualdad de los promedios del índice en todas las condiciones de tratamiento (use el p valor así como el cociente F calculado de la tabla para concluir)?
R: Encontrando un F calculado de 0.8891 siendo muy inferior al tabulado de 2.8, aceptando la Ho,entonces no se rechaza la Hipotesis Nula con lo que concluimos que los tratamientos tienen similitud en sus medias. Ademas al ser menor a 1 indica que la variabilidad entre las repeticiones es mayor que la variabilidad entre tratamientos
¿Vale la pena comparar las medias de tratamientos a posteriori del ANOVA (prueba de Tukey)?
R:No usamos la prueba de Tukey pues no encontramos una diferencia signifcativa y al ser esta prueba util para generar comparaciones entre los pares de los tratamientos no resulta siendo util al tener similitud entre estos.
- De ser afirmativa su respuesta aplique este procedimiento y concluya (seleccione el o los tratamientos con mayores contenidos de clorofila). Para esto genere unos datos para cada tratamiento que tengan una varianza como la dada en el ejercicio y una media que desde el primer al último tratamiento cambie a razón de 10 unidades del índice, empezando con 403.
set.seed(2020) #Semilla
data <- c(rnorm(8, 403, sqrt(833)),
rnorm(8, 413, sqrt(833)),
rnorm(8, 423, sqrt(833)),
rnorm(8, 433, sqrt(833)),
rnorm(8, 443, sqrt(833)),
rnorm(8, 453, sqrt(833)),
rnorm(8, 463, sqrt(833)),
rnorm(8, 473, sqrt(833)),
rnorm(8, 483, sqrt(833))) #Generación de datos con rnormdata1 <- as.data.frame(matrix(data, 8, 9, byrow = F))
colnames(data1) <- c("Data1", "Data2", "Data3", "Data4", "Data5","Data6", "Data7", "Data8", "Data9") #Generación de etiquetas
tratamiento <- gl(9, 8, 72, c('Tratameinto1', 'Tratameinto2', 'Tratameinto3', 'Tratameinto4', 'Tratameinto5', 'Tratameinto6','Tratameinto7', 'Tratameinto8', 'Tratameinto9')) # Generación de etiquetasboxplot(data ~ tratamiento) # Comparación medias de tratamientos
anova <- aov(data ~ tratamiento) # Desarrollo del ANOVA
summary(anova)## Df Sum Sq Mean Sq F value Pr(>F)
## tratamiento 8 67900 8488 8.359 1.18e-07 ***
## Residuals 63 63966 1015
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Revise los tres supuestos del ANOVA discutidos en clase y con este valide el uso del ANOVA en este estudio (con los datos generados). De no cumplirse el supuesto de homocedasticidad utilice el ANOVA de Welch discutido en clase y vuelva a concluir.
Revision de supuestos
- Prueba de normalidad del vector de resiudales del modelo.
aov_residuals = residuals(object = anova)
hist(aov_residuals)
shapiro.test(x = aov_residuals) # Prueba de Shapiro-Wilk para normalidad##
## Shapiro-Wilk normality test
##
## data: aov_residuals
## W = 0.98276, p-value = 0.4291Se observa que los residuales tienen una distribucion normal ya que p-value (0.4291) es mayor a 0.05.
plot(anova,2)
- Prueba de homogeneidad de varianzas
bartlett.test(aov_residuals ~ tratamiento)##
## Bartlett test of homogeneity of variances
##
## data: aov_residuals by tratamiento
## Bartlett's K-squared = 22.007, df = 8, p-value = 0.004903En este caso encontramos aplicando la prueba de bartlett.test de la comparación de los residuales que los grupos son heteroscedáticos debido a que el p-valor es menor a 0,05 por lo que no tenemos homocedasticidad, refiriéndose a que la varianza no es constante. Mostrando una alta dispersión en los residuales de los datos.
Por lo cual usaremos el ANOVA de Welch compara dos medias para ver si son iguales. Es una alternativa al ANOVA clásico y se puede utilizar incluso si sus datos violan el supuesto de homogeneidad de las variaciones.
anova_welch = oneway.test(data ~ tratamiento)
anova_welch##
## One-way analysis of means (not assuming equal variances)
##
## data: data and tratamiento
## F = 10.723, num df = 8.000, denom df = 25.981, p-value = 1.583e-06Para nuestros resultados, el p valor muy bajo indica que estos resultados son estadísticamente significativos. Nuestra prueba de muestra proporciona suficientes pruebas para concluir que las medias de todos los grupos no son iguales en los tratamientos.
plot(anova,1)
- Independencia de los residuales.
resid = anova$residuals
plot(resid, pch = 20)
No se observa un patron en los residuales por lo cual se concluye que son independientes.
Segundo Punto
Antes de hilar el algodón, éste debe ser procesado para eliminar las materias extrañas y la humedad. El limpiador de pelusas más común es el limpiador de pelusas tipo sierra de batería controlada. Aunque el limpiador de pelusas de motor de sierra (M1) es uno de los más efectivos, también es uno de los limpiadores que causa más daño a las fibras de algodón. Un investigador del algodón diseñó un estudio para comparar cuatro alternativas de limpieza de las fibras de algodón: M2, M3, M4 y M5. Los métodos M2 y M3 son mecánicos, mientras que los métodos M4 y M5 son una combinación mecánica y química. El investigador quiso tener en cuenta el impacto de los diferentes cultivadores en el proceso y, por lo tanto, obtuvo fardos de algodón de seis diferentes granjas algodoneras. Las granjas fueron consideradas como bloques en el estudio. Después de una limpieza preliminar del algodón, los seis fardos fueron mezclados a fondo, y luego fue procesada una igual cantidad de algodón por cada uno de los cinco métodos de limpieza de pelusas. Las pérdidas en peso (en kg) después de la limpieza las fibras de algodón se dan en la siguiente tabla. Durante el procesamiento de las muestras de algodón, lasa mediciones de la granja 1 procesada por el limpiador M1 se perdieron.
- Realice el ANOVA para este diseño recordando que es un caso desbalanceado. Concluya sobre el resultado de la tabla del ANOVA obtenida. (¿Afecta el orden de colocación de los efectos del modelo dentro del software R?, Verifique si la tabla del ANOVA cambia): en R: lm() lm(respuesta~ bloque+método) lm(respuesta~ método + bloque)
library(readxl)
Algodon <- read_excel("C:/Users/yisus/Downloads/Algodon.xlsx") # Se reorganizó la tabla
#View(Algodon)
head(Algodon)## # A tibble: 6 x 3
## Metodo Granjero Perdida
## <chr> <dbl> <chr>
## 1 M1 1 <NA>
## 2 M1 2 6.75
## 3 M1 3 13.05
## 4 M1 4 10.26
## 5 M1 5 8.01
## 6 M1 6 8.42tail(Algodon)## # A tibble: 6 x 3
## Metodo Granjero Perdida
## <chr> <dbl> <chr>
## 1 M5 1 9.27
## 2 M5 2 4.39
## 3 M5 3 13.44
## 4 M5 4 9.13
## 5 M5 5 9.20
## 6 M5 6 7.13names(Algodon) # Etiquetas## [1] "Metodo" "Granjero" "Perdida"# Reorganización de los datos
Talgodon = as.data.frame(Algodon)
Metodof = factor(Talgodon$Metodo)
Granjerof = factor(Talgodon$Granjero)
Perdidaf = as.vector(Talgodon$Perdida)
Perdidaf1 = as.numeric(Perdidaf)
str(Talgodon)## 'data.frame': 30 obs. of 3 variables:
## $ Metodo : chr "M1" "M1" "M1" "M1" ...
## $ Granjero: num 1 2 3 4 5 6 1 2 3 4 ...
## $ Perdida : chr NA "6.75" "13.05" "10.26" ...Usamos el ANOVA Tipo III ya que permite calcular la suma de cuadrados ajustando cada tratamiento teniendo en cuenta cualquier otro tratamiento que no lo contenga y de forma independiente de cualquier tratamiento que lo contenga. Estas sumas de cuadrados no se alteran por las variaciones del tamaño muestral de las casillas, de modo que son útiles para modelos no equilibrados.
library(car)## Warning: package 'car' was built under R version 4.0.3## Loading required package: carData## Warning: package 'carData' was built under R version 4.0.3Modelo = lm(Perdidaf1 ~ Granjerof + Metodof) # Modelo lineal
ANOVA = Anova(Modelo, type = 3) # Desarrollo del ANOVA desbalanceado mediante el uso de anova(... ,type = 3).
ANOVA## Anova Table (Type III tests)
##
## Response: Perdidaf1
## Sum Sq Df F value Pr(>F)
## (Intercept) 198.62 1 143.8751 2.610e-10 ***
## Granjerof 139.66 5 20.2331 5.163e-07 ***
## Metodof 49.12 4 8.8951 0.0003186 ***
## Residuals 26.23 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Modelo2 = lm(Perdidaf1 ~ Metodof + Granjerof) # Modelo lineal
ANOVA2 = Anova(Modelo2, type = 3) # Desarrollo del ANOVA desbalanceado mediante el uso de anova(... ,type = 3).
ANOVA2## Anova Table (Type III tests)
##
## Response: Perdidaf1
## Sum Sq Df F value Pr(>F)
## (Intercept) 198.62 1 143.8751 2.610e-10 ***
## Metodof 49.12 4 8.8951 0.0003186 ***
## Granjerof 139.66 5 20.2331 5.163e-07 ***
## Residuals 26.23 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1R: El orden de colocación de los efectos del modelo no afecta el resultados final.
- Estimar el valor de la observación usando el promedio de los datos para los cinco granjeros del mismo método M1 y luego realice el análisis de varianza para probar las diferencias en las pérdidas medias de peso para los cinco métodos de limpiado de las fibras de algodón. Compare este resultado con el caso desbalanceado (de ser posible).
Obtención del promedio de los datos para los cinco granjeros del mismo método M1
library(dplyr)## Warning: package 'dplyr' was built under R version 4.0.3##
## Attaching package: 'dplyr'## The following object is masked from 'package:car':
##
## recode## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionTalgodon2 = filter(Talgodon, Metodo == "M1") # Usa la libreria "dplyr" con la funcion filter para aislar los granjeros del metodo M1
Talgodon2 = na.omit(Talgodon2) # Se elimino el dato faltante para que no interfiriera en el calculo del promedio
Talgodon2## Metodo Granjero Perdida
## 2 M1 2 6.75
## 3 M1 3 13.05
## 4 M1 4 10.26
## 5 M1 5 8.01
## 6 M1 6 8.42# calculo del promedio
PromAOV = (Talgodon2$Perdida)
PromAOV2 = as.numeric(PromAOV)
Prom = mean(PromAOV2)
Prom## [1] 9.298Talgodon3 = Talgodon %>% mutate(Perdida = replace(Perdida, which(is.na (Perdida)), Prom)) # Remplazo del valor faltante por el promedio obtenido
Talgodon3## Metodo Granjero Perdida
## 1 M1 1 9.298
## 2 M1 2 6.75
## 3 M1 3 13.05
## 4 M1 4 10.26
## 5 M1 5 8.01
## 6 M1 6 8.42
## 7 M2 1 5.54
## 8 M2 2 3.53
## 9 M2 3 11.20
## 10 M2 4 7.21
## 11 M2 5 3.24
## 12 M2 6 6.45
## 13 M3 1 7.67
## 14 M3 2 4.15
## 15 M3 3 9.79
## 16 M3 4 8.27
## 17 M3 5 6.75
## 18 M3 6 5.50
## 19 M4 1 7.89
## 20 M4 2 1.97
## 21 M4 3 8.97
## 22 M4 4 6.12
## 23 M4 5 4.22
## 24 M4 6 7.84
## 25 M5 1 9.27
## 26 M5 2 4.39
## 27 M5 3 13.44
## 28 M5 4 9.13
## 29 M5 5 9.20
## 30 M5 6 7.13Metodof3 = factor(Talgodon3$Metodo)
Granjerof3 = factor(Talgodon3$Granjero)
Perdida3 = as.vector(Talgodon3$Perdida)
Perdidaf3 = as.numeric(Perdida3)Modelof3 = lm(Perdidaf3 ~ Granjerof3 + Metodof3) # Modelo lineal
ANOVAf3 = aov(Modelof3) # Desarrollo del ANOVA
summary(ANOVAf3)## Df Sum Sq Mean Sq F value Pr(>F)
## Granjerof3 5 139.36 27.873 21.015 2.45e-07 ***
## Metodof3 4 51.15 12.787 9.641 0.000164 ***
## Residuals 20 26.53 1.326
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Modelo = lm(Perdidaf1 ~ Granjerof + Metodof)
ANOVA = Anova(Modelo, type = 3)
ANOVA## Anova Table (Type III tests)
##
## Response: Perdidaf1
## Sum Sq Df F value Pr(>F)
## (Intercept) 198.62 1 143.8751 2.610e-10 ***
## Granjerof 139.66 5 20.2331 5.163e-07 ***
## Metodof 49.12 4 8.8951 0.0003186 ***
## Residuals 26.23 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Por medio de la comparación del ANOVA con promedio y el ANOVA desbalanceado se comprobó la efectividad del ANOVA tipo 3 al dar unos valores muy cercanos, mostrando pequeñas diferencias en los decimales, sin embargo la utilización del ANOVA tipo 3 dará resultados más adecuados.
- Probemos nuestro ingenio y capacidad de trabajo en R: Proponga un promedio basado en la medida de M1 con los datos disponibles, pero sumando una cantidad delta de modo que el coeficiente de variación entre las seis mediciones sea menor al 20%.
# Formula de la Covarianza con variables.
Funcionglob = function(st) {
x = st[1]
y = st[2]
((sqrt(((x-9.298-y)^2 + (6.75-9.298-y)^2 + (13.05-9.298-y)^2 + (10.26-9.298-y)^2 + (8.01-9.298-y)^2 + (8.42-9.298-y)^2)/5 ))/9.298 + y)
}# Funcion optimize para hallar los valores que permitan que la covarianza sea inferior a 0.2
Optimi = optimize(Funcionglob, lower = 0, upper = 2)## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive value
## Warning in optimize(Funcionglob, lower = 0, upper = 2): NA/Inf replaced by
## maximum positive valueOptimi## $minimum
## [1] 1.999959
##
## $objective
## [1] NATercer Punto
Use la función de R para generar de la distribución uniforme unos datos de carbono orgánico del suelo medido a 5 cm y 10 cm de profundidad. Suponga que la medida de la capa superior osciló entre 3.0 y 3.3+0.1 y de la capa inferior osciló entre 2 y 2.0+0.2. Use expand.grid para generar una ventana de observación de 0 a 100 m para la longitud y de 0 a 200 m para la latitud. Genere 50 datos en cada capa. Use la función sort.int de R para ordenar los datos de cada capa con la opción partial=25+3 dentro de la propia función sort.int. Una vez cree los datos realice algún diagrama de color (preferiblemente 3D) que permita visualizar las medidas de carbono en cada capa generadas por computadora. Compare si se encuentran diferencias en la media de carbono entre capas utilizando un nivel de confianza del 95%.
set.seed(2020)
Suelo5 = runif(50,3,3.4) # Generación de datos uniformes para la capa superior
Suelo5## [1] 3.258761 3.157690 3.247401 3.190756 3.054439 3.026954 3.051661 3.157247
## [9] 3.001033 3.248082 3.305766 3.297534 3.330466 3.169092 3.163715 3.215877
## [17] 3.384289 3.261423 3.218686 3.106425 3.078718 3.031148 3.327357 3.376962
## [25] 3.353690 3.066351 3.142041 3.299238 3.180380 3.222350 3.385629 3.028588
## [33] 3.382324 3.379191 3.000475 3.146912 3.004416 3.373543 3.209299 3.088795
## [41] 3.394064 3.090132 3.345575 3.353873 3.249936 3.228621 3.188338 3.280267
## [49] 3.319174 3.066298set.seed(2020)
Suelo10 = runif(50,2,2.2) # Generación de datos uniformes para la capa inferior
Suelo10## [1] 2.129381 2.078845 2.123700 2.095378 2.027219 2.013477 2.025831 2.078624
## [9] 2.000517 2.124041 2.152883 2.148767 2.165233 2.084546 2.081858 2.107939
## [17] 2.192144 2.130711 2.109343 2.053213 2.039359 2.015574 2.163679 2.188481
## [25] 2.176845 2.033176 2.071020 2.149619 2.090190 2.111175 2.192815 2.014294
## [33] 2.191162 2.189596 2.000238 2.073456 2.002208 2.186772 2.104649 2.044397
## [41] 2.197032 2.045066 2.172788 2.176937 2.124968 2.114311 2.094169 2.140133
## [49] 2.159587 2.033149xy = expand.grid(x=0:sqrt(100), y=0:sqrt(200)) # Ventana de Observación
capa_su = sort.int(Suelo5,partial = 28) # Se ordeno los datos
capa_in = sort.int(Suelo10, partial = 28)CO = c(capa_in, capa_su)
plot(xy, col = CO, pch=20)
Gra = plot(xy,cex= CO,pch=20)
names(xy)## [1] "x" "y"library(scatterplot3d) # Grafica 3D
xn = xy$x[1:100]
yn = xy$y[1:100]
s3d = scatterplot3d(xn, yn, CO, color = CO, pch = 20, xlab = "x", ylab = "y", zlab = "Carbono Orgánico")
library(plot3D)
scatter3D(xn, yn, CO, ticktype = "detailed", pch = 20, cex =0.8, clab = c("Carbono Orgánico", ""), cex = 0.2, zlab = "CO", bty = "g")
Con el fin de comparar si se encuentran diferencias en la media de carbono entre capas utilizando un nivel de confianza del 95% se aplicara T-test para muestras pareada.
t.test(Suelo10, Suelo5, paired = T, mu = 0, conf.level = 0.95)##
## Paired t-test
##
## data: Suelo10 and Suelo5
## t = -126.82, df = 49, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.12191 -1.08691
## sample estimates:
## mean of the differences
## -1.10441\[H_o=\mu Suelo10-\mu Suelo5 =0 \\ H_a=\mu Suelo10-\mu Suelo5 < 0\]
Al arrojarnos un p-valor menor a 0.05 (2.2e-16) se rechaza Ho, además el intervalo de confianza de la diferencia va desde -1.12 a -1.08 al no estar presente el valor 0, confirmamos el rechazo de la hipotesis Ho. Entonces como conclusión obtuvimos que en el capa superior existio un mayor contenido de carbono organico con respecto a la capa inferior.
Cuarto Punto
El siguiente diseño se corresponde con un factorial completo (32) en arreglo completamente al azar.
D = expand.grid(F1= c(3.25, 3.75, 4.25), F2 = c(4,5,6)) # CREA EL DISEÑO 3^2
D = rbind(D,D) # Crea la estructura para dos repeticiones por tratamiento
set.seed(2020)
D = D[order(sample(1:18)),] # Aleatoriza la estructura
class(D)## [1] "data.frame"D$biomasa = sort.int(rnorm(18, 3, 0.3), partial = 9) # Crea la respuesta
D## F1 F2 biomasa
## 2 3.75 4 2.708826
## 10 3.25 4 2.772692
## 16 3.25 6 2.143359
## 4 3.25 5 2.560519
## 13 3.25 5 2.708666
## 6 4.25 5 2.773705
## 17 3.75 6 2.770350
## 8 3.75 6 2.832470
## 14 3.75 5 2.898280
## 5 3.75 5 3.359619
## 9 4.25 6 3.054099
## 1 3.25 4 3.157896
## 12 4.25 4 3.487669
## 11 3.75 4 3.451547
## 15 4.25 5 3.016111
## 7 3.25 6 3.042156
## 3 4.25 4 3.200552
## 18 4.25 6 2.989329- Escriba (completamente especificado) el modelo del diseño.
El diseño 3^2 tiene dos factores, cada uno con tres niveles. Están presentes 3^2 = 9 combinaciones de tratamientos.
\[y_{ijk}=\mu+\tau_i+\delta_j+(\tau\delta)_{ij}+\epsilon_{ijk\\i:1\cdots 3\\j:1,3\\k:1\cdots 2}\\Condiciones~Laterales~Respectivas\]
Tenemos que \(Y_{ijk}\) es la \(ijk-esima\) biomasa para el \(i-esimo\) nivel del factor 1 (F1) y el \(j-esimo\) nivel del factor 2 (F2) y la \(k-esima\) repetición.
*\(\mu\) como media general de biomasa.
*\(\tau_i\) es el efecto del \(i-esimo\) nivel del F1.
*\(\delta_j\) es el efecto del \(j-esimo\) nivel del F2.
*\((\tau\delta)_{ij}\) es la interacción del \(i-esimo\) nivel del F1 con el \(j-esimo\) nivel del F2.
*\(\epsilon_{ijk}\) es el error residual para el \(i-esimo\) nivel del F1, el \(j-esimo\) nivel del F2 y la \(k-esima\) repetición.
** F1: Dosis de un insecticida **
** F2: Número de aplicaciones **
Hipótesis
Hipótesis 1: El efecto de la interacción entre los dos factores, es nulo para todo i,j.
\[H_{1}: (\tau\delta)_{ij}=0~\\ \forall i,j\]
Hipótesis 2: El efecto del factor 1 (Dosis de insecticida) es nulo.
\[H_{12}: \tau_i=0~\\ i=1\cdots 3\]
Hipótesis 3: El efecto del factor 2 (Aplicaciones) es nulo.
\[H_{03}: \delta_j=0~\\ j= 1,3\]
- Realice el Anova para este diseño y de ser necesario realice la prueba de comparaciones de medias para los efectos principales de F1: dosis de un insecticida que se sospecha tiene un efecto de disminución del crecimiento (biomasa) y F2: número de aplicaciones durante el desarrollo del cultivo.
D1 = as.data.frame(D)
biomasa1 = as.vector(D$biomasa)
Dosis = as.factor(D$F1)
Aplicaciones = as.factor(D$F2)str(D)## 'data.frame': 18 obs. of 3 variables:
## $ F1 : num 3.75 3.25 3.25 3.25 3.25 4.25 3.75 3.75 3.75 3.75 ...
## $ F2 : num 4 4 6 5 5 5 6 6 5 5 ...
## $ biomasa: num 2.71 2.77 2.14 2.56 2.71 ...
## - attr(*, "out.attrs")=List of 2
## ..$ dim : Named int [1:2] 3 3
## .. ..- attr(*, "names")= chr [1:2] "F1" "F2"
## ..$ dimnames:List of 2
## .. ..$ F1: chr [1:3] "F1=3.25" "F1=3.75" "F1=4.25"
## .. ..$ F2: chr [1:3] "F2=4" "F2=5" "F2=6"aov4 = aov(biomasa1 ~ Dosis * Aplicaciones)
summary(aov4)## Df Sum Sq Mean Sq F value Pr(>F)
## Dosis 2 0.4161 0.20804 1.979 0.194
## Aplicaciones 2 0.3426 0.17128 1.630 0.249
## Dosis:Aplicaciones 4 0.1635 0.04087 0.389 0.812
## Residuals 9 0.9459 0.10510Según el ANOVA no existe interacción entre las dosis del insecticida y las aplicaciones, ya que el p-valor (0.812) es mayor a 0.05, por lo tanto no se rechaza la hipótesis 1, tampoco las hipótesis 2 y 3 se rechazan a causa de que el efecto del Factor 1 y del Factor 2 son nulos.
Por otro lado se observa que los valores F para cada fuente de variación son superiores a 1, lo que indica que la principal fuente de variación son los tratamientos, es decir que las diferencia son entre tratamientos, no intratratamientos.
aov4si = aov(biomasa1~Dosis+Aplicaciones) # ANOVA sin interacción para corroborar
summary(aov4si)## Df Sum Sq Mean Sq F value Pr(>F)
## Dosis 2 0.4161 0.20804 2.438 0.126
## Aplicaciones 2 0.3426 0.17128 2.007 0.174
## Residuals 13 1.1094 0.08534Con el ANOVA sin interacción se llega al mismo resultado, confirmando lo anteriormente dicho, en la que la variacion se da entre tratamientos.Ademas se comprobo que no existe interacción entre los factores, se puede aplicar la prueba de comparación de medias, para encontrar cuales son los causantes de las diferencias.
TukeyHSD(aov4)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = biomasa1 ~ Dosis * Aplicaciones)
##
## $Dosis
## diff lwr upr p adj
## 3.75-3.25 0.27263391 -0.2499630 0.7952308 0.3549307
## 4.25-3.25 0.35602957 -0.1665674 0.8786265 0.1934891
## 4.25-3.75 0.08339567 -0.4392013 0.6059926 0.8975539
##
## $Aplicaciones
## diff lwr upr p adj
## 5-4 -0.24371366 -0.7663106 0.2788833 0.4289015
## 6-4 -0.32456967 -0.8471666 0.1980273 0.2454279
## 6-5 -0.08085601 -0.6034529 0.4417409 0.9033230
##
## $`Dosis:Aplicaciones`
## diff lwr upr p adj
## 3.75:4-3.25:4 0.11489233 -1.1676566 1.3974413 0.9999813
## 4.25:4-3.25:4 0.37881648 -0.9037325 1.6613655 0.9452628
## 3.25:5-3.25:4 -0.33070147 -1.6132505 0.9518475 0.9734324
## 3.75:5-3.25:4 0.16365519 -1.1188938 1.4462042 0.9997380
## 4.25:5-3.25:4 -0.07038587 -1.3529349 1.2121631 0.9999996
## 3.25:6-3.25:4 -0.37253651 -1.6550855 0.9100125 0.9497229
## 3.75:6-3.25:4 -0.16388379 -1.4464328 1.1186652 0.9997353
## 4.25:6-3.25:4 0.05642012 -1.2261289 1.3389691 0.9999999
## 4.25:4-3.75:4 0.26392415 -1.0186248 1.5464731 0.9930776
## 3.25:5-3.75:4 -0.44559380 -1.7281428 0.8369552 0.8826430
## 3.75:5-3.75:4 0.04876285 -1.2337861 1.3313118 1.0000000
## 4.25:5-3.75:4 -0.18527820 -1.4678272 1.0972708 0.9993609
## 3.25:6-3.75:4 -0.48742885 -1.7699778 0.7951201 0.8301399
## 3.75:6-3.75:4 -0.27877612 -1.5613251 1.0037729 0.9902626
## 4.25:6-3.75:4 -0.05847221 -1.3410212 1.2240768 0.9999999
## 3.25:5-4.25:4 -0.70951795 -1.9920669 0.5730310 0.4788847
## 3.75:5-4.25:4 -0.21516130 -1.4977103 1.0673877 0.9981920
## 4.25:5-4.25:4 -0.44920235 -1.7317513 0.8333466 0.8784883
## 3.25:6-4.25:4 -0.75135300 -2.0339020 0.5311960 0.4171999
## 3.75:6-4.25:4 -0.54270027 -1.8252493 0.7398487 0.7487731
## 4.25:6-4.25:4 -0.32239636 -1.6049453 0.9601526 0.9769758
## 3.75:5-3.25:5 0.49435666 -0.7881923 1.7769056 0.8206001
## 4.25:5-3.25:5 0.26031560 -1.0222334 1.5428646 0.9936553
## 3.25:6-3.25:5 -0.04183504 -1.3243840 1.2407139 1.0000000
## 3.75:6-3.25:5 0.16681768 -1.1157313 1.4493667 0.9996989
## 4.25:6-3.25:5 0.38712159 -0.8954274 1.6696706 0.9389879
## 4.25:5-3.75:5 -0.23404106 -1.5165900 1.0485079 0.9968189
## 3.25:6-3.75:5 -0.53619170 -1.8187407 0.7463573 0.7588975
## 3.75:6-3.75:5 -0.32753897 -1.6100880 0.9550100 0.9748241
## 4.25:6-3.75:5 -0.10723507 -1.3897840 1.1753139 0.9999890
## 3.25:6-4.25:5 -0.30215064 -1.5846996 0.9803983 0.9841808
## 3.75:6-4.25:5 -0.09349792 -1.3760469 1.1890511 0.9999962
## 4.25:6-4.25:5 0.12680599 -1.1557430 1.4093550 0.9999606
## 3.75:6-3.25:6 0.20865273 -1.0738963 1.4912017 0.9985347
## 4.25:6-3.25:6 0.42895663 -0.8535923 1.7115056 0.9008106
## 4.25:6-3.75:6 0.22030391 -1.0622451 1.5028529 0.9978778Puesto que los p-valores ajustados en todos los casos nos arroja un valor mayor a 0.05, se puede deducir que utilizar cualquiera de las dosis equivale lo mismo, pues la evidencia estadística comprueba un mismo efecto sobre la biomasa.
- Use los resultados del ANOVA y el gráfico de interacción (ggplot2) para visualizar si existe o no interacción entre los factores.
R = aggregate(x=list(media_bio=D1$biomasa), by=list(Dosis = Dosis, Aplicaciones = Aplicaciones),
FUN=mean, na.rm=TRUE)
R## Dosis Aplicaciones media_bio
## 1 3.25 4 2.965294
## 2 3.75 4 3.080186
## 3 4.25 4 3.344110
## 4 3.25 5 2.634592
## 5 3.75 5 3.128949
## 6 4.25 5 2.894908
## 7 3.25 6 2.592757
## 8 3.75 6 2.801410
## 9 4.25 6 3.021714library(ggplot2)
R %>%
ggplot(aes(Dosis,media_bio,color= Aplicaciones))+
geom_line(aes(group = Aplicaciones))
En el gráfico se aprecia que entre las aplicaciones 4 y 6 no va ha existir interacción alguna. Mientras que la aplicación 5 tendrá una minima interacción con las aplicaciones 4 y 6 cuando cambian las dosis, pero no resulta siendo significativa, por esto mismo en la prueba anova nos arroja una interacción nula. En las aplicaciones 4 y 6 la medida de biomasa parece tener un comportamiento similar cuando variamos las dosis, cosa que no se ve en el tratamiento 5. Como conlusión el intercepto que hace la aplicación 5 es despreciable, mostrando poca interacción como fue corroborado con el ANOVA Test.
- El investigador quiso colocar como covariable el contenido de arcilla(expansible) en el suelo utilizado en cada unidad experimental. Genere unos datos con la distribución uniforme cuya medida oscile entre 0.20 y 0.40 , ordene estas medidas en forma decreciente y meta dentro del análisis esta variable. Especifique nuevamente el modelo y realice el análisis de covarianza respectivo ¿se justifica el uso de la covariable? Construya nuevamente el gráfico de interacción y compare con el caso sin covariable (discuta el resultado). Revise en internet los supuestos que deben tener las covariables para ser utilizadas en el modelo. ¿Se está incumpliendo en nuestros datos alguno de los supuestos necesarios? Revise los supuestos sobre los residuales tanto del ANOVA como del ANCOVA ¿qué puede percibir? ¿recomendaría el uso de arcillas para minimizar el efecto sobre el contenido de biomasa que puede ocasionar el uso del insecticida?
A1=sort.default(runif(18,0.20,0.40),decreasing=TRUE)# Generación Data uniforme del contenido de arcilla
A1## [1] 0.3985124 0.3768873 0.3726635 0.3345261 0.3302072 0.3228664 0.3227439
## [8] 0.3153687 0.3111702 0.3060593 0.2940596 0.2848798 0.2775186 0.2590587
## [15] 0.2575276 0.2457057 0.2334734 0.2114003D2 = cbind(D1, A1)
D3 = as.data.frame(D2)- Modelo
\[y_{ijk}=\mu+\tau_i+\delta_j+(\tau\delta)_{ij}+\beta(x_{ijk}-\bar{x}) +\epsilon_{ijk\\i:1\cdots 3\\j:1,3\\k:1\cdots 2}\\\] Exploración:
*\(Y_{ijk}\) Representa la \(ijk-esima\) observación de biomasa para el \(i-esimo\) nivel del factor 1 (F1) y el \(j-esimo\) nivel del factor 2 (F2) y la \(k-esima\) repetición.
*\(\tau_i\) = Efecto del \(i-esimo\) nivel del F1.
*\(\delta_j\) Efecto del \(j-esimo\) nivel del F2.
*\(\mu\) Media de biomasa.
*\((\tau\delta)_{ij}\) Es la interacción del \(i-esimo\) nivel del F1 con el \(j-esimo\) nivel del F2.
*\(\bar{x}\) Media de los valores de \(x_{ijk}\)
*\(x_{ijk}\)Medida de la covariable que se hace para \(Y_{ijk}\)
*\(x_{ijk}\) Medida de la covariable que se hace para \(Y_{ijk}\)
*\(\epsilon_{ijk}\) Error residual de \(i-esimo\) nivel del F1, el \(j-esimo\) nivel del F2 y la \(k-esima\) repetición.
*\(\beta\) Coeficiente de regresión que relaciona \(Y_{ijk}\) con la covariable $ x_{ijk}$
ancova1=aov(biomasa1~A1+(Dosis+Aplicaciones),data=D3)
summary(ancova1)## Df Sum Sq Mean Sq F value Pr(>F)
## A1 1 0.8130 0.8130 22.233 0.000501 ***
## Dosis 2 0.1884 0.0942 2.577 0.117222
## Aplicaciones 2 0.4279 0.2139 5.850 0.016846 *
## Residuals 12 0.4388 0.0366
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#Sin InteracciónSegun los datos obtenidos, la covariable arcilla tendria un efecto significativo sobre la biomasa ya que el Pr() es menor a 0.05, por otro lado se observa tambien un posible efecto de las aplicaciones aunque sería reducido.
ancova2=aov(biomasa1~A1+Dosis*Aplicaciones,data=D3)
summary(ancova2)## Df Sum Sq Mean Sq F value Pr(>F)
## A1 1 0.8130 0.8130 17.123 0.00326 **
## Dosis 2 0.1884 0.0942 1.984 0.19961
## Aplicaciones 2 0.4279 0.2139 4.506 0.04891 *
## Dosis:Aplicaciones 4 0.0589 0.0147 0.310 0.86330
## Residuals 8 0.3798 0.0475
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#Con interacciónSe compruba que no existe interacción entre dosis y aplicaciones.
Grap = data.frame(D3$biomasa, D3$A1)
Grap## D3.biomasa D3.A1
## 1 2.708826 0.3985124
## 2 2.772692 0.3768873
## 3 2.143359 0.3726635
## 4 2.560519 0.3345261
## 5 2.708666 0.3302072
## 6 2.773705 0.3228664
## 7 2.770350 0.3227439
## 8 2.832470 0.3153687
## 9 2.898280 0.3111702
## 10 3.359619 0.3060593
## 11 3.054099 0.2940596
## 12 3.157896 0.2848798
## 13 3.487669 0.2775186
## 14 3.451547 0.2590587
## 15 3.016111 0.2575276
## 16 3.042156 0.2457057
## 17 3.200552 0.2334734
## 18 2.989329 0.2114003- Grafica relacion lineal entre Respuesta y Covariable
library(ggplot2)
ggplot(data = Grap, aes(x =D3.A1, y =D3.biomasa) ) + geom_point() + xlab("Arcilla") + ylab("Biomasa") + geom_smooth(method = "lm", se = FALSE)## `geom_smooth()` using formula 'y ~ x'
La gráfica de línea ajustada indica una relación lineal entre la arcilla y la biomasa.
Prueba de Corelación
cor.test(D3$biomasa, D3$A1)##
## Pearson's product-moment correlation
##
## data: D3$biomasa and D3$A1
## t = -3.5112, df = 16, p-value = 0.002894
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.8612947 -0.2786494
## sample estimates:
## cor
## -0.6596974El valor p obtenido es de 0.002894 que es menor a 0.05. De modo que se puede concluir que Arcilla y Biomasa estan significativamente correlacionsados con con un coeficiente de correlación de -0.65.
- ¿Se justifica el uso de la covariable?
Si se justifica el uso de covariable ya que se observa una fuerte correlacion entre la Arcilla y la Biomasa, de modo que el ANCOVA permitira aumentar la potencia estadística del diseño experimental pues reducira la variabilidad y el error.
- Construya nuevamente el gráfico de interacción y compare con el caso sin covariable (discuta el resultado).
Gnew = aggregate(x=list(media_bio2=D3$biomasa), by=list(Dosis = Dosis, Aplicaciones = Aplicaciones),
FUN=mean, na.rm=TRUE)
R## Dosis Aplicaciones media_bio
## 1 3.25 4 2.965294
## 2 3.75 4 3.080186
## 3 4.25 4 3.344110
## 4 3.25 5 2.634592
## 5 3.75 5 3.128949
## 6 4.25 5 2.894908
## 7 3.25 6 2.592757
## 8 3.75 6 2.801410
## 9 4.25 6 3.021714library(ggplot2)
R %>%
ggplot(aes(Dosis,Gnew$media_bio,color= Aplicaciones))+
geom_line(aes(group = Aplicaciones))
Aunque en el grafico puede llegar a parecer que existe una interaccion entre los tratamientos, las pruebas ANOVA mostraron que no existia una interacción significativa, de modo que el gráfico solo nos estaría diciendo que el tratamiento de 4 aplicaciones con dosis de 4.25 es mejor en terminos de biomasa. Por otro lado con el ANCOVA se demostró que habia una fuerte corelacion entre Arcilla y Biomasa, ya que su correlación nos arrojó un valor de -0.65 se dice que su relación es signidicativa, sin embargo es una relación inversamente proporcional, es decir por cada unidad que aumente una de las variables la otra disminuira.
Revise los supuestos sobre los residuales tanto del ANOVA como del ANCOVA.
Probamos Los supuestos de ANCOVA:
- Supuesto de linealidad
library(tidyverse)## Warning: package 'tidyverse' was built under R version 4.0.3## -- Attaching packages -------- tidyverse 1.3.0 --## v tibble 3.0.3 v purrr 0.3.4
## v tidyr 1.1.2 v stringr 1.4.0
## v readr 1.4.0 v forcats 0.5.0## Warning: package 'tidyr' was built under R version 4.0.3## Warning: package 'readr' was built under R version 4.0.3## Warning: package 'purrr' was built under R version 4.0.3## Warning: package 'forcats' was built under R version 4.0.3## -- Conflicts ----------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
## x dplyr::recode() masks car::recode()
## x purrr::some() masks car::some()library(ggpubr)## Warning: package 'ggpubr' was built under R version 4.0.3library(rstatix)## Warning: package 'rstatix' was built under R version 4.0.3##
## Attaching package: 'rstatix'## The following object is masked from 'package:stats':
##
## filterlibrary(broom)## Warning: package 'broom' was built under R version 4.0.3ggscatter( data = D3, x = "A1" , y = "biomasa", facet.by =c("F1", "F2"), short.panel.labs = FALSE) + stat_smooth(method = "loess", span = 1)## `geom_smooth()` using formula 'y ~ x'
Hubo una relación lineal entre la covariable Arcilla y la variable Biomasa para cada grupo, según se evaluó mediante la inspección visual de un diagrama de dispersión.
2.Homogeneidad de las pendientes de regresión.
D3 %>%
anova_test( biomasa ~ A1 + F1 + F2 +
F1*F2 + A1*F1 +
A1*F2 + A1*F2*F1)## Coefficient covariances computed by hccm()## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 A1 1 10 15.173 0.003 * 0.603
## 2 F1 1 10 0.807 0.390 0.075
## 3 F2 1 10 11.658 0.007 * 0.538
## 4 F1:F2 1 10 0.456 0.515 0.044
## 5 A1:F1 1 10 5.183 0.046 * 0.341
## 6 A1:F2 1 10 0.012 0.914 0.001
## 7 A1:F1:F2 1 10 0.402 0.540 0.039Existio homogeneidad en la prueba de pendientes de regresión ya que los terminos de interacción entre las variables de la covariable arcilla y de los factores aplicación y dosis no fueron estadísticamente significativos ya que obtuvimos un p-valor mayor a 0.005 (A1,F2,A1:F1). Comprobando que no existe una interacción entre la covariable y los dos factores, afirmando el segundo supuesto.
3.Normalidad de los residuos
residA=ancova1$residuals
normaA=shapiro.test(residA); normaA##
## Shapiro-Wilk normality test
##
## data: residA
## W = 0.93707, p-value = 0.258ifelse(0.05>normaA$p.value, "No_normal","Residuales_Normales")## [1] "Residuales_Normales"El supuesto se cumple ya que tenemos una normalidad en los residuos.
4.Igualdad de varianzas entre tratamientos
BarVar=bartlett.test(residA~interaction(Aplicaciones,Dosis)); BarVar##
## Bartlett test of homogeneity of variances
##
## data: residA by interaction(Aplicaciones, Dosis)
## Bartlett's K-squared = 9.2079, df = 8, p-value = 0.3251ifelse(0.05>BarVar$p.value,"Varianzas_Desiguales","Varianzas_Iguales")## [1] "Varianzas_Iguales"El segundo supuesto fue corroborado ya que por medio de la prueba de bartlett.test() nos arroja que tenemos varinzas iguales
- Independecia de los residuos
aov4_residuals1 = residuals(object = ancova1)
plot(aov4_residuals1, pch = 20)
No se observa un patron en los residuales por lo cual se concluye que son independientes.
Probamos los supuestos del ANOVA:
- Prueba de normalidad del vector de resiudales del modelo.
aov4_residuals = residuals(object = aov4)
hist(aov4_residuals)
shap_aov4 =shapiro.test(aov4_residuals); shap_aov4 # Prueba de Shapiro-Wilk para normalidad##
## Shapiro-Wilk normality test
##
## data: aov4_residuals
## W = 0.99029, p-value = 0.999ifelse(0.05>shap_aov4$p.value, "No_normal","Residuales_Normales")## [1] "Residuales_Normales"Se observa que los residuales tienen una distribucion normal ya que p-value (0.999) es mayor a 0.05.
plot(anova,2)
- Prueba de homogeneidad de varianzas
Bar_aov4 = bartlett.test(aov4_residuals ~ interaction(Dosis, Aplicaciones)) ; Bar_aov4##
## Bartlett test of homogeneity of variances
##
## data: aov4_residuals by interaction(Dosis, Aplicaciones)
## Bartlett's K-squared = 7.5942, df = 8, p-value = 0.4741ifelse(0.05>Bar_aov4$p.value,"Varianzas_Desiguales","Varianzas_Iguales")## [1] "Varianzas_Iguales"El segundo supuesto fue comprobado ya que por medio de la prueba de bartlett.test() nos arroja que tenemos varinzas iguales
- Independencia de los residuales.
plot(aov4_residuals, pch = 20)
No se observa un patron en los residuales por lo cual se concluye que son independientes.
¿Se está incumpliendo en nuestros datos alguno de los supuestos necesarios?
En conclusión se confirma que los tres supuestos que deben tener las covariables para ser utilizadas en el modelo.
¿Recomendaría el uso de arcillas para minimizar el efecto sobre el contenido de biomasa que puede ocasionar el uso del insecticida?
No se recomendaria para minimizar el efecto sobre el contenido de biomasa que puede ocasionar el uso del insecticida ya que la covariable arcilla tiene una correlación de -0,65, lo que indica que su contenido es inversamente proporcional al de la biomasa, por lo que si se agrega arcilla tambien se estaria reduciendo la biomasa generando un impacto negativo.
Quinto Punto
Existe un tipo de diseño anidado (factorial incompleta) conocido como anidado escalonado (staggered nested design) y ocurre tal como se muestra en la imagen, donde se tienen dos fincas sembradas con variedades de papa solo que la finca A permite que se desarrollen las dos variedades mientras que la altitud de la finca B solo permite el desarrollo de una de ellas. Además, se tienen dos parcelas con la variedad 1 en la primera finca y solo una en el resto de las fincas.
Para desarollar este ejercicio tomamos en cuenta que corresponde a un diseño anidado, factorial incompleto.
Definimos Terminos:
Tomamos A: como niveles del factor superior. Finca A y B. a: 20 y a-1: 19
El diseño anidado escalonado se utiliza para estimar las fuentes de variabilidad estimando la varianza atribuible a las parcelas, a las parcelas anidadas en las fincas, y a la variedad dentro de la finca.
Tenemos un diseño anidado escalonado cada uno de 4 etapas, a cada uno de ellas se le asigna una letra de la A a la D.
Una respuesta puede ser el rendimiento de la variedad de papa.
library(daewr)## Warning: package 'daewr' was built under R version 4.0.3## Registered S3 method overwritten by 'DoE.base':
## method from
## factorize.factor conf.designdata("polymer")
polymer## lot box prep test strength
## 1 1 1 1 1 9.76
## 2 1 1 1 2 9.24
## 3 1 1 2 1 11.91
## 4 1 2 1 1 9.02
## 5 2 1 1 1 10.65
## 6 2 1 1 2 7.77
## 7 2 1 2 1 10.00
## 8 2 2 1 1 13.69
## 9 3 1 1 1 6.50
## 10 3 1 1 2 6.26
## 11 3 1 2 1 8.02
## 12 3 2 1 1 7.95
## 13 4 1 1 1 8.08
## 14 4 1 1 2 5.28
## 15 4 1 2 1 9.15
## 16 4 2 1 1 7.46
## 17 5 1 1 1 7.84
## 18 5 1 1 2 5.91
## 19 5 1 2 1 7.43
## 20 5 2 1 1 6.11
## 21 6 1 1 1 9.00
## 22 6 1 1 2 8.38
## 23 6 1 2 1 7.01
## 24 6 2 1 1 8.58
## 25 7 1 1 1 12.81
## 26 7 1 1 2 13.58
## 27 7 1 2 1 11.13
## 28 7 2 1 1 10.00
## 29 8 1 1 1 10.62
## 30 8 1 1 2 11.71
## 31 8 1 2 1 14.07
## 32 8 2 1 1 14.56
## 33 9 1 1 1 4.88
## 34 9 1 1 2 4.96
## 35 9 1 2 1 4.08
## 36 9 2 1 1 4.76
## 37 10 1 1 1 9.38
## 38 10 1 1 2 8.02
## 39 10 1 2 1 6.73
## 40 10 2 1 1 6.99
## 41 11 1 1 1 5.91
## 42 11 1 1 2 5.79
## 43 11 1 2 1 6.59
## 44 11 2 1 1 6.55
## 45 12 1 1 1 7.19
## 46 12 1 1 2 7.22
## 47 12 1 2 1 5.77
## 48 12 2 1 1 8.33
## 49 13 1 1 1 7.93
## 50 13 1 1 2 6.48
## 51 13 1 2 1 8.12
## 52 13 2 1 1 7.43
## 53 14 1 1 1 3.70
## 54 14 1 1 2 2.86
## 55 14 1 2 1 3.95
## 56 14 2 1 1 5.92
## 57 15 1 1 1 4.64
## 58 15 1 1 2 5.70
## 59 15 1 2 1 5.96
## 60 15 2 1 1 5.88
## 61 16 1 1 1 5.94
## 62 16 1 1 2 6.28
## 63 16 1 2 1 4.18
## 64 16 2 1 1 5.24
## 65 17 1 1 1 9.50
## 66 17 1 1 2 8.00
## 67 17 1 2 1 11.25
## 68 17 2 1 1 11.14
## 69 18 1 1 1 10.93
## 70 18 1 1 2 12.16
## 71 18 1 2 1 9.51
## 72 18 2 1 1 12.71
## 73 19 1 1 1 11.95
## 74 19 1 1 2 10.58
## 75 19 1 2 1 16.79
## 76 19 2 1 1 13.08
## 77 20 1 1 1 4.34
## 78 20 1 1 2 5.45
## 79 20 1 2 1 7.51
## 80 20 2 1 1 5.21
## 81 21 1 1 1 7.60
## 82 21 1 1 2 6.72
## 83 21 1 2 1 6.51
## 84 21 2 1 1 6.35
## 85 22 1 1 1 5.12
## 86 22 1 1 2 5.85
## 87 22 1 2 1 6.31
## 88 22 2 1 1 8.74
## 89 23 1 1 1 5.28
## 90 23 1 1 2 5.73
## 91 23 1 2 1 4.53
## 92 23 2 1 1 5.07
## 93 24 1 1 1 5.44
## 94 24 1 1 2 5.38
## 95 24 1 2 1 4.35
## 96 24 2 1 1 7.04
## 97 25 1 1 1 3.50
## 98 25 1 1 2 3.88
## 99 25 1 2 1 2.57
## 100 25 2 1 1 3.76
## 101 26 1 1 1 4.80
## 102 26 1 1 2 4.46
## 103 26 1 2 1 3.48
## 104 26 2 1 1 3.18
## 105 27 1 1 1 5.35
## 106 27 1 1 2 6.39
## 107 27 1 2 1 4.38
## 108 27 2 1 1 5.50
## 109 28 1 1 1 3.09
## 110 28 1 1 2 3.19
## 111 28 1 2 1 3.79
## 112 28 2 1 1 2.59
## 113 29 1 1 1 5.30
## 114 29 1 1 2 4.72
## 115 29 1 2 1 4.39
## 116 29 2 1 1 6.13
## 117 30 1 1 1 7.09
## 118 30 1 1 2 7.82
## 119 30 1 2 1 5.96
## 120 30 2 1 1 7.14Se creó el data.frame con los primeros 80 datos, cambiando el nombre de las variables
Parcelas =gl(20,4,80)
Finca = rep(c(1, 1, 1, 2), times = 20)
Variedad = c(1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1)
test = c(1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1,1,2,1,1)
Rendimiento = c(10.76, 7.14, 10.11, 8.12, 10.15, 7.27, 9.00, 11.29, 13.50, 4.26, 8.12, 7.95, 8.08, 5.28, 9.15,7.46, 7.84, 5.91, 7.43, 6.11, 9.00, 8.38, 7.01, 8.58, 11.81, 12.88, 9.23, 10.00, 14.62, 12.31, 14.07, 14.56, 14.88, 8.96, 4.08, 4.76, 9.38, 8.02, 6.73, 6.99, 5.91, 5.79, 6.59, 6.55, 11.19, 5.31, 7.77, 6.23, 7.93, 6.48, 3.12, 5.3, 8.70, 12.6, 3.25, 3.92, 4.74, 5.60, 5.96, 5.88, 5.94, 6.28, 4.18, 5.24, 9.80, 11.00, 8.5, 11.14, 5.93, 12.16, 9.51, 12.71, 11.95, 10.58, 16.79, 13.08, 4.34, 5.45, 7.51, 4.11)
P <- data.frame(Parcelas, Finca, Variedad, test, Rendimiento)
P## Parcelas Finca Variedad test Rendimiento
## 1 1 1 1 1 10.76
## 2 1 1 1 2 7.14
## 3 1 1 2 1 10.11
## 4 1 2 1 1 8.12
## 5 2 1 1 1 10.15
## 6 2 1 1 2 7.27
## 7 2 1 2 1 9.00
## 8 2 2 1 1 11.29
## 9 3 1 1 1 13.50
## 10 3 1 1 2 4.26
## 11 3 1 2 1 8.12
## 12 3 2 1 1 7.95
## 13 4 1 1 1 8.08
## 14 4 1 1 2 5.28
## 15 4 1 2 1 9.15
## 16 4 2 1 1 7.46
## 17 5 1 1 1 7.84
## 18 5 1 1 2 5.91
## 19 5 1 2 1 7.43
## 20 5 2 1 1 6.11
## 21 6 1 1 1 9.00
## 22 6 1 1 2 8.38
## 23 6 1 2 1 7.01
## 24 6 2 1 1 8.58
## 25 7 1 1 1 11.81
## 26 7 1 1 2 12.88
## 27 7 1 2 1 9.23
## 28 7 2 1 1 10.00
## 29 8 1 1 1 14.62
## 30 8 1 1 2 12.31
## 31 8 1 2 1 14.07
## 32 8 2 1 1 14.56
## 33 9 1 1 1 14.88
## 34 9 1 1 2 8.96
## 35 9 1 2 1 4.08
## 36 9 2 1 1 4.76
## 37 10 1 1 1 9.38
## 38 10 1 1 2 8.02
## 39 10 1 2 1 6.73
## 40 10 2 1 1 6.99
## 41 11 1 1 1 5.91
## 42 11 1 1 2 5.79
## 43 11 1 2 1 6.59
## 44 11 2 1 1 6.55
## 45 12 1 1 1 11.19
## 46 12 1 1 2 5.31
## 47 12 1 2 1 7.77
## 48 12 2 1 1 6.23
## 49 13 1 1 1 7.93
## 50 13 1 1 2 6.48
## 51 13 1 2 1 3.12
## 52 13 2 1 1 5.30
## 53 14 1 1 1 8.70
## 54 14 1 1 2 12.60
## 55 14 1 2 1 3.25
## 56 14 2 1 1 3.92
## 57 15 1 1 1 4.74
## 58 15 1 1 2 5.60
## 59 15 1 2 1 5.96
## 60 15 2 1 1 5.88
## 61 16 1 1 1 5.94
## 62 16 1 1 2 6.28
## 63 16 1 2 1 4.18
## 64 16 2 1 1 5.24
## 65 17 1 1 1 9.80
## 66 17 1 1 2 11.00
## 67 17 1 2 1 8.50
## 68 17 2 1 1 11.14
## 69 18 1 1 1 5.93
## 70 18 1 1 2 12.16
## 71 18 1 2 1 9.51
## 72 18 2 1 1 12.71
## 73 19 1 1 1 11.95
## 74 19 1 1 2 10.58
## 75 19 1 2 1 16.79
## 76 19 2 1 1 13.08
## 77 20 1 1 1 4.34
## 78 20 1 1 2 5.45
## 79 20 1 2 1 7.51
## 80 20 2 1 1 4.11# Creamos el árbol para el diseño anidado escalonado.
library(collapsibleTree) ## Warning: package 'collapsibleTree' was built under R version 4.0.3collapsibleTree(P,hierarchy=c("Finca","Variedad","test"))Mostrando dos sistemas de finca el inicial subdivido en 2 variedades la primera que contiene dos repeticiones y la ultima con solo una l y la ultimo en una variedad con una sola repetición como se observa en la grafica generada por collapsible Tree
Estimamos los componentes de la varianza por medio del método REML, para ello se usa la función (lmer) del paquete lm4.
library(lme4)## Warning: package 'lme4' was built under R version 4.0.3## Loading required package: Matrix##
## Attaching package: 'Matrix'## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack## Registered S3 methods overwritten by 'lme4':
## method from
## cooks.distance.influence.merMod car
## influence.merMod car
## dfbeta.influence.merMod car
## dfbetas.influence.merMod car##
## Attaching package: 'lme4'## The following object is masked from 'package:daewr':
##
## cakemodelo3 <- lmer(Rendimiento ~ 1 + (1|Parcelas) + (1|Parcelas:Finca)+
(1|Parcelas:Finca:Variedad), data = P)## boundary (singular) fit: see ?isSingularsummary(modelo3)## Linear mixed model fit by REML ['lmerMod']
## Formula: Rendimiento ~ 1 + (1 | Parcelas) + (1 | Parcelas:Finca) + (1 |
## Parcelas:Finca:Variedad)
## Data: P
##
## REML criterion at convergence: 388.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.80461 -0.51179 -0.06751 0.56309 2.88979
##
## Random effects:
## Groups Name Variance Std.Dev.
## Parcelas:Finca:Variedad (Intercept) 0.000 0.000
## Parcelas:Finca (Intercept) 0.000 0.000
## Parcelas (Intercept) 4.478 2.116
## Residual 5.334 2.309
## Number of obs: 80, groups:
## Parcelas:Finca:Variedad, 60; Parcelas:Finca, 40; Parcelas, 20
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 8.3275 0.5391 15.45
## convergence code: 0
## boundary (singular) fit: see ?isSingularEl valor de REML para el componente de varianza para la finca dentro de la parcela es igual a 0, por lo tanto, no hay varianza atribuible a las parcelas anidadas en las fincas.
#Calculamos el porcentaje de la variación de las Parcelas
VariacionP= (100*7.0186)/(7.0186+1.2361+0.8743)# Sacar el porcentaje
VariacionP## [1] 76.88246#Del 100% de la variacón total el 76% se debe a las parcelasMostrando unas fuertes relaciones entre la variaciones de los tratamientos, que puede deberse a la varianza entre parcelas.
#Porcentaje de la variación de las variedades dentro de la finca
Variabilidad_variedad= (100*1.2361)/(7.0186+1.2361+0.8743)
Variabilidad_variedad## [1] 13.54037Como conclusión la mayor variabilidad la causan las parcelas, los resultados nos arrojan una variación del 76,88% como fue explicada anteriormente, el valor de 13,54% interpretamos que posiblemente corresponda a la variabilidad de las variedades dentro de la finca. Mostrando que la variabilidad dentro de parcela o de finca a finca es muy pequeña casí insignifcante, en este caso 0.
Entonces, si el objetivo es disminuir la variabilidad en el rendimiento del producto , deberian intentar disminuir de alguna forma la variación de los factores que cambian entre parcelas, como puede ser la parte química del suelo u otros factores vinculados a este.
Sexto Punto
En el enlace https://cran.r-project.org/web/packages/asbio/asbio.pdf se tienen unos datos de potasio de muestras de suelos medidas en 8 diferentes laboratorios. Compare descriptivamente (medidas, tablas y gráficos) para representar los datos. ¿qué prueba me recomendaría para comparar la medida que usted seleccione? Proponga una solución. Sabiendo que son muestras mezcladas de una misma finca, ¿ Se perciben diferencias en las medidas como consecuencia probable de los laboratorios? Sugerencia: Use el enfoque no paramétrico considerado en clase y su respectiva prueba de comparación por pares (Nemenyi)
Datos de potasio en 8 laboratorios diferentes.
library(asbio)## Warning: package 'asbio' was built under R version 4.0.3## Loading required package: tcltk##
## Attaching package: 'asbio'## The following object is masked from 'package:broom':
##
## bootstrapdata(K)
K## K lab
## 1 296 B
## 2 260 B
## 3 341 B
## 4 359 B
## 5 323 B
## 6 321 B
## 7 287 B
## 8 413 B
## 9 335 B
## 10 315 D
## 11 330 D
## 12 326 D
## 13 354 D
## 14 266 D
## 15 348 D
## 16 343 D
## 17 284 D
## 18 324 D
## 19 351 E
## 20 302 E
## 21 395 E
## 22 357 E
## 23 400 E
## 24 187 E
## 25 376 E
## 26 283 E
## 27 198 E
## 28 327 F
## 29 354 F
## 30 308 F
## 31 274 F
## 32 324 F
## 33 305 F
## 34 347 F
## 35 297 F
## 36 305 F
## 37 326 G
## 38 301 G
## 39 316 G
## 40 312 G
## 41 297 G
## 42 280 G
## 43 300 G
## 44 319 G
## 45 286 G
## 46 218 H
## 47 280 H
## 48 241 H
## 49 226 H
## 50 243 H
## 51 199 H
## 52 205 H
## 53 225 H
## 54 227 H
## 55 338 I
## 56 303 I
## 57 341 I
## 58 311 I
## 59 355 I
## 60 269 I
## 61 284 I
## 62 279 I
## 63 339 I
## 64 359 J
## 65 318 J
## 66 313 J
## 67 352 J
## 68 334 J
## 69 356 J
## 70 342 J
## 71 299 J
## 72 353 Jdata(K)Iniciamos el analisis Grafico
lattice::bwplot(K$K~K$lab)
A partir de la generación gráfica encontramos unos mayores valores de potasio en el apartado E, aunque muestra igualmente una mayor dispersión entre datos, es decir una mayor variabilidad entre datos teniendo unos atípicos más notables en comparación con los demás. El lab G creemos que son los datos más homogéneos y consistentes.Como apartado adicional encontramos en el H una media muy inferior a comparación de los demás, un menor rendimiento .
Realizamos un plot de densidad, útil en observar la forma de la distribución.
#Gráfico de densidad
library(ggplot2)
ggplot(data=K, aes(x= K, fill=lab)) +
geom_density(alpha=0.8)
A partir de la gráfica de densidad encontramos que la distribución de los datos en el laboratorio E presenta un mayor grosor de densidad, pues sus datos no están encasillados en un intervalo, como ya lo explicamos anteriormente cuenta con una mayor variabilidad, desde 200-400. Parece ultimo el lab G tenemos datos más concentrados una rango exacto de 275 a 320.Por último Parece que el intervalo de la densidad de F el más uniforme con respecto a los demás
Realizamos multiples pruebas estadisticas para poder inferir mejor
mean<-tapply(K$K, K$lab,mean);mean## B D E F G H I J
## 326.1111 321.1111 316.5556 315.6667 304.1111 229.3333 313.2222 336.2222des = tapply(K$K, list(K$lab), sd)
des## B D E F G H I J
## 44.40564 29.25510 80.35097 25.04995 15.37404 23.89037 31.39577 21.60890median<- tapply(K$K, K$lab,median);median## B D E F G H I J
## 323 326 351 308 301 226 311 342medias <- tapply(K$K, list(K$lab), FUN = mean)
medias## B D E F G H I J
## 326.1111 321.1111 316.5556 315.6667 304.1111 229.3333 313.2222 336.2222#Cálculo del coeficiente de variación para verificar si en verdad el laboratorio J presenta los mayores valores de potasioComenzando el análisis encontramos que el laboratorio J presenta valores más homogéneos, lo que indica más confiabilidad a los resultados de este laboratorio siendo el más idóneo si se desea elegir . Además el laboratorio E es el que mayor grado de dispersión presenta con respecto a la media, un valor poco recomendable para el tratamientos. los demas datos presentan datos similares de estos no es importante generar un analisis exhaustivo
coef_var = 100 * des/medias
coef_var## B D E F G H I J
## 13.616722 9.110586 25.382896 7.935570 5.055402 10.417315 10.023481 6.426969#Ayuda a descartar los promediosCon medias distintas usamos Kruskal- Wallis para comprobar si la variable respuesta es la misma en todas las poblaciones Es decir si las variabes pertenecen a una misma distribución, enn pocas palabras por medio de rangos
Prueba de Krustall-Wallis :Análisis de varianza de rangos.
\[H_0: \mu_{rangoB}=\mu_{rangoD}=\mu_{rangoE}=...=\mu_{rangoJ}\]
K$K## [1] 296 260 341 359 323 321 287 413 335 315 330 326 354 266 348 343 284 324 351
## [20] 302 395 357 400 187 376 283 198 327 354 308 274 324 305 347 297 305 326 301
## [39] 316 312 297 280 300 319 286 218 280 241 226 243 199 205 225 227 338 303 341
## [58] 311 355 269 284 279 339 359 318 313 352 334 356 342 299 353which.min(K$K)## [1] 24K$rangos <-rank(K$K)
tapply(K$rangos, K$lab, mean)## B D E F G H I J
## 42.444444 42.444444 42.777778 38.277778 30.611111 7.611111 36.888889 50.944444#Tenemos rangos promedio son muy similares en los laboratorio E,B y D. boxplot(K$rangos~K$lab,xlab='Laboratorio', ylab='Rangos', main='Rangos Vs Lab PLot')
De todos los resultados el lab H si difiere en mayor grado.
kruskal.test(K$K~K$lab) ##
## Kruskal-Wallis rank sum test
##
## data: K$K by K$lab
## Kruskal-Wallis chi-squared = 24.482, df = 7, p-value = 0.000937Con Kuskall wallis obtuvimos p_valor obtenido de 0.000937 (<0.05) rechaza la hipótesis nula de igualdad de rangos entre los laboratorios, y comprueba que los rangos de los promedios asignados difieren en los 8 laboratorios.
#Nemenyi comparación por par de medias de labs
PMCMR::posthoc.kruskal.nemenyi.test(K$K~K$lab) ## Warning in posthoc.kruskal.nemenyi.test.default(c(296, 260, 341, 359, 323, :
## Ties are present, p-values are not corrected.##
## Pairwise comparisons using Tukey and Kramer (Nemenyi) test
## with Tukey-Dist approximation for independent samples
##
## data: K$K by K$lab
##
## B D E F G H I
## D 1.0000 - - - - - -
## E 1.0000 1.0000 - - - - -
## F 0.9999 0.9999 0.9998 - - - -
## G 0.9324 0.9324 0.9222 0.9943 - - -
## H 0.0098 0.0098 0.0087 0.0397 0.2764 - -
## I 0.9993 0.9993 0.9989 1.0000 0.9984 0.0600 -
## J 0.9893 0.9893 0.9916 0.9051 0.4405 0.0003 0.8461
##
## P value adjustment method: noneSeptimo punto
Strip plot design. Diseño de parcelas divididas en bloques completos.
En un estudio llevado a cabo por agrónomos para determinar si existían diferencias importantes en el rendimiento de cafe que respondieran a la fertilización con N entre en las diferentes variedades de plantas, los principales tratamientos de la parcela eran tres variedades y 4 niveles de fertilizacion. El estudio se replicó con cuatro bloqueos y 3 parcelas, cada una dividida en 4 subparcelas y los datos recogidos para el experimento se muestran:
#Daata y planteo del problema: https://fliphtml5.com/waov/gjll/basic
set.seed(2020)
Bloques = gl(4, 14, 50)
VariedadN = gl(3, 4,50, c("V1","V2","V3"))
TratN = gl(4,1,50, c(0.2, 0.5, 1.4, 1.9))
AlturaP = rnorm(50, 98, 0.5)
Cafe = data.frame(Bloques, VariedadN, TratN, AlturaP)
Cafe## Bloques VariedadN TratN AlturaP
## 1 1 V1 0.2 98.18849
## 2 1 V1 0.5 98.15077
## 3 1 V1 1.4 97.45099
## 4 1 V1 1.9 97.43480
## 5 1 V2 0.2 96.60173
## 6 1 V2 0.5 98.36029
## 7 1 V2 1.4 98.46956
## 8 1 V2 1.9 97.88531
## 9 1 V3 0.2 98.87957
## 10 1 V3 0.5 98.05868
## 11 1 V3 1.4 97.57344
## 12 1 V3 1.9 98.45463
## 13 1 V1 0.2 98.59819
## 14 1 V1 0.5 97.81421
## 15 2 V1 1.4 97.93837
## 16 2 V1 1.9 98.90002
## 17 2 V2 0.2 98.85200
## 18 2 V2 0.5 96.48062
## 19 2 V2 1.4 96.85551
## 20 2 V2 1.9 98.02915
## 21 2 V3 0.2 99.08718
## 22 2 V3 0.5 98.54909
## 23 2 V3 1.4 98.15911
## 24 2 V3 1.9 97.96343
## 25 2 V1 0.2 98.41713
## 26 2 V1 0.5 98.09938
## 27 2 V1 1.4 98.64892
## 28 2 V1 1.9 98.46836
## 29 3 V2 0.2 97.92628
## 30 3 V2 0.5 98.05522
## 31 3 V2 1.4 97.59375
## 32 3 V2 1.9 97.62815
## 33 3 V3 0.2 98.54767
## 34 3 V3 0.5 99.21769
## 35 3 V3 1.4 98.19406
## 36 3 V3 1.9 98.14531
## 37 3 V1 0.2 97.85720
## 38 3 V1 0.5 98.03801
## 39 3 V1 1.4 97.71985
## 40 3 V1 1.9 98.22359
## 41 3 V2 0.2 98.45425
## 42 3 V2 0.5 97.74747
## 43 4 V2 1.4 97.84950
## 44 4 V2 1.9 97.63698
## 45 4 V3 0.2 97.40996
## 46 4 V3 0.5 98.12654
## 47 4 V3 1.4 97.81464
## 48 4 V3 1.9 98.01109
## 49 4 V1 0.2 98.33002
## 50 4 V1 0.5 98.24440str(Cafe)## 'data.frame': 50 obs. of 4 variables:
## $ Bloques : Factor w/ 4 levels "1","2","3","4": 1 1 1 1 1 1 1 1 1 1 ...
## $ VariedadN: Factor w/ 3 levels "V1","V2","V3": 1 1 1 1 2 2 2 2 3 3 ...
## $ TratN : Factor w/ 4 levels "0.2","0.5","1.4",..: 1 2 3 4 1 2 3 4 1 2 ...
## $ AlturaP : num 98.2 98.2 97.5 97.4 96.6 ...Las variedades de Nitrogeno fueron asignadas a las parcelas principales y los tratamientos a las subparcelas
Cafe$VariedadN = as.factor(Cafe$VariedadN)
Cafe$TratN = as.factor(Cafe$TratN)
Cafe$Bloques = as.numeric(Cafe$Bloques)
str(Cafe)## 'data.frame': 50 obs. of 4 variables:
## $ Bloques : num 1 1 1 1 1 1 1 1 1 1 ...
## $ VariedadN: Factor w/ 3 levels "V1","V2","V3": 1 1 1 1 2 2 2 2 3 3 ...
## $ TratN : Factor w/ 4 levels "0.2","0.5","1.4",..: 1 2 3 4 1 2 3 4 1 2 ...
## $ AlturaP : num 98.2 98.2 97.5 97.4 96.6 ...- Modelo del diseño :
\[ y_{ijk}=\mu+\alpha_i+\gamma_k+\eta_{ik}+\beta_j+\ (\alpha\beta)_{ij} +\epsilon_{ijk\\i:1\cdots 3\\j:1,4\\k:1\cdots 4} \] En el cual
*\(Y_{ijk}\) La altura de la planta de cafe
*\(\mu\) es la media general
*\(\alpha_i\) Efecto fijo de la variedad
*\(\eta_{ik}\) El error en parcela principal
*\(\beta_j\) Efecto producido po el Nitrogeno
*\((\alpha\beta)_{ij}\) Interacción entre ambos factores
*\(\epsilon_{ijk}\) Error en la subparcela
*\(\gamma_{ik}\) Efecto fijo del bloqueo
library(agricolae)## Warning: package 'agricolae' was built under R version 4.0.3modeloInicial<-with(Cafe,strip.plot(Bloques,VariedadN,TratN,AlturaP))##
## ANALYSIS STRIP PLOT: AlturaP
## Class level information
##
## VariedadN : V1 V2 V3
## TratN : 0.2 0.5 1.4 1.9
## Bloques : 1 2 3 4
##
## Number of observations: 50
##
## model Y: AlturaP ~ Bloques + VariedadN + Ea + TratN + Eb + TratN:VariedadN + Ec
##
## Analysis of Variance Table
##
## Response: AlturaP
## Df Sum Sq Mean Sq F value Pr(>F)
## Bloques 3 0.4028 0.13427 0.3853 0.7648
## VariedadN 2 2.2424 1.12121 3.7965 0.0860 .
## Ea 6 1.7720 0.29533 0.8474 0.5488
## TratN 3 0.9569 0.31896 1.3974 0.3057
## Eb 9 2.0543 0.22825 0.6550 0.7384
## TratN:VariedadN 6 0.7490 0.12483 0.3582 0.8966
## Ec 20 6.9700 0.34850
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 0.6 %, cv(b) = 0.5 %, cv(c) = 0.6 %, Mean = 98.06281modeloInicial## $data
## Bloques VariedadN TratN AlturaP
## 1 1 V1 0.2 98.18849
## 2 1 V1 0.5 98.15077
## 3 1 V1 1.4 97.45099
## 4 1 V1 1.9 97.43480
## 5 1 V2 0.2 96.60173
## 6 1 V2 0.5 98.36029
## 7 1 V2 1.4 98.46956
## 8 1 V2 1.9 97.88531
## 9 1 V3 0.2 98.87957
## 10 1 V3 0.5 98.05868
## 11 1 V3 1.4 97.57344
## 12 1 V3 1.9 98.45463
## 13 1 V1 0.2 98.59819
## 14 1 V1 0.5 97.81421
## 15 2 V1 1.4 97.93837
## 16 2 V1 1.9 98.90002
## 17 2 V2 0.2 98.85200
## 18 2 V2 0.5 96.48062
## 19 2 V2 1.4 96.85551
## 20 2 V2 1.9 98.02915
## 21 2 V3 0.2 99.08718
## 22 2 V3 0.5 98.54909
## 23 2 V3 1.4 98.15911
## 24 2 V3 1.9 97.96343
## 25 2 V1 0.2 98.41713
## 26 2 V1 0.5 98.09938
## 27 2 V1 1.4 98.64892
## 28 2 V1 1.9 98.46836
## 29 3 V2 0.2 97.92628
## 30 3 V2 0.5 98.05522
## 31 3 V2 1.4 97.59375
## 32 3 V2 1.9 97.62815
## 33 3 V3 0.2 98.54767
## 34 3 V3 0.5 99.21769
## 35 3 V3 1.4 98.19406
## 36 3 V3 1.9 98.14531
## 37 3 V1 0.2 97.85720
## 38 3 V1 0.5 98.03801
## 39 3 V1 1.4 97.71985
## 40 3 V1 1.9 98.22359
## 41 3 V2 0.2 98.45425
## 42 3 V2 0.5 97.74747
## 43 4 V2 1.4 97.84950
## 44 4 V2 1.9 97.63698
## 45 4 V3 0.2 97.40996
## 46 4 V3 0.5 98.12654
## 47 4 V3 1.4 97.81464
## 48 4 V3 1.9 98.01109
## 49 4 V1 0.2 98.33002
## 50 4 V1 0.5 98.24440
##
## $ANOVA
## Analysis of Variance Table
##
## Response: AlturaP
## Df Sum Sq Mean Sq F value Pr(>F)
## Bloques 3 0.4028 0.13427 0.3853 0.7648
## VariedadN 2 2.2424 1.12121 3.7965 0.0860 .
## Ea 6 1.7720 0.29533 0.8474 0.5488
## TratN 3 0.9569 0.31896 1.3974 0.3057
## Eb 9 2.0543 0.22825 0.6550 0.7384
## TratN:VariedadN 6 0.7490 0.12483 0.3582 0.8966
## Ec 20 6.9700 0.34850
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $gl.a
## [1] 6
##
## $gl.b
## [1] 9
##
## $gl.c
## [1] 20
##
## $Ea
## [1] 0.2953316
##
## $Eb
## [1] 0.2282548
##
## $Ec
## [1] 0.3484981#ANOVA- con la librería agricolae
library(agricolae)
model1 = with(Cafe, strip.plot(Bloques, VariedadN, TratN, AlturaP))##
## ANALYSIS STRIP PLOT: AlturaP
## Class level information
##
## VariedadN : V1 V2 V3
## TratN : 0.2 0.5 1.4 1.9
## Bloques : 1 2 3 4
##
## Number of observations: 50
##
## model Y: AlturaP ~ Bloques + VariedadN + Ea + TratN + Eb + TratN:VariedadN + Ec
##
## Analysis of Variance Table
##
## Response: AlturaP
## Df Sum Sq Mean Sq F value Pr(>F)
## Bloques 3 0.4028 0.13427 0.3853 0.7648
## VariedadN 2 2.2424 1.12121 3.7965 0.0860 .
## Ea 6 1.7720 0.29533 0.8474 0.5488
## TratN 3 0.9569 0.31896 1.3974 0.3057
## Eb 9 2.0543 0.22825 0.6550 0.7384
## TratN:VariedadN 6 0.7490 0.12483 0.3582 0.8966
## Ec 20 6.9700 0.34850
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 0.6 %, cv(b) = 0.5 %, cv(c) = 0.6 %, Mean = 98.06281Analisis Split plpot: Basado en el calculo anterior se tiene que el término de interacción es significativo en un nivel de probabilidad del 5%, de modo que cada dosis de Nitrogeno difiere entre sí y por tanto, causa diferencias en los la altura de nuestro experimento.
modelo <- lmer(AlturaP~Bloques+VariedadN*TratN+ (1|Bloques:VariedadN), data = Cafe)
anova(modelo)## Analysis of Variance Table
## npar Sum Sq Mean Sq F value
## Bloques 1 0.00884 0.00884 0.0287
## VariedadN 2 1.86980 0.93490 3.0329
## TratN 3 0.92456 0.30819 0.9998
## VariedadN:TratN 6 0.51445 0.08574 0.2782Añadimos Covariables:
En este caso, se escogió como covariable la Temperatura, medida como (C°). Pues esta esta fuertemente relacionada con la altura de la mayoria de plantas de cafe #Fuente:https://www.canna.es/influencia_temperatura_ambiental_en_las_plantas#:~:text=Junto%20con%20los%20niveles%20de,corto%20como%20a%20largo%20plazo.
temp <- sort.default(runif(50, 10, 20), decreasing = TRUE)
temp2 <- data.frame(Bloques,VariedadN,TratN,temp, AlturaP)
temp2## Bloques VariedadN TratN temp AlturaP
## 1 1 V1 0.2 19.99317 98.18849
## 2 1 V1 0.5 19.71871 98.15077
## 3 1 V1 1.4 19.50119 97.45099
## 4 1 V1 1.9 19.36078 97.43480
## 5 1 V2 0.2 19.30595 96.60173
## 6 1 V2 0.5 19.26770 98.36029
## 7 1 V2 1.4 18.97273 98.46956
## 8 1 V2 1.9 18.41939 97.88531
## 9 1 V3 0.2 18.36661 98.87957
## 10 1 V3 0.5 18.30271 98.05868
## 11 1 V3 1.4 18.21975 97.57344
## 12 1 V3 1.9 17.86937 98.45463
## 13 1 V1 0.2 17.85995 98.59819
## 14 1 V1 0.5 17.26200 97.81421
## 15 2 V1 1.4 17.18649 97.93837
## 16 2 V1 1.9 16.83804 98.90002
## 17 2 V2 0.2 16.82981 98.85200
## 18 2 V2 0.5 16.81611 96.48062
## 19 2 V2 1.4 16.43804 96.85551
## 20 2 V2 1.9 16.30077 98.02915
## 21 2 V3 0.2 16.07848 99.08718
## 22 2 V3 0.5 15.98999 98.54909
## 23 2 V3 1.4 15.87435 98.15911
## 24 2 V3 1.9 15.66811 97.96343
## 25 2 V1 0.2 15.48244 98.41713
## 26 2 V1 0.5 15.47265 98.09938
## 27 2 V1 1.4 15.37035 98.64892
## 28 2 V1 1.9 15.18330 98.46836
## 29 3 V2 0.2 15.17052 97.92628
## 30 3 V2 0.5 15.00602 98.05522
## 31 3 V2 1.4 14.65515 97.59375
## 32 3 V2 1.9 14.26204 97.62815
## 33 3 V3 0.2 14.25129 98.54767
## 34 3 V3 0.5 14.18582 99.21769
## 35 3 V3 1.4 14.02943 98.19406
## 36 3 V3 1.9 13.75112 98.14531
## 37 3 V1 0.2 13.68747 97.85720
## 38 3 V1 0.5 12.85383 98.03801
## 39 3 V1 1.4 12.66088 97.71985
## 40 3 V1 1.9 12.50232 98.22359
## 41 3 V2 0.2 12.48112 98.45425
## 42 3 V2 0.5 12.23560 97.74747
## 43 4 V2 1.4 11.82961 97.84950
## 44 4 V2 1.9 11.63166 97.63698
## 45 4 V3 0.2 10.92044 97.40996
## 46 4 V3 0.5 10.65376 98.12654
## 47 4 V3 1.4 10.54986 97.81464
## 48 4 V3 1.9 10.42996 98.01109
## 49 4 V1 0.2 10.28819 98.33002
## 50 4 V1 0.5 10.28754 98.24440length(temp)## [1] 50length(Bloques)## [1] 50pt=anova(lm(AlturaP~VariedadN+Bloques+temp));pt## Analysis of Variance Table
##
## Response: AlturaP
## Df Sum Sq Mean Sq F value Pr(>F)
## VariedadN 2 2.0531 1.02655 3.5454 0.03756 *
## Bloques 3 0.5922 0.19738 0.6817 0.56808
## temp 1 0.0518 0.05179 0.1789 0.67447
## Residuals 43 12.4503 0.28954
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Encontrando en este caso que las variedades son responsables de causar diferencias las temperaturas.
plot(temp2$temp,temp2$AlturaP)
Grafico Disperso sin correlación, por ende no existe la necesidad de realizar analisis de covarianza
Para corroborar lo anteriormente dicho, ejecutamos un cor test
cor.test(temp2$temp,temp2$AlturaP)##
## Pearson's product-moment correlation
##
## data: temp2$temp and temp2$AlturaP
## t = -0.16106, df = 48, p-value = 0.8727
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.2996503 0.2567677
## sample estimates:
## cor
## -0.02324103Arrojandonos un valor de -0.02 mostrando con correlación casí existente
library(agricolae)
# First factor variable (varieties)
LSD_A = LSD.test(y = AlturaP,
trt = VariedadN,
DFerror = modeloInicial$gl.a,
MSerror = modeloInicial$Ea,
alpha = 0.05,
p.adj = "bonferroni",
group = TRUE,
console = TRUE)##
## Study: AlturaP ~ VariedadN
##
## LSD t Test for AlturaP
## P value adjustment method: bonferroni
##
## Mean Square Error: 0.2953316
##
## VariedadN, means and individual ( 95 %) CI
##
## AlturaP std r LCL UCL Min Max
## V1 98.14015 0.3966811 18 97.82672 98.45358 97.43480 98.90002
## V2 97.77661 0.6634064 16 97.44417 98.10905 96.48062 98.85200
## V3 98.26201 0.5044871 16 97.92957 98.59445 97.40996 99.21769
##
## Alpha: 0.05 ; DF Error: 6
## Critical Value of t: 3.287455
##
## Groups according to probability of means differences and alpha level( 0.05 )
##
## Treatments with the same letter are not significantly different.
##
## AlturaP groups
## V3 98.26201 a
## V1 98.14015 a
## V2 97.77661 a#Intento adicional de usar algunos test encontrados en internet LSDEn el primer caso mostrando que los tratamientos con la misma letra no son significativamente diferentes.
Usamos LSD.test() como prueba adicional para calcular as diferencias significativas entre los valores medios. Esta prueba crea múltiples comparaciones de tratamientos por medio del LSD y una agrupación de tratamientos. En el cual el argumento trt representa el vector de los tratamientos aplicados a cada unidad experimental, en nuestro caso las variedades aplciadas de Nitrogeno. El valor de DFerror puede especificarse utilizando model\(gl.a, model\)gl.b y model$gl.c para las variedades. anteriormente habia sido extraido estos valores
Octavo Punto
Realice un resumen con la nota que aparece en las siguientes direcciones sobre: El uso de los diseños en parcelas divididas: http://207.67.83.164/quality-progress/2007/10/laboratory/when-should-you-consider-a-split-plot-design.html
When Should You Consider A Split-Plot Design?
Un diseño de parcelas divididas es una opción muy viable cuando tenemos una condición en la cual algunos factores tienen niveles muy difíciles de cambiar con respecto a otros, por lo que el Diseño de parcelas divididas, intenta que los factores difíciles de cambiar se reajusten con menos frecuencia que los factores fáciles de cambiar. Para entender este diseño debemos conocer que la mayoría de los términos empleados se originaron en entornos agrícolas. Parcelas de tierra: Los factores difíciles de cambiar se denominan factores de parcela completa, mientras que los factores más fáciles de cambiar se denominan factores de subparcela.
Un ejemplo adecuado para exponer lo dicho puede ser un experimento con varias parcelas adyacentes, en las cuales todos reciben el mismo fertilizante (Factor) y luego el tipo de semilla particular (factor de subparcela) se asigna al azar y se aplica a las parcelas dentro de ese grupo mayor. Donde el número total de los experimentos (N) es igual al número de subparcelas, además en el cual el número de parcelas enteras (WP) corresponde al número de veces que se reinician los factores de las parcelas enteras. Existen una cantidad de ventajas asociadas a la implementación de este diseño, ya que algunos muy puntuales, en caso de no efectuar el diseño de parcelas divididas, se podría caer en el error de ejecutar el problema bajo uno completamente al azar, cuyos resultados podrían declarar erróneamente un factor significativo cuando en realidad no lo es.
En el desarrollo del diseño de experimentos la aleatorización puede estar sujeta a diversas modificaciones, todo esto dependiendo del investigador, en el que este mismo puede elegir una forma de aleatorización más estructurada con el fin de reajustar todos sus valores, generando muchas veces puede tener mejor precisión. Sin realizar dicho reajuste todo el factor de la parcela puede estar correlacionadas entre sí, ya que se esperaría que fueran más similares que las observaciones que implicaban reajustar todos los factores.
Una buena justificación para la ejecución de este experimentó es tener en cuenta el costo de ejecución, pues si los factores difíciles de cambiar son caros o lleva mucho tiempo cambiarlos, entonces puede ser mejor seleccionar un diseño que contenga menos parcelas enteras. Además, en muchas ocasiones los factores difíciles de cambiar son extremadamente difíciles o costosos de cambiar en relación con la reubicación de los factores de la subparcela. En definitiva, podría representar un ahorro potencialmente grande.
Por último, para la elección del diseño correcto existen múltiples aspectos, como el equilibrio entre observaciones por parcela entera y el número de niveles de cada factor, además una buena estimación de los términos de error de la parcela entera y de las subparcelas. Aunque generalmente lo óptimo puede variar según el experimento.
Short communication: On recognizing the proper experimental unit in animal studies in the dairy sciences((https://www.sciencedirect.com/science/article/pii/S002203021630621X)
Un diseño adecuado de experimentos, combinado con la aplicación correcta de métodos es una herramienta útil para producir resultados científicos significativos que sean a la vez reproducibles. Partiendo del supuesto en el cual debemos considerar una “unidad experimental” como la entidad más pequeña que se asigna de forma independiente de todas las demás unidades a un tratamiento concreto, además de ser en la cual se aplica un tratamiento de forma independiente. Con este concepto es pertinente aclarar el concepto de “unidad de observación”, también conocida como unidad de muestreo, que se define como la entidad física sobre la que se mide un resultado de interés en un experimento. Esto como pie de inicio para poder entender ampliamente el experimento que será analizado basado en la ciencia lechera.
Partiendo de un ejemplo con Dietas o antibióticos aplicadas a un corral con x individuos de vacas o en otros casos usando cada vaca como una unidad, en la cual intentan explicar múltiples conceptos, por medio de la ejemplificación de algunos diseños aplicados con los corrales o vacas. Logrando explicar una estructura de diseño anidado, la cual se ejemplifica como el corral está anidado dentro de un tratamiento y la vaca individual está anidada dentro de un corral, creando así una estructura jerárquica en los datos, siendo esta última definida como la configuración de datos en la que las observaciones no son mutuamente independientes sino al revés. tienen alguna correlación dada por el diseño experimental, mostrando que en los casos en los que existe correlación no es posible separar los efectos independientemente. Evidenciando que siempre que las unidades de observación se encuentren anidadas dentro de una unidad experimental, se les denominará submuestras, para indicar que estos datos están correlacionados y por lo tanto no construye una observación entera. Estas estructuras jerárquicas muchas veces se les conoce como andantes o bloqueadoras en el diseño experimental. Para poder identificar la estructura jerárquica de los datos; es decir, cuando se produce una replicación independiente y cuando no, se puede lograr por medio del ejercicio llamado ¿Qué haría Fisher?" (WWFD) es una estrategia general que otorga un proceso de construcción de modelos mixtos lineales generalizados, logrando explorar las diferencias entre efectos fijos y aleatorios. La importancia de los diseños cuadrados latinos en las ciencias lecheras en las cuales, cuando se repite más de un cuadrado se le conoce como latinos replicados o rectángulos latinos, en las cuales el periodo como factor es de suma importancia en donde se le debería tratar a como un elemento del diseño experimental o de la estructura del tratamiento, puesto que tiene consecuencias para la inferencia posterior porque determina cómo se define el error experimental para el tratamiento de interés.
En rasgos generales se menciona que diversas entidades físicas pueden servir como unidad experimental o de observación dentro de un mismo experimento. Mostrando que en la ciencia lechera puede ser capaz de aplicar lo aprendido en una unidad experimental, generando que la aplicación correcta de esto nos ayude a reducir posibles errores. Además de mostrar la importancia de identificar correctamente la unidad experimental en estudios determinados, en la cual para saber con plena exactitud debemos establecer las siguientes preguntas: ¿cuál es el efecto de un tratamiento en un resultado de interés concreto? es decir, el diseño experimental y la estructura de datos correspondiente. Lo más importante que cabe destacar es que la especificación de la unidad experimental en un diseño determinado no es una cuestión de opinión, sino que está determinada por la forma en que se estableció el experimento y de la forma en que se recogieron los datos y el alcance previsto de la inferencia.
Unidad experimental y replicación https://online.stat.psu.edu/stat502/lesson/6/6.1-0
Una unidad experimental es donde se aplica el tratamiento de un experimento, por ejemplo, se quiere evaluar las lesiones que el agua de un arroyo contaminado causa a unos peces, para esto se usan dos acuarios, a los cuales se les asigna aleatoriamente un tratamiento, es decir agua contaminada o agua sin contaminar. A cada acuario se le agregan 50 peces y pasados 30 días se capturan 10 para evaluar el número de lesiones, este es un diseño de factor único con dos niveles de tratamiento por lo que sería adecuado aplicarle un análisis de varianza de una vía. En este caso la unidad experimental sería el acuario que contiene 50 peces, no serían los peces pues a estos no se les aplica individualmente el tratamiento.
Por otro lado, una réplica se entiende como cada vez que se aplica todos los niveles de tratamiento, es decir que en el ejemplo una réplica sería cada vez que se toman dos acuarios, se les agrega aleatoriamente agua contaminada o agua sin contaminar y se introducen a cada uno 50 peces que permanecen durante 30 días. Por ejemplo, al aplicarse solo una vez el conjunto de tratamientos se habla de un estudio no replicado.
Las unidades de muestreo en el ejemplo fueron los peces que fueron capturados y contados para evaluar sus lesiones. Es importante no tratar las unidades de muestreo como si fueran unidades experimentales pues esto causa pseudo - replicación e infla el Error d.f., de modo que se reduce el MSE y se produce un estadístico F incorrecto.
Fundamentos del diseño experimental: pautas para diseñar experimentos exitosos https://acsess.onlinelibrary.wiley.com/doi/full/10.2134/agronj2013.0114
Para los campos de la biología, los experimentos comparativos o manipulativos son los más usados, este tipo de experimentos permiten comparar o contrastar dos o más prácticas o sistemas, durante estos experimentos se sigue la teoría general de la investigación científica. Se comienza con una pregunta y/o hipótesis, y de esta se realizan modelos basados en el tema a investigar, este a su vez se traduce en un modelo estadístico desarrollado juntamente con un diseño estadístico. El diseño estadístico reúne el diseño del tratamiento y el diseño experimental, y proporciona las reglas y procedimientos que permiten realizar el experimento.
Una vez se obtengan los datos, se pasa a realizar el análisis estadístico previamente planeado, de modo que se obtienen interpretaciones de los resultados, creando inferencias y conclusiones sobre la pregunta o hipótesis inicial, adicionalmente debe haber un ciclo de retroalimentación cuando se completa el experimento, compuesto por respuestas científicas a preguntas biológicas e información estadística que servirán para el planteamiento de futuros experimentos mejores y más eficientes. En el caso de los experimentos que tienen que ver con el campo de la biología, un resultado negativo es un fracaso por el costo de la inversión tanto capital como emocional. Un resultado negativo puede significar que el experimento está mal diseñado por lo que no tienen poder suficiente para detectar diferencias entre dos medias, o que los tratamientos están mal diseñados y no refleja adecuadamente la pregunta o hipótesis inicial, o que no se supervisó adecuadamente los protocolos de tratamiento y la recopilación de datos, o que verdaderamente no existen diferencias entre tratamientos. Sin embargo, si el investigador formula una pregunta o hipótesis válida y piensa bien sobre el diseño de los tratamientos los diseños deberían arrojar diferencias entre tratamientos. Teniendo en cuenta todo lo que puede indicar un resultado negativo los editores experimentados generalmente rechazan artículos con resultados negativos en especial cuando estos se basan en una falta de significancia estadística. Con el fin de evitar fallas en la experimentación comparativa es necesario responder a la pregunta “¿Qué hace que un experimento sea mejor o más eficiente?”, revisando cuatro pilares del diseño experimental (replicación, aleatorización, bloqueo y unidades experimentales).
Replicación
- Forma y escala de replicación La replicación tiene cuatro funciones dentro de la experimentación comparativa, primero la replicación de un mecanismo para estimar el error experimental, esencial para proporcionar pruebas de hipótesis válidas e intervalos de confianza, lo que permite estimar la variación de fondo en las estimaciones de los efectos del tratamiento. En segundo y tercer lugar, la replicación aumenta la precisión y el alcance de la inferencia del experimento, es decir que la replicación conlleva un rango más amplio de observaciones, lo que aumenta el rango de condiciones encontradas durante el curso del experimento. Por último, la replicación también afecta el control del error, pues permite al investigador controlar la magnitud del error experimental y regular el nivel de precisión o potencia para el experimento.
Cuando se realiza replicación es fundamental que las observaciones repetidas ocurran a una escala espacial y temporal que coincidan con la aplicación de los tratamientos, esto por la variabilidad inherente que existe en los sistemas biológicos y para evitar confusiones de tratamientos con otros factores que puedan variar entre las unidades experimentales, los diseños experimentales deben dar una interpretación inequívoca de los resultados. Para esta situación en la que hay confusión completa de tratamientos con unidades experimentales se anidó el término pseudorreplicación. En el caso de que un factor de tratamiento no se replique a una escala adecuada, se creará una confusión por la falta de replicación del tratamiento que se extenderá a todo el experimento, provocando que un efecto simple se atribuya incorrectamente, por ejemplo, a variable ambiental no repetidas.
Para que una replicación sea exitosa primero se debe definir explícitamente la unidad experimental, después de definir la unidad experimental se diseñan los niveles adicionales de replicación experimental. Durante esto se debe tomar en cuenta en primer lugar las escalas que son necesarias para lograr las metas experimentales, esto depende de la inferencia deseada y de cómo se recopilarán los datos. En segundo lugar se toma en cuenta los recursos que se usarán para la replicación de todos los niveles, esto depende más que todo de las condiciones del experimento ya que muchas variables no pueden medirse en unidades experimentales completas.
- Número de réplicas
La replicación experimental ocurre en cuatro niveles básicos, en primer lugar la unidad experimental, segundo lugar la replicación de todo el experimento, tercero, muestreo en uno o más niveles dentro de las unidades experimentales y cuarto, medidas repetidas. La teoría clásica dicta que la asignación de recursos se basa en estimaciones precisas de los componentes de varianza de experimentos anteriores. El número de veces que se replica es fundamental ya que incide directamente en la precisión y la capacidad de detectar diferencias entre tratamientos, así mismo si se utiliza el número de replicaciones adecuadas se podrá percibir mejor las diferencias entre tratamientos y así obtener un diseño de experimentos que tenga un amplio potencial estadístico. Para desarrollar un análisis de la potencia que tendrá un experimento, Gubur et al. (2012) brindó un método de distribución de probabilidad con los siguientes pasos: Obtener una estimación de la variabilidad del error experimental a partir de experimentos previos realizados; identificar las distribución de la variable; determinar el valor de P que se utilizará para establecer los límites de detección y la diferencia mínima de tratamientos a detectar; Elegir una estructura de diseño experimental que incluya los arreglos de bloqueo y los niveles de replicación, y por último, crear un conjuntos de datos ejemplares que coincidan con la estructura del diseño.
- Experimentos no repetidos
Existen ocasiones en que existe la necesidad de dedicar todos los recursos a múltiples tratamientos pero no a replicas u observación independiente de los tratamientos. A menudo los investigadores que desarrollan este tipo de experimentos tienen tres opciones, primero, realizar los experimentos en múltiples fincas, tomando las fincas como bloques, segundo, usar diseños de parcelas de control en los que los tratamientos están intercalados y tercero, utilizar una combinación de los dos enfoques. En este tema también entra en juego el costo de las unidades experimentales que puede favorecer el uso de tratamientos adicionales en vez de réplicas.
Aleatorización
La correcta realización de experimentos tiene que ver con el principio de aleatorización en dos niveles, primero para una cuidadosa definición de los materiales e instalaciones experimentales se exige que estos estén muestreados para garantizar que estén correctamente representados. El segundo nivel involucra la asignación de tratamientos a unidades experimentales. La aleatorización permitirá una estimación no sesgada de las medias del tratamiento y los errores experimentales y tener una precaución contra una perturbación que pueda surgir, especialmente en el campo de la biología la aleatorización garantiza que las fuentes de variación desconocida o inesperada produzcan sesgo, confusión o eliminación de pruebas de hipótesis válidas.
Bloqueo
El bloqueo se usa con el objetivo de aumentar la precisión, crea grupos de unidades experimentales que sean más homogéneas, el bloqueo también permite utilizar distintos tamaños de unidades experimentales, que facilitaran la aplicación de un experimento que requiera un área experimental más grande para un determinado tratamiento. Además el bloqueo funciona también como una precaución ante posibles perturbaciones que puedan surgir, por ejemplo el diseño RCBD se utiliza para recuperarse con exito de las perturbaciones posteriores al establecimiento de los experimentos, este diseño es el de bloqueo simple y tienen la restricción de que el tamaño de los bloques debe ser igual al número de tratamientos. El diseño en parcelas divididas existe como una restricción específica de la aleatorización que se coloca en algunos experimentos factoriales.
La aleatorización para un diseño de parcelas divididas es flexible y versátil, aplicable a diseños de tratamiento factorial con dos o más factores. Usualmente se usa este diseño cuando el factor de las parcela completa requiere unidades experimentales más grandes que los factores de la subparcela. Añadiendo a los bloques existe una correlación directa y positiva entre la complejidad y la eficiencia del diseño experimental, los diseños más simples como CRD y RCBD son los menos eficientes.
Tamaño de las unidades experimentales.
El concepto de tamaño óptimo de las unidades experimentales se basa en la ley de Smith que se deriva de la observación general de una relación asintótica negativa entre varianzas y el tamaño de la parcela, la ley de Smith se aplica a una amplia gama de especies y condiciones ambientales. La naturaleza asintótica de la relación arroja que para tamaños de parcela inicialmente pequeños, se espera que cualquier cambio tenga un gran efecto en la varianza media de la unidad, por el contrario para parcelas grandes los cambios grandes pueden tener un impacto poco o nulo en la varianza media de la unidad. Es de agregar que la definición de que si una parcela es “grande” o “pequeña” es relativa y puede variar entre especies, tipos de suelo y muchos otros factores. A la hora de decidir si se debe aumentar el tamaño de las parcelas se debe tener la ubicación en la curva en la que esta la unidad experimental.
Noveno Punto
Seleccionar un artículo científico de una revista de agronomía donde se haya utilizado un diseño en parcelas divididas. Hacer las críticas constructivas sobre:
La mención de la estructura factorial. Se usó un diseño experimental completamente al azar con cuatro repeticiones, y un arreglo de tratamientos en parcelas divididas. En las parcelas divididas, las parcelas grandes fueron las cinco fechas de siembra y las subparcelas estuvieron constituidas por las 28 variedades.Se establecieron 28 unidades experimentales que se conformaron de cuatro surcos con 5 m de largo × 0.3 m de ancho, y de cada una se realizaron cuatro repeticiones.
La razón de colocar cada factor en la parcela principal o en la subparcela. Las parcelas principales fueron las cinco fechas de siembra y las subparcelas fueron los 28 genotipos de avena, la elección de elegir los genotipos como parcela subprincipal es debido a que estos pueden asociarse fácilmente a cada uno de los factores a unidades experimentales de mayor tamaño (Fechas de siembra ), generando un mayor ahorro en la producción del experimento.
La revisión de supuestos para el análisis de varianza. De acuerdo al Anova, en las cinco fechas de siembra, hubo diferencias altamente significativas, esto con un valor de confianza( p= 0.05), por lo cual la fecha de siembra influye en la variable respuesta.
La tabla del análisis de varianza. No es presentada en el artículo, pero hubiese sido adecuado la creación de esta, ya que hubiese facilitado el análisis e interpretación por parte de los lectores.
El uso de muchos análisis de varianzas en lugar de uno solo multivariante. Solo se ejecutó un ANOVA con el fin de determinar dos análisis, uno entre variedades-fechas de siembra, y otro simple con fechas de siembra.
El método de comparaciones de medias después del Anova. (DMS) Comparaciones de medias empleando la prueba de diferencia mínima significativa, que identifican que pares de medias son diferentes.
La interacción de factores. Interacciones entre cada una de las variedades y en todas las fechas de siembra( Factores), en estos se observaron que, en la mayoría de las variedades de avena, la fecha 1 y 2 podrían considerarse como las mejores.
La presencia de bloques. No presenta bloques , aunque hubiese sido adecuado la presencia de algunos para poder ver la variación en la variable de respuesta que no es causada por los factores. Un bloque que podría ser alguna variable climática.
El balance o desbalance. Balanceado, ya que no presenta la ausencia de algún dato.
La definición clara de la unidad experimental. Las a 28 genotipos-variedades de planta de avena, sembradas en cada fecha particular con distancia de 7 días.
Software utilizado y librería específica (en caso de ser R). SAS (Versión 9.1), Curve Expert 2 para hacer el modelo de progreso de la epidemia de tipo exponencial.
Otros aspectos que considere de interés:
La variabilidad entre fechas de siembra también puede deberse a la influencia de las siembras anteriores causada por un protocolo deficiente o un mal tratamiento.
Podrá tomar en cuenta factores como el clima ya que saber la mejor fecha de siembra solo para 2013-2014 no nos da una respuesta precisa a la pregunta de cuándo es mejor sembrar ya que el clima no es constante todos los años.
No se definió con claridad la unidad experimental, lo que puede conllevar a algunas confusiones entre los lectores.
La explicación del modelo (ABCPE ) se omite solo mostrando su utilidad.