parcial
Felipe Palacio, Nashly salazar
9/11/2020
Ejercicios basados en la biblioteca agricolae
Ejercicio 11 Diseño en franjas divididas
Lo primero que hago es descargar datos de un archivo publico sobre el maíz este sera un ejercicio bifactorial en donde se toman 4 diferentes bloques en donde en cada bloque se realiza 3 aplicaciones de nitrogeno (0, 50 y 100) y dos metodos de riego uno alto y uno bajo en donde se busca hallar si esto interviene en el rendimiento del maiz
library(readr)
Datos <- read_delim("C:/Users/felip/OneDrive/Escritorio/Diseño y experimentos/Parcial/Maiz.txt",
"\t", escape_double = FALSE, trim_ws = TRUE)##
## -- Column specification --------------------------------------------------------
## cols(
## Bloque = col_double(),
## Nitrog = col_double(),
## Riego = col_character(),
## Rend = col_double()
## )View(Datos)Ahora visualizo los primeros registros de tabla
print(head(Datos))## # A tibble: 6 x 4
## Bloque Nitrog Riego Rend
## <dbl> <dbl> <chr> <dbl>
## 1 1 0 Bajo 55
## 2 1 100 Bajo 69
## 3 1 50 Bajo 62
## 4 1 0 Alto 71
## 5 1 100 Alto 78
## 6 1 50 Alto 77Ahora es importante seleccionar los datos para realizar el modelo y el analisis de la varianza y establecemos vectores para cada variable en donde los vectores seran; nitrogeno, riego, bloques, rendimiento este sera la variable respuesta
Luego vamos a establecer un vector numerico para hacer mas facil el manejo en r
nitrogeno <- factor(Datos$Nitrog)
riego <- factor(Datos$Riego)
bloque <-factor(Datos$Bloque)
rend <-as.vector(Datos$Rend)
Re1 <-as.numeric(rend)library(agricolae)vamos a usar un Análisis de varianza usando la función Strip-Plot, este tipo de analisis se usa cuando se tienen 2 factores en el diseño y la presición deseada para medir el efecto entre las interacciones de los dos factores es mayor que el de midiendo el efecto principal de cualquiera de los dos factores esto se logra gracias a que con el con el uso de tres tamaños de parcela, existen varios graficos; gráfico de bandas verticales para el primer factor también llamado factor vertical, gráfico de franjas horizontales para el segundofactor también llamado factor horizontal, gráfico de intersección para la interacciónentre los dos factores la vertical parcela de franjas y franjas horizontales son siempre perpendiculares entre sí sin embargo, no hay relación entre sus tamaños a diferencia del caso de la trama principaly subparcela del diseño de parcela dividida en este caso la gráfica de interacción es de Por supuesto, el más pequeño desde la parcela de viaje.diseñar los grados de precisión asociado con los principales efectos de ambos Se sacrifican factores paramejorar la precisión de la interacciónefecto en el diseño de parcela divididalos niveles del segundo factor o factor B sonindependientemente aleatorizados dentro delsubparcelas de cada nivel del primer factor o factor a en el diseño del diagrama de bandaslos niveles de un factor se asignan a tiras de parcelas en una dirección y el nivel del segundo factor a la tira es perpendicular a la primera tirase realiza una aleatorización separada paracada bloque para cada factor a y B este diseño facilita el físico operaciones y aumenta la precisión para la estimación del efecto de interacción
fran.anova<-strip.plot(bloque, nitrogeno, riego, Re1)##
## ANALYSIS STRIP PLOT: Re1
## Class level information
##
## nitrogeno : 0 100 50
## riego : Bajo Alto
## bloque : 1 2 3 4
##
## Number of observations: 24
##
## model Y: Re1 ~ bloque + nitrogeno + Ea + riego + Eb + riego:nitrogeno + Ec
##
## Analysis of Variance Table
##
## Response: Re1
## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 3 123.46 41.15 28.7670 0.0005860 ***
## nitrogeno 2 339.08 169.54 60.1330 0.0001073 ***
## Ea 6 16.92 2.82 1.9709 0.2147416
## riego 1 570.38 570.38 52.1817 0.0054703 **
## Eb 3 32.79 10.93 7.6408 0.0179415 *
## riego:nitrogeno 2 94.75 47.38 33.1165 0.0005731 ***
## Ec 6 8.58 1.43
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 2.3 %, cv(b) = 4.6 %, cv(c) = 1.7 %, Mean = 72.45833fran.anova## $data
## bloque nitrogeno riego Re1
## 1 1 0 Bajo 55
## 2 1 100 Bajo 69
## 3 1 50 Bajo 62
## 4 1 0 Alto 71
## 5 1 100 Alto 78
## 6 1 50 Alto 77
## 7 2 100 Alto 81
## 8 2 0 Alto 77
## 9 2 50 Alto 79
## 10 2 100 Bajo 77
## 11 2 0 Bajo 63
## 12 2 50 Bajo 66
## 13 3 50 Bajo 70
## 14 3 100 Bajo 79
## 15 3 0 Bajo 63
## 16 3 50 Alto 78
## 17 3 100 Alto 80
## 18 3 0 Alto 77
## 19 4 100 Bajo 76
## 20 4 50 Bajo 66
## 21 4 0 Bajo 65
## 22 4 100 Alto 79
## 23 4 50 Alto 76
## 24 4 0 Alto 75
##
## $ANOVA
## Analysis of Variance Table
##
## Response: Re1
## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 3 123.46 41.15 28.7670 0.0005860 ***
## nitrogeno 2 339.08 169.54 60.1330 0.0001073 ***
## Ea 6 16.92 2.82 1.9709 0.2147416
## riego 1 570.38 570.38 52.1817 0.0054703 **
## Eb 3 32.79 10.93 7.6408 0.0179415 *
## riego:nitrogeno 2 94.75 47.38 33.1165 0.0005731 ***
## Ec 6 8.58 1.43
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $gl.a
## [1] 6
##
## $gl.b
## [1] 3
##
## $gl.c
## [1] 6
##
## $Ea
## [1] 2.819444
##
## $Eb
## [1] 10.93056
##
## $Ec
## [1] 1.430556Se observa que si hay una interacción entre riego y nitrogeno asi que ahora aplicamos turkey
Prueba de Tukey
Para Dosis de Nitrógeno
FactorA<-HSD.test(Re1,nitrogeno,fran.anova$gl.a,fran.anova$Ea);FactorA## $statistics
## MSerror Df Mean CV MSD
## 2.819444 6 72.45833 2.31736 2.576001
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey nitrogeno 3 4.339195 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## 0 68.250 7.995534 8 55 77 63.00 68.0 75.50
## 100 77.375 3.739270 8 69 81 76.75 78.5 79.25
## 50 71.750 6.562883 8 62 79 66.00 73.0 77.25
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## 100 77.375 a
## 50 71.750 b
## 0 68.250 c
##
## attr(,"class")
## [1] "group"Para Niveles de Riego
FactorB<-HSD.test(Re1,riego,fran.anova$gl.b,fran.anova$Eb);FactorB## $statistics
## MSerror Df Mean CV MSD
## 10.93056 3 72.45833 4.562814 4.295429
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey riego 2 4.500659 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## Alto 77.33333 2.605356 12 71 81 76.75 77.5 79.0
## Bajo 67.58333 7.012435 12 55 79 63.00 66.0 71.5
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## Alto 77.33333 a
## Bajo 67.58333 b
##
## attr(,"class")
## [1] "group"Para interacción de Dosis de Nitrógeno y Niveles de Riego
FactorAB<-HSD.test(Re1,riego:nitrogeno,fran.anova$gl.c,fran.anova$Ec);FactorAB## $statistics
## MSerror Df Mean CV MSD
## 1.430556 6 72.45833 1.650684 3.365919
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey riego:nitrogeno 6 5.628353 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## Alto:0 75.00 2.828427 4 71 77 74.00 76.0 77.00
## Alto:100 79.50 1.290994 4 78 81 78.75 79.5 80.25
## Alto:50 77.50 1.290994 4 76 79 76.75 77.5 78.25
## Bajo:0 61.50 4.434712 4 55 65 61.00 63.0 63.50
## Bajo:100 75.25 4.349329 4 69 79 74.25 76.5 77.50
## Bajo:50 66.00 3.265986 4 62 70 65.00 66.0 67.00
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## Alto:100 79.50 a
## Alto:50 77.50 ab
## Bajo:100 75.25 b
## Alto:0 75.00 b
## Bajo:50 66.00 c
## Bajo:0 61.50 d
##
## attr(,"class")
## [1] "group"Ahora realizamos pruebas de comparación múltiple de medias y Gráficas de barras
Prueba HSD de Tukey para Dosis de Nitrógeno
out1<-with(fran.anova,HSD.test(Re1,nitrogeno,fran.anova$gl.a,fran.anova$Ea));out1## $statistics
## MSerror Df Mean CV MSD
## 2.819444 6 72.45833 2.31736 2.576001
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey nitrogeno 3 4.339195 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## 0 68.250 7.995534 8 55 77 63.00 68.0 75.50
## 100 77.375 3.739270 8 69 81 76.75 78.5 79.25
## 50 71.750 6.562883 8 62 79 66.00 73.0 77.25
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## 100 77.375 a
## 50 71.750 b
## 0 68.250 c
##
## attr(,"class")
## [1] "group"Se vizualiza como gráfica de barras para Dosis de Nitrógeno
bar.err(out1$means, variation ="SE", ylim=c(0,90), xlab="Dosis de Nitrógeno", ylab="Rendimiento")
Prueba HSD de Tukey para Niveles de Riego
out2<-with(fran.anova,HSD.test(Re1,riego,fran.anova$gl.b,fran.anova$Eb));out2## $statistics
## MSerror Df Mean CV MSD
## 10.93056 3 72.45833 4.562814 4.295429
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey riego 2 4.500659 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## Alto 77.33333 2.605356 12 71 81 76.75 77.5 79.0
## Bajo 67.58333 7.012435 12 55 79 63.00 66.0 71.5
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## Alto 77.33333 a
## Bajo 67.58333 b
##
## attr(,"class")
## [1] "group"Gráfica de barras para Niveles de Riego
bar.err(out2$means, variation ="SE", ylim=c(0,90), xlab="Riego", ylab="Rendimiento")
Prueba HSD de Tukey para la interacción
out3<-with(fran.anova,HSD.test(Re1,riego:nitrogeno,fran.anova$gl.c,fran.anova$Ec));out3## $statistics
## MSerror Df Mean CV MSD
## 1.430556 6 72.45833 1.650684 3.365919
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey riego:nitrogeno 6 5.628353 0.05
##
## $means
## Re1 std r Min Max Q25 Q50 Q75
## Alto:0 75.00 2.828427 4 71 77 74.00 76.0 77.00
## Alto:100 79.50 1.290994 4 78 81 78.75 79.5 80.25
## Alto:50 77.50 1.290994 4 76 79 76.75 77.5 78.25
## Bajo:0 61.50 4.434712 4 55 65 61.00 63.0 63.50
## Bajo:100 75.25 4.349329 4 69 79 74.25 76.5 77.50
## Bajo:50 66.00 3.265986 4 62 70 65.00 66.0 67.00
##
## $comparison
## NULL
##
## $groups
## Re1 groups
## Alto:100 79.50 a
## Alto:50 77.50 ab
## Bajo:100 75.25 b
## Alto:0 75.00 b
## Bajo:50 66.00 c
## Bajo:0 61.50 d
##
## attr(,"class")
## [1] "group"Gráfica de barras para la interacción
bar.err(out3$means, variation ="SE", ylim=c(0,90), xlab="Riego por Nitrógeno", ylab="Rendimiento")
Gráfico de interacción
interaction.plot(nitrogeno,riego,Re1, fixed=T, xlab="Dosis de Nitrógeno", ylab="Rendimiento", pch = c(20,21),type = "b") 
Esta ultima grafica nos ayuda a dilusidar una parte importante de el diseño y es que no hay interacion entre el riego y el nitrogeno en funcion del rendimiento
De esto se puede concluir que la aplicacion de 100 g de nitrogeno con un sistema de riego bajo ve un gran amento del rendimiento con respecto al mismo modelo en un sistema de alto riego, es mas si se analiza se puede ver como un sistema de alto riego apenas si presenta mejora de rendimiento con respecto a el aumento de la dosis de nitrogeno, esto puede deberse a la forma en la que se asimila el nitrogeno, se puede intuir que en el suelo el nitrogeno esta luego de su aplicacion y tras pasar por un proceso de inorganico a organico en forma de nitrato, el cual es sumamente lavable y se pierde por procesos de lixiviación lo cual explica el porque no se presenta un aumento considerable en el riego alto, ya que puede que se haya lavado y no haya sido absorvido por la planta, por ello la recomendacion preliminar (si no se tiene en cuenta el valor del nitrogeno) es mejor un alto riego para obtener mejores rendimientos en maíz
Por ultimo podemos hacer el analisis de los efectos simples
Efectos simples
para ello primero se crea un vector de datos A
A<-fran.anova$dataSe crea una lista a con los niveles de Dosis de Nitrógeno
a<-nlevels(A$Ni)Se crea una lista b con los Niveles de riego
b<-nlevels(A$Riego)Se crea una lista r con los niveles de Bloques
r<-nlevels(A$BLOF)Se crean vectores de datos para los errores experimentales; a, b y c
Ea<-fran.anova$Ea; Eb<-fran.anova$Eb; Ec<-fran.anova$EcSe crean vectores de datos para los grados de libertad asociados a los errores a, b y c respectivamente
gla<-fran.anova$gl.a; glb<-fran.anova$gl.b; glc<-fran.anova$gl.cMedias de las interacciones
B <-tapply.stat(A[,4],A[,2:3],mean);B## nitrogeno riego A[, 4]
## 1 0 Alto 75.00
## 2 0 Bajo 61.50
## 3 100 Alto 79.50
## 4 100 Bajo 75.25
## 5 50 Alto 77.50
## 6 50 Bajo 66.00Errores estándar de las medias de las interacciones
std<-tapply.stat(A[,4],A[,2:3],function(x) sd(x)/sqrt(length(x)));std## nitrogeno riego A[, 4]
## 1 0 Alto 1.4142136
## 2 0 Bajo 2.2173558
## 3 100 Alto 0.6454972
## 4 100 Bajo 2.1746647
## 5 50 Alto 0.6454972
## 6 50 Bajo 1.6329932Medias de las interacciones y sus errores estándares
B<-data.frame(B[,1:2],Re1=B[,3],std=std[,3]);B## nitrogeno riego Re1 std
## 1 0 Alto 75.00 1.4142136
## 2 0 Bajo 61.50 2.2173558
## 3 100 Alto 79.50 0.6454972
## 4 100 Bajo 75.25 2.1746647
## 5 50 Alto 77.50 0.6454972
## 6 50 Bajo 66.00 1.6329932Ejercicios basados en la biblioteca agricolae
Ejercicio 2 Diseño Carolino
Diseñados por Comstock y Robinson en 1948 en Carolina del Norte.Relacionado con la técnica de apareamiento entre progenitores (DISEÑOS GENÉTICOS) y puede ser de 3 formas el diseño:
Diseño Carolina I
También conocido domo diseño jerárquico o anidado. Acá cada macho se aparea con un grupo de hembras, las cuales estas tienen la restricción de que cada una solo participa en un solo cruzamiento.
Diseño Carolina II
También conocido como diseño cruzado o factorial. Un ejemplo es que en términos genéticos el apareamiento se da entre el cruce de un grupo de progenitores machos con un grupo de hembras y se dan en todas las combinaciones posibles.
Diseño Carolina III
Se crearon con el fin de estimar el grado de dominancia de los genes que controlan los rasgos de estudio. Este diseño tiene la finalidad de estimar la varianza aditiva y la de dominancia.
\(NOTA\)
Para este caso donde usaremos carolina no me centrare en un ejemplo de interes agronomico ya que en terminos geneticos las plantas son mucho mas complejas que los animales y no se si se puede aplicar de igual modo.
en este montaje quieren mirar la dominancia de un gen en especifico y la varianza que se puede generar entre hembras de su misma especie, diferente raza
Para este primer montaje se toma un macho y se cruza con 4 hembras y se busca la varianza del gen de interes con respecto a sus crias
library(agricolae)
data(DC)
carolina1 <- DC$carolina1
# str(carolina1)
output<-carolina(model=1,carolina1)## Response(y): yield
##
## Analysis of Variance Table
##
## Response: y
## Df Sum Sq Mean Sq F value Pr(>F)
## set 1 0.5339 0.5339 7.2120 0.0099144 **
## set:replication 2 2.9894 1.4947 20.1914 4.335e-07 ***
## set:male 4 22.1711 5.5428 74.8743 < 2.2e-16 ***
## set:male:female 6 4.8250 0.8042 10.8630 1.311e-07 ***
## set:replication:male:female 10 3.2072 0.3207 4.3325 0.0002462 ***
## Residuals 48 3.5533 0.0740
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## CV: 8.286715 Mean: 3.283333se observa una interacción desde el nivel mas bajo donde el cruce de las hembras con el macho y sus repeticiones muestran una interaccion entre ellas, esto debe ser asi ya que se esta haciendo un cruce de 4 hembras con el mismo macho y ahi es donde se refleja la interacción
output[][-1]## $var.m
## [1] 0.3948843
##
## $var.f
## [1] 0.08057407
##
## $var.A
## [1] 1.579537
##
## $var.D
## [1] -1.257241En este primero se analiza la varianza de cada una de las hembras con las que se cruza el unico macho se observa que la varianza es baja en todas las hembras, lo que me esta hablando de una homogeneidad en los resultados esto puede hablar de una posible dominancia de los genes en el macho
Para el segundo analisis se toman machos y hembras de diferentes razas y se cruzan entre ellos para medir la variabilidad entre las especies
carolina2 <- DC$carolina2
# str(carolina2)
majes<-subset(carolina2,carolina2[,1]==1)
majes<-majes[,c(2,5,4,3,6:8)]
output<-carolina(model=2,majes[,c(1:4,6)])## Response(y): yield
##
## Analysis of Variance Table
##
## Response: y
## Df Sum Sq Mean Sq F value Pr(>F)
## set 1 847836 847836 45.6296 1.097e-09 ***
## set:replication 4 144345 36086 1.9421 0.109652
## set:male 8 861053 107632 5.7926 5.032e-06 ***
## set:female 8 527023 65878 3.5455 0.001227 **
## set:male:female 32 807267 25227 1.3577 0.129527
## Residuals 96 1783762 18581
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## CV: 19.08779 Mean: 714.1301se observa que esta vez la interacción entre el cruce de machos y hembras de diferente especie no presentan una interación entre ellos
output[][-1]## $var.m
## [1] 2746.815
##
## $var.f
## [1] 1355.024
##
## $var.mf
## [1] 2215.415
##
## $var.Am
## [1] 10987.26
##
## $var.Af
## [1] 5420.096
##
## $var.D
## [1] 8861.659En este caso observamos variablilidades muy altas esto debido que el analisis se hace en base a las crias de las hembras sin tener en cuenta el macho que las monto para la variación m, f, D, tambien se realizo un analisis sobre la variabilidad de crias comparando las crias de ambos grupos mf en esta se observa que la varianza es mayor que en las demas.
Este mide la dominancia de los genes
carolina3 <- DC$carolina3
# str(carolina3)
output<-carolina(model=3,carolina3)## Response(y): yield
##
## Analysis of Variance Table
##
## Response: y
## Df Sum Sq Mean Sq F value Pr(>F)
## set 3 2.795 0.93167 1.2784 0.300965
## set:replication 4 3.205 0.80125 1.0995 0.376215
## set:female 4 1.930 0.48250 0.6621 0.623525
## set:male 12 20.970 1.74750 2.3979 0.027770 *
## set:female:male 12 27.965 2.33042 3.1978 0.005493 **
## Residuals 28 20.405 0.72875
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## CV: 21.95932 Mean: 3.8875output[][-1]## $var.mi
## [1] 0.8008333
##
## $var.m
## [1] 0.2546875
##
## $var.A
## [1] 1.01875
##
## $var.D
## [1] 1.601667Ejercicio de youden
para este montaje se esta interesado en estudiar el rendimiento (en kg/parcela) de cuatro tipos de semillas de avena, y se necesita eliminar estadísticamente el efecto del relieve, en este caso la pendiente, y del pH del suelo (los cuales actúan como bloques).en este solo se presentan tres gradientes de pendiente. Para analizar el experimento se usó un cuadrado de Youden con cuatro filas, los pH (pH1, pH2, pH3, pH4); tres columnas, los tres gradientes de la pendiente (P1, P2, P3) y cuatro letras latinas que serían los tipos de semilla (G,D,L,E).Se aleatorizó.
library(readr)
avena <- read_delim("C:/Users/felip/OneDrive/Escritorio/Diseño y experimentos/Parcial/avena.txt",
"\t", escape_double = FALSE, trim_ws = TRUE)## Warning: Missing column names filled in: 'X1' [1]##
## -- Column specification --------------------------------------------------------
## cols(
## X1 = col_double(),
## Rendimiento = col_double(),
## pH = col_character(),
## Pendiente = col_character(),
## Semilla = col_character()
## )View(avena)Youden = data.frame(avena)
Youden## X1 Rendimiento pH Pendiente Semilla
## 1 1 1210 pH1 P1 G
## 2 2 1335 pH2 P1 D
## 3 3 1442 pH3 P1 L
## 4 4 1129 pH4 P1 E
## 5 5 1450 pH1 P2 L
## 6 6 1337 pH2 P2 G
## 7 7 1120 pH3 P2 E
## 8 8 1240 pH4 P2 D
## 9 9 1116 pH1 P3 E
## 10 10 1463 pH2 P3 L
## 11 11 1256 pH3 P3 D
## 12 12 1320 pH4 P3 GAhora se transformará la columna de los tratamientos y bloques en factores para realizar los cálculos.
Youden$Rendimiento=factor(Youden$Rendimiento)
Youden$Rendimiento## [1] 1210 1335 1442 1129 1450 1337 1120 1240 1116 1463 1256 1320
## Levels: 1116 1120 1129 1210 1240 1256 1320 1335 1337 1442 1450 1463Youden$pH = factor(Youden$pH)
Youden$pH## [1] pH1 pH2 pH3 pH4 pH1 pH2 pH3 pH4 pH1 pH2 pH3 pH4
## Levels: pH1 pH2 pH3 pH4Youden$Pendiente = factor(Youden$Pendiente)
Youden$Pendiente## [1] P1 P1 P1 P1 P2 P2 P2 P2 P3 P3 P3 P3
## Levels: P1 P2 P3Youden$Semilla = factor(Youden$Semilla)
Youden$Semilla## [1] G D L E L G E D E L D G
## Levels: D E G LPara cada factor se realizará su respectivo Análisis de Varianza.
Factor Principal:Semilla
Primero se introducen los factores bloques y después los tratamientos.
modelo <- aov(Rendimiento~ pH + Pendiente+ Semilla, data = avena)
summary(modelo)## Df Sum Sq Mean Sq F value Pr(>F)
## pH 3 37855 12618 5.887 0.0897 .
## Pendiente 2 212 106 0.049 0.9525
## Semilla 3 134102 44701 20.853 0.0164 *
## Residuals 3 6431 2144
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Como se puede observar en la tabla el p.valor de Semilla 0.0164 dio menor a 0.05 o al 5% lo que significa que este factor principal o los efectos de los tratamientos (tipo de semilla) son significativos y por ello se puede decir que hay efecto en ella.
Factor Bloque:pH
Para evaluar el efecto del primer bloque, la suma de cuadrados de bloques debe ajustarse por los tratamientos, por lo tanto primero se colocan los tratamientos y después los factores bloques.
mod2=aov(Rendimiento~Semilla+Pendiente+pH, data=avena)
summary(mod2)## Df Sum Sq Mean Sq F value Pr(>F)
## Semilla 3 163606 54535 25.441 0.0123 *
## Pendiente 2 212 106 0.049 0.9525
## pH 3 8350 2783 1.299 0.4176
## Residuals 3 6431 2144
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Como se puede observar en la tabla el p.valor del pH dio 0.4176 mayor a 0.05 o al 5%, lo que significa que el efecto del factor bloque de pH no es significativo. Esto no es de gran relevancia debido a que lo que interesa es el efecto de los tratamientos (semilla) y no de los bloques.Sin embargo se hizo su análisis para que las respuestas futuras puedan ser analizadas
Factor Bloque:Pendiente
Luego, para evaluar el efecto del segundo bloque, la suma de cuadrados de bloques debe ajustarse por los tratamientos, por lo tanto primero se colocan los tratamientos y después los bloques.
mod3=aov(Rendimiento~Semilla+pH+Pendiente, data=avena)
summary(mod3)## Df Sum Sq Mean Sq F value Pr(>F)
## Semilla 3 163606 54535 25.441 0.0123 *
## pH 3 8350 2783 1.299 0.4176
## Pendiente 2 212 106 0.049 0.9525
## Residuals 3 6431 2144
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Como se puede observar en la tabla el p.valor de la pendiente dio 0.9525, mayor al nivel de significancia del 5%, lo que se deduce que el efecto del factor bloque pendiente no es significativo.Esto no es de gran relevancia debido a que lo que interesa es el efecto de los tratamientos (semilla) y no de los bloques.Sin embargo se hizo su análisis para más adelante sacar su conclusión.
Ahora un análisis general
modC=lm(aov(Rendimiento~Semilla+pH+Pendiente, data=avena))
summary(modC)##
## Call:
## lm(formula = aov(Rendimiento ~ Semilla + pH + Pendiente, data = avena))
##
## Residuals:
## 1 2 3 4 5 6 7 8 9 10
## -36.750 30.625 8.125 -2.000 28.625 6.125 -6.125 -28.625 8.125 -36.750
## 11 12
## -2.000 30.625
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1228.00 46.30 26.523 0.000118 ***
## SemillaE -129.88 40.10 -3.239 0.047884 *
## SemillaG 18.75 40.10 0.468 0.671924
## SemillaL 185.62 40.10 4.630 0.018982 *
## pHpH2 76.38 40.10 1.905 0.152906
## pHpH3 20.25 40.10 0.505 0.648297
## pHpH4 32.88 40.10 0.820 0.472345
## PendienteP2 7.75 32.74 0.237 0.828113
## PendienteP3 9.75 32.74 0.298 0.785277
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 46.3 on 3 degrees of freedom
## Multiple R-squared: 0.964, Adjusted R-squared: 0.868
## F-statistic: 10.04 on 8 and 3 DF, p-value: 0.04201Como se puede observar, y con lo presentado anteriormente, se deben tener en cuenta los p.valores del factor semilla (tratamientos) los cuales dieron significativos (el rendimiento de las semillas de avena se comportó diferente) y no los p.valores de los bloques. Esto último lo que indica es que valió la pena bloquear por pH y pendiente. Los p.valores de las semillas E y L dieron significativos.
Revisión de supuestos
Prueba de residuales
modC$residuals## 1 2 3 4 5 6 7 8 9 10
## -36.750 30.625 8.125 -2.000 28.625 6.125 -6.125 -28.625 8.125 -36.750
## 11 12
## -2.000 30.625shapiro.test(modC$residuals)##
## Shapiro-Wilk normality test
##
## data: modC$residuals
## W = 0.90208, p-value = 0.1687plot(modC$residuals)
Como se puede observar los residuales son normales y no hay dependencia espacial ni temporal debido a que el p.valor dio mayor al nivel de significancia del 5% ademas se observa en el grafico una dispersion de estos lo que comprueba su independencia
Igualdad de varianzas
bartlett.test(modC$residuals~avena$Semilla)##
## Bartlett test of homogeneity of variances
##
## data: modC$residuals by avena$Semilla
## Bartlett's K-squared = 3.1685, df = 3, p-value = 0.3664Como se puede observar el p.valor dio mayor al nivel de significancia del 5%, de ahí que las varianzas son estadísticamente iguales. Se cumple el supuesto de homocedasticidad.
Pruebas de comparaciones multiples
Sirve para probar las diferencias existentes entre las medias de los tratamientos
mod1 = aov(avena$Rendimiento ~ avena$Semilla)
mod1## Call:
## aov(formula = avena$Rendimiento ~ avena$Semilla)
##
## Terms:
## avena$Semilla Residuals
## Sum of Squares 163606.33 14993.33
## Deg. of Freedom 3 8
##
## Residual standard error: 43.29165
## Estimated effects may be unbalancedsmod1 = summary(mod1)
tukey =TukeyHSD(mod1)
tukey## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = avena$Rendimiento ~ avena$Semilla)
##
## $`avena$Semilla`
## diff lwr upr p adj
## E-D -155.3333 -268.52841 -42.13826 0.0099008
## G-D 12.0000 -101.19508 125.19508 0.9855486
## L-D 174.6667 61.47159 287.86174 0.0049679
## G-E 167.3333 54.13826 280.52841 0.0064248
## L-E 330.0000 216.80492 443.19508 0.0000657
## L-G 162.6667 49.47159 275.86174 0.0075885plot(tukey, col="red", las=1,cex.axis=0.5, cex.lab=0.5, cex=0.5)
Se observa que todas las comparaciones (5), menos la comparación entre semillas G y D, fueron significativas lo que indica que hay diferencias en el rendimiento entre esos tipos de semillas comparadas.
Ejercicio de parcelas divididas
para este ejercicio tomamos un ejercicio trabajado previamente del cual se le incluyo un codigo de la biblioteca agricolae para su montaje
library(readr)
datos <- read_delim("~/JEP/oats (2).csv",
";", escape_double = FALSE, trim_ws = TRUE)##
## -- Column specification --------------------------------------------------------
## cols(
## observation = col_double(),
## replicate = col_character(),
## variety = col_character(),
## nitro = col_double(),
## yield = col_double()
## )View(datos)
sp.dato <- within(datos, nitroF <- factor(nitro))
head(sp.dato)## # A tibble: 6 x 6
## observation replicate variety nitro yield nitroF
## <dbl> <chr> <chr> <dbl> <dbl> <fct>
## 1 1 I Victory 0 111 0
## 2 2 I Victory 0.2 130 0.2
## 3 3 I Victory 0.4 157 0.4
## 4 4 I Victory 0.6 174 0.6
## 5 5 I Golden.rain 0 117 0
## 6 6 I Golden.rain 0.2 114 0.2llamo a las variables para poder hacer mi grafico de arbol y hacer un mejor analisis sobre si es completo o incompleto
yield = datos$yield
variety= datos$variety
nitroF= datos$nitrodfsp<-data.frame(yield = sp.dato$yield,variety =sp.dato$variety, nitroF = sp.dato$nitroF)
dfsp## yield variety nitroF
## 1 111 Victory 0
## 2 130 Victory 0.2
## 3 157 Victory 0.4
## 4 174 Victory 0.6
## 5 117 Golden.rain 0
## 6 114 Golden.rain 0.2
## 7 161 Golden.rain 0.4
## 8 141 Golden.rain 0.6
## 9 105 Marvellous 0
## 10 140 Marvellous 0.2
## 11 118 Marvellous 0.4
## 12 156 Marvellous 0.6
## 13 61 Victory 0
## 14 91 Victory 0.2
## 15 97 Victory 0.4
## 16 100 Victory 0.6
## 17 70 Golden.rain 0
## 18 108 Golden.rain 0.2
## 19 126 Golden.rain 0.4
## 20 149 Golden.rain 0.6
## 21 96 Marvellous 0
## 22 124 Marvellous 0.2
## 23 121 Marvellous 0.4
## 24 144 Marvellous 0.6
## 25 68 Victory 0
## 26 64 Victory 0.2
## 27 112 Victory 0.4
## 28 86 Victory 0.6
## 29 60 Golden.rain 0
## 30 102 Golden.rain 0.2
## 31 89 Golden.rain 0.4
## 32 96 Golden.rain 0.6
## 33 89 Marvellous 0
## 34 129 Marvellous 0.2
## 35 132 Marvellous 0.4
## 36 124 Marvellous 0.6
## 37 74 Victory 0
## 38 89 Victory 0.2
## 39 81 Victory 0.4
## 40 122 Victory 0.6
## 41 64 Golden.rain 0
## 42 103 Golden.rain 0.2
## 43 132 Golden.rain 0.4
## 44 133 Golden.rain 0.6
## 45 70 Marvellous 0
## 46 89 Marvellous 0.2
## 47 104 Marvellous 0.4
## 48 117 Marvellous 0.6
## 49 62 Victory 0
## 50 90 Victory 0.2
## 51 100 Victory 0.4
## 52 116 Victory 0.6
## 53 80 Golden.rain 0
## 54 82 Golden.rain 0.2
## 55 94 Golden.rain 0.4
## 56 126 Golden.rain 0.6
## 57 63 Marvellous 0
## 58 70 Marvellous 0.2
## 59 109 Marvellous 0.4
## 60 99 Marvellous 0.6
## 61 53 Victory 0
## 62 74 Victory 0.2
## 63 118 Victory 0.4
## 64 113 Victory 0.6
## 65 89 Golden.rain 0
## 66 82 Golden.rain 0.2
## 67 86 Golden.rain 0.4
## 68 104 Golden.rain 0.6
## 69 97 Marvellous 0
## 70 99 Marvellous 0.2
## 71 119 Marvellous 0.4
## 72 121 Marvellous 0.6library(collapsibleTree)
collapsibleTree::collapsibleTree(dfsp,c("variety","nitroF"))vemos que se usan la mismas cantidades de nitrogeno en las tres variedades, pero estas posiblemente se veran afectadas por la capacidad de la variedad a absorber el nitrogeno, por ello me atrevo a decir que estamos hablando de un modelo 4^3 en bloques completamente aleatorizado factorial incompleto
library(lattice) # Can only list one package at a time
library(car)## Loading required package: carDatalibrary(agricolae)
with(sp.dato, xyplot(yield ~ nitroF | variety))
## En esta grafica podemos ver que se toman las tres variedades y se evalua su rendimiento en campo con la aplicación de 4 dosis diferentes de Nitrogeno, en donde por la grafica se puede ver que la variedad Golden rain optine su mejor rendimiento con la aplicacion de 0.4, en Marvvellous se ve al aplicar O.6 y en la victory se ve al aplicar en 0.6, se observa que victory obtuvo los mejores rendimientos comparadas con las tres, pero a su vez fue la de rendimientos mas bajos cuando no se aplico nitrogeno
Primero definimos los tratamientos y generamos datos
Trt 1 : Variedad
Victory (V)
Golden rain (Gr)
Marvellous (M)
Trt 2: cantidad de nitrogeno
0.0 (a)
0.2 (b)
0.4 (c)
0.6 (d)
trt1<-c("V","Gr","M")
trt2<-c("a","b","c", "d")luego usamos la función design.split eta permite encontrar el plan experimental para este diseño
library(agricolae)
diseño <-design.split(trt1,trt2,r=3,seed=2020, randomization = TRUE, first = TRUE)
diseño## $parameters
## $parameters$design
## [1] "split"
##
## $parameters[[2]]
## [1] TRUE
##
## $parameters$trt1
## [1] "V" "Gr" "M"
##
## $parameters$applied
## [1] "rcbd"
##
## $parameters$r
## [1] 3
##
## $parameters$serie
## [1] 2
##
## $parameters$seed
## [1] 2020
##
## $parameters$kinds
## [1] "Super-Duper"
##
##
## $book
## plots splots block trt1 trt2
## 1 101 1 1 Gr a
## 2 101 2 1 Gr c
## 3 101 3 1 Gr b
## 4 101 4 1 Gr d
## 5 102 1 1 M d
## 6 102 2 1 M b
## 7 102 3 1 M c
## 8 102 4 1 M a
## 9 103 1 1 V d
## 10 103 2 1 V c
## 11 103 3 1 V a
## 12 103 4 1 V b
## 13 104 1 2 Gr d
## 14 104 2 2 Gr a
## 15 104 3 2 Gr b
## 16 104 4 2 Gr c
## 17 105 1 2 V d
## 18 105 2 2 V a
## 19 105 3 2 V b
## 20 105 4 2 V c
## 21 106 1 2 M b
## 22 106 2 2 M a
## 23 106 3 2 M c
## 24 106 4 2 M d
## 25 107 1 3 Gr a
## 26 107 2 3 Gr c
## 27 107 3 3 Gr b
## 28 107 4 3 Gr d
## 29 108 1 3 M c
## 30 108 2 3 M d
## 31 108 3 3 M a
## 32 108 4 3 M b
## 33 109 1 3 V b
## 34 109 2 3 V c
## 35 109 3 3 V d
## 36 109 4 3 V aAhora contamos con un modelo experimental de parcela divididas de dos factores
libro <- diseño$book
print(libro)## plots splots block trt1 trt2
## 1 101 1 1 Gr a
## 2 101 2 1 Gr c
## 3 101 3 1 Gr b
## 4 101 4 1 Gr d
## 5 102 1 1 M d
## 6 102 2 1 M b
## 7 102 3 1 M c
## 8 102 4 1 M a
## 9 103 1 1 V d
## 10 103 2 1 V c
## 11 103 3 1 V a
## 12 103 4 1 V b
## 13 104 1 2 Gr d
## 14 104 2 2 Gr a
## 15 104 3 2 Gr b
## 16 104 4 2 Gr c
## 17 105 1 2 V d
## 18 105 2 2 V a
## 19 105 3 2 V b
## 20 105 4 2 V c
## 21 106 1 2 M b
## 22 106 2 2 M a
## 23 106 3 2 M c
## 24 106 4 2 M d
## 25 107 1 3 Gr a
## 26 107 2 3 Gr c
## 27 107 3 3 Gr b
## 28 107 4 3 Gr d
## 29 108 1 3 M c
## 30 108 2 3 M d
## 31 108 3 3 M a
## 32 108 4 3 M b
## 33 109 1 3 V b
## 34 109 2 3 V c
## 35 109 3 3 V d
## 36 109 4 3 V aObservamos el diseño de parcelas y sus respectivos tratamientos de forma grafica en matrices
p<-libro$trt1[seq(1,36,3)]
q<-NULL
for(i in 1:12)
q <- c(q,paste(libro$trt2[3*(i-1)+1],libro$trt2[3*(i-1)+2], libro$trt2[3*(i-1)+3]))
print(t(matrix(p,c(4,3))))## [,1] [,2] [,3] [,4]
## [1,] "Gr" "Gr" "M" "V"
## [2,] "Gr" "Gr" "V" "M"
## [3,] "Gr" "Gr" "M" "V"print(t(matrix(q,c(4,3))))## [,1] [,2] [,3] [,4]
## [1,] "a c b" "d d b" "c a d" "c a b"
## [2,] "d a b" "c d a" "b c b" "a c d"
## [3,] "a c b" "d c d" "a b b" "c d a"res.bad <- lm(yield ~ variety * nitroF, data = sp.dato)
anova(res.bad)## Analysis of Variance Table
##
## Response: yield
## Df Sum Sq Mean Sq F value Pr(>F)
## variety 2 1786.4 893.2 1.7949 0.1750
## nitroF 3 20020.5 6673.5 13.4108 8.367e-07 ***
## variety:nitroF 6 321.7 53.6 0.1078 0.9952
## Residuals 60 29857.3 497.6
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1A raíz del analisis se puede decir que no hay interacción entre la variedad y el nitrogeno claro esta que esta respuesta se debe a que el anialisis esta mal desarrollado
A continuación se realiza el anova correctamente
res.good <- aov(yield ~ variety * nitroF + Error(replicate:variety), data = sp.dato)## Warning in aov(yield ~ variety * nitroF + Error(replicate:variety), data =
## sp.dato): Error() model is singularsummary(res.good)##
## Error: replicate:variety
## Df Sum Sq Mean Sq F value Pr(>F)
## variety 2 1786 893.2 0.612 0.555
## Residuals 15 21889 1459.2
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## nitroF 3 20020 6673 37.686 2.46e-12 ***
## variety:nitroF 6 322 54 0.303 0.932
## Residuals 45 7969 177
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1pero nuevamente se observa como no hay interaccion en variedad y nitrogeno
res.good2 <- aov(yield ~ variety * nitroF + replicate:variety, data = sp.dato)
summary(res.good2)## Df Sum Sq Mean Sq F value Pr(>F)
## variety 2 1786 893 5.044 0.0106 *
## nitroF 3 20020 6673 37.686 2.46e-12 ***
## variety:nitroF 6 322 54 0.303 0.9322
## variety:replicate 15 21889 1459 8.240 1.61e-08 ***
## Residuals 45 7969 177
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Haciendo un anova que incluye esta vez a las replicas vemos como para este caso si se observa una interacción entre las variedades y las replicas
plot(res.good2, 1)
los residuales son independientes, mostrando que no hay dependencia ni espacial ni temporal
plot(res.good2, 2)
la distribución de datos es normal
library(ggplot2)
library(dplyr)##
## 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, union sp.dato %>%
group_by(variety,nitroF) %>%
summarise(medsp.dato = mean(yield)) -> medias## `summarise()` regrouping output by 'variety' (override with `.groups` argument) ggplot(medias)+
aes (x= nitroF, y = medsp.dato, color = variety) +
geom_line(aes(group=variety)) +
geom_point()
library(ggplot2)
library(dplyr)
sp.dato %>%
group_by(variety,nitroF) %>%
summarise(medsp.dato = mean(yield)) -> medias## `summarise()` regrouping output by 'variety' (override with `.groups` argument) ggplot(medias)+
aes (x= variety, y = medsp.dato, color = nitroF) +
geom_line(aes(group=nitroF)) +
geom_point()
pru = with(sp.dato, HSD.test(yield, nitroF, DFerror = 45, MSerror = 117))
pru## $statistics
## MSerror Df Mean CV MSD
## 117 45 103.9722 10.40341 9.618527
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey nitroF 4 3.772697 0.05
##
## $means
## yield std r Min Max Q25 Q50 Q75
## 0 79.38889 19.39417 18 53 117 63.25 72.0 94.25
## 0.2 98.88889 21.84407 18 64 140 83.75 95.0 112.50
## 0.4 114.22222 22.31738 18 81 161 97.75 115.0 124.75
## 0.6 123.38889 22.99908 18 86 174 106.25 121.5 139.00
##
## $comparison
## NULL
##
## $groups
## yield groups
## 0.6 123.38889 a
## 0.4 114.22222 a
## 0.2 98.88889 b
## 0 79.38889 c
##
## attr(,"class")
## [1] "group"Diseño en Cuadrados Latinos
Los diseños en cuadrados latinos utilizan dos variables de bloque para reducir el error experimental.
Se quiere evaluar la productividad de cuatro variedades de aguacates, A, B, C y D. Para ello se decide realizar el ensayo en un terreno que posee un gradiente de pendiente de oriente a occidente y además, diferencias en la disponibilidad de Nitrógeno de norte a sur, para controlar los efectos de la pendiente y la disponibilidad de Nitrógeno, se utilizó un diseño de cuadrado latino, los datos corresponden a la producción en kg/parcela.
library(readr)
cuadrado_lat<- read_delim("C:/Users/felip/Downloads/AGUACATE.csv",
";", escape_double = FALSE, trim_ws = TRUE)##
## -- Column specification --------------------------------------------------------
## cols(
## productividad = col_double(),
## aguacate = col_character(),
## nitrogeno = col_character(),
## pendiente = col_character()
## )View(cuadrado_lat)El análisis de la productividad de las variedades de aguacate corresponde al análisis de un factor con 4 niveles. Dado que en el estudio intervienen dos fuentes de variación: la Disponibilidad de Nitrógeno y la Pendiente, se consideran dos factores de bloque, cada uno de ellos con 4 niveles.
• Variable respuesta: Productividad
• Factor: Variedad de aguacate. Es un factor de efectos fijos ya que desde el principio se establecen los niveles concretos que se van a analizar.
• Bloques: Disponibilidad de Nitrógeno y Pendiente, ambos con 4 niveles y ambos de efectos fijos.
• Tamaño del experimento: Número total de observaciones (42).
cuadrado_lat$latina <- factor(cuadrado_lat$aguacate)
cuadrado_lat$latina## [1] D A C B A B D C C D B A B C A D
## Levels: A B C Dcuadrado_lat$Bloque1 <- factor(cuadrado_lat$nitrogeno)
cuadrado_lat$Bloque1## [1] N1 N1 N1 N1 N2 N2 N2 N2 N3 N3 N3 N3 N4 N4 N4 N4
## Levels: N1 N2 N3 N4cuadrado_lat$Bloque2 = factor(cuadrado_lat$pendiente)
cuadrado_lat$Bloque2## [1] P1 P2 P3 P4 P1 P2 P3 P4 P1 P2 P3 P4 P1 P2 P3 P4
## Levels: P1 P2 P3 P4mod1 <- aov(productividad~ aguacate + nitrogeno + pendiente, data = cuadrado_lat )
mod1## Call:
## aov(formula = productividad ~ aguacate + nitrogeno + pendiente,
## data = cuadrado_lat)
##
## Terms:
## aguacate nitrogeno pendiente Residuals
## Sum of Squares 5556.25 92518.75 52556.25 112.50
## Deg. of Freedom 3 3 3 6
##
## Residual standard error: 4.330127
## Estimated effects may be unbalancedsummary(mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## aguacate 3 5556 1852 98.78 1.70e-05 ***
## nitrogeno 3 92519 30840 1644.78 3.92e-09 ***
## pendiente 3 52556 17519 934.33 2.13e-08 ***
## Residuals 6 113 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(agricolae)Observando los valores de los p-valores, 1.70e-05, 3.92e-09 y 2.13e-08; menores respectivamente que el nivel de significación del 5%, deducimos que los tres efectos son significativos. Tanto las variedades de aguacates utilizadas, como la pendiente del terreno y la disponibilidad de nitrógeno influyen en la productividad de los aguacates.
Los supuestos que han de verificarse en un diseño de cuadrados latinos son Normalidad, Homocedasticidad e Independencia además del supuesto de aditividad entre filas, columnas y tratamientos (es decir, que no haya interacciones entre los mismos).
Hipótesis de normalidad
Comprobamos la hipótesis de normalidad mediante el análisis de la normalidad de los residuos. Para ello, hacemos uso del test de Shapiro-Wilks:
shapiro.test(mod1$residuals)##
## Shapiro-Wilk normality test
##
## data: mod1$residuals
## W = 0.89854, p-value = 0.07616Observamos el contraste de Shapiro-Wilk que es adecuado cuando las muestras son pequeñas (n<50) y es una alternativa más potente que el test de Kolmogorov-Smirnov. El p-valor (0.07616) es mayor que el nivel de significación del 5%, aceptándose la hipótesis de normalidad.
Homogeneidad de varianzas
En este caso hacemos uso del Test de Barlett para contrastar la igualdad entre varianzas del factor.
bartlett.test(cuadrado_lat$productividad, cuadrado_lat$aguacate)##
## Bartlett test of homogeneity of variances
##
## data: cuadrado_lat$productividad and cuadrado_lat$aguacate
## Bartlett's K-squared = 3.355, df = 3, p-value = 0.3401bartlett.test(cuadrado_lat$productividad, cuadrado_lat$nitrogeno)##
## Bartlett test of homogeneity of variances
##
## data: cuadrado_lat$productividad and cuadrado_lat$nitrogeno
## Bartlett's K-squared = 0.56609, df = 3, p-value = 0.9041bartlett.test(cuadrado_lat$productividad, cuadrado_lat$pendiente)##
## Bartlett test of homogeneity of variances
##
## data: cuadrado_lat$productividad and cuadrado_lat$pendiente
## Bartlett's K-squared = 0.4168, df = 3, p-value = 0.9368El p-valor es 0.3401 por lo tanto no se puede rechazar la hipótesis de homogeneidad de las varianzas y se concluye que la cuatro variedades tienen varianzas homogéneas.
Los p-valores son mayores que 0.05, por lo tanto no se puede rechazar la hipótesis de homogeneidad de las varianzas en ninguno de los factores bloques.
Independencia de los residuos
Para ello tenemos que representar un gráfico de los residuos tipificados frente a los pronosticados. En R obtenemos varios gráficos a la vez que están incluidos en la estimación del modelo.
Para verlos de forma correcta hacemos uso de las siguientes órdenes:
layout(matrix(c(1,2,3,4),2,2)) # para que salgan en la misma pantalla
plot(mod1)## hat values (leverages) are all = 0.625
## and there are no factor predictors; no plot no. 5
Nos fijamos en el primer gráfico que representa los valores ajustados frente a los residuos y observamos que no hay ninguna tendencia sistemática. Concluimos que no hay sospechas para que se incumpla la hipótesis de independencia.
Analizamos si Se obtiene la misma producción con las cuatro variedades de aguacate, en caso negativo, analizar mediante el procedimiento de Tukey y Newman-Keuls, con qué variedad de aguacate hay mayor producción.
TukeyHSD(mod1, "aguacate", ordered = TRUE)## Tukey multiple comparisons of means
## 95% family-wise confidence level
## factor levels have been ordered
##
## Fit: aov(formula = productividad ~ aguacate + nitrogeno + pendiente, data = cuadrado_lat)
##
## $aguacate
## diff lwr upr p adj
## A-B 33.75 23.1507168 44.34928 0.0001393
## D-B 38.75 28.1507168 49.34928 0.0000629
## C-B 50.00 39.4007168 60.59928 0.0000134
## D-A 5.00 -5.5992832 15.59928 0.4289199
## C-A 16.25 5.6507168 26.84928 0.0072734
## C-D 11.25 0.6507168 21.84928 0.0392225plot(TukeyHSD(mod1, "aguacate"))
### Así por ejemplo, si comparamos la concentración media correspondiente a las variedades A Y B, tenemos una diferencia entre ambas medias de 33.75, un p-valor (Sig.) de 0.0001393 significativo. Por lo tanto, las producciones medias de las variedades A y B pueden considerarse distintas estadísticamente y un intervalo de confianza con un límite inferior 23.1507168 y un límite superior 44.34928 y por lo tanto no contiene al cero de lo que también deducimos que hay diferencias significativas entre estas dos variedades de aguacate. ### Como se puede observar, todos los intervalos de confianza construidos para las diferencias entre las producciones medias de las variedades no contienen al 0, excepto el correspondiente a la pareja de variedades de aguacates A y D. Lo que significa que todas las producciones medias pueden considerarse distintas estadísticamente excepto las producciones medias correspondientes a las variedades A y D. Se deduce que únicamente no se observan diferencias significativas entre las producciones de las variedades de aguacates A y D (P-valor = 0.4289199).
Newman_Keuls <- SNK.test(mod1,"aguacate", console=TRUE)##
## Study: mod1 ~ "aguacate"
##
## Student Newman Keuls Test
## for productividad
##
## Mean Square Error: 18.75
##
## aguacate, means
##
## productividad std r Min Max
## A 811.25 68.84463 4 730 880
## B 777.50 143.38177 4 595 945
## C 827.50 141.50972 4 700 950
## D 816.25 55.43389 4 760 885
##
## Alpha: 0.05 ; DF Error: 6
##
## Critical Range
## 2 3 4
## 7.492106 9.394633 10.599283
##
## Means with the same letter are not significantly different.
##
## productividad groups
## C 827.50 a
## D 816.25 b
## A 811.25 b
## B 777.50 cEl contraste de Newman-keuls muestra una tabla de Subconjuntos homogéneos. En nuestro estudio sobre las producciones de aguacates se observan que hay tres subgrupos homogéneos, al primer subgrupo pertenece la Variedad C, con una producción de 827.50 media kg/parcela, el segundo las variedades D y A, con producciones medias de 816.25 y 811.25 Kg(parcela, respectivamente y el tercero la Variedad B, con una producción media de 777.50 kg/parcela. Por lo tanto, se observa que la producción media mayor se obtiene con la Variedad C (827.5 Kg/ parcela) y la menor con la Variedad B (777.50 Kg/parcela).
Diseño en Cuadrado greco latino
El estudio se realiza con cuatro tipos de semillas diferentes que se plantan en cuatro tipos de terreno, se les aplican cuatro tipos de abonos y cuatro tipos de insecticidas. La asignación de los tratamientos a las plantas se realiza de forma aleatoria. Para controlar estas posibles fuentes de variabilidad se decide plantear un diseño por cuadrados greco-latinos como el que se muestra en la siguiente tabla, donde las letras griegas corresponden a los cuatro tipos de semilla y las latinas a los abonos.
• Variable respuesta: Crecimiento
• Factor: Tipo_abonos que tiene cuatro niveles. Es un factor de efectos fijos ya que viene decidido que niveles concretos se van a utilizar.
• Bloques: Insecticidas, Terrenos y Semillas, cada uno con cuatro niveles y de efectos fijos.
• Tamaño del experimento: Número total de observaciones (16).
library(readr)
tomate <- read_delim("C:/Users/felip/Downloads/greco latino.csv",
";", escape_double = FALSE, trim_ws = TRUE)##
## -- Column specification --------------------------------------------------------
## cols(
## Crecimiento = col_double(),
## Abono = col_character(),
## Semilla = col_character(),
## Insecticida = col_character(),
## Terreno = col_character()
## )View(tomate)tomate$Abono <- factor(tomate$Abono)
tomate$Abono## [1] C B A D B C D A D A B C A D C B
## Levels: A B C Dtomate$Semilla <- factor(tomate$Semilla)
tomate$Semilla## [1] beta alfa delta gamma gamma delta alfa beta delta gamma beta alfa
## [13] alfa beta gamma delta
## Levels: alfa beta delta gammatomate$Insecticida <- factor(tomate$Insecticida)
tomate$Insecticida## [1] Insecticida1 Insecticida2 Insecticida3 Insecticida4 Insecticida1
## [6] Insecticida2 Insecticida3 Insecticida4 Insecticida1 Insecticida2
## [11] Insecticida3 Insecticida4 Insecticida1 Insecticida2 Insecticida3
## [16] Insecticida4
## Levels: Insecticida1 Insecticida2 Insecticida3 Insecticida4tomate$Terreno <- factor(tomate$Terreno)
tomate$Terreno## [1] T.1 T.1 T.1 T.1 T.2 T.2 T.2 T.2 T.3 T.3 T.3 T.3 T.4 T.4 T.4 T.4
## Levels: T.1 T.2 T.3 T.4Para afirmar que el crecimiento de las plantas es el mismo para los cuatro tipos de abonos, junto con los distintos insecticidas. Se realizo una tabla anova
mod1 <- aov(Crecimiento~ Abono + Semilla + Insecticida + Terreno, data = tomate )
mod1## Call:
## aov(formula = Crecimiento ~ Abono + Semilla + Insecticida + Terreno,
## data = tomate)
##
## Terms:
## Abono Semilla Insecticida Terreno Residuals
## Sum of Squares 42.25 31.25 20.75 64.25 1.25
## Deg. of Freedom 3 3 3 3 3
##
## Residual standard error: 0.6454972
## Estimated effects may be unbalancedsummary(mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## Abono 3 42.25 14.083 33.8 0.00820 **
## Semilla 3 31.25 10.417 25.0 0.01266 *
## Insecticida 3 20.75 6.917 16.6 0.02260 *
## Terreno 3 64.25 21.417 51.4 0.00445 **
## Residuals 3 1.25 0.417
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Como podemos ver son significativos los efectos de todos los factores
Ahora analizaremos si el crecimiento de las plantas es el mismo para los cuatro tipos de abonos y con los distintos insecticidas, vamos a realizar contrastes de comparaciones múltiples, por ejemplo de Duncan para los tipos de abonos y Newman Keuls para los insecticidas.
library(agricolae)
(duncan=duncan.test(mod1, "Abono" , group = T))## $statistics
## MSerror Df Mean CV
## 0.4166667 3 9.375 6.885304
##
## $parameters
## test name.t ntr alpha
## Duncan Abono 4 0.05
##
## $duncan
## Table CriticalRange
## 2 4.500659 1.452581
## 3 4.515652 1.457420
## 4 4.472854 1.443607
##
## $means
## Crecimiento std r Min Max Q25 Q50 Q75
## A 10.00 3.464102 4 5 13 9.50 11.0 11.50
## B 8.00 3.162278 4 5 12 5.75 7.5 9.75
## C 7.75 1.707825 4 6 10 6.75 7.5 8.50
## D 11.75 3.774917 4 7 16 10.00 12.0 13.75
##
## $comparison
## NULL
##
## $groups
## Crecimiento groups
## D 11.75 a
## A 10.00 b
## B 8.00 c
## C 7.75 c
##
## attr(,"class")
## [1] "group"El mayor crecimiento de las plantas se produce con el Abono D siendo la altura que alcanza de 11.75 y la menor altura (7.75) la alcanza con el tipo de abono C
Newman_Keuls <- SNK.test(mod1,"Insecticida", console=TRUE, main=" Comparación de Newman-Keuls para el factor tipo de insecticida")##
## Study: Comparación de Newman-Keuls para el factor tipo de insecticida
##
## Student Newman Keuls Test
## for Crecimiento
##
## Mean Square Error: 0.4166667
##
## Insecticida, means
##
## Crecimiento std r Min Max
## Insecticida1 7.5 2.380476 4 6 11
## Insecticida2 9.5 3.109126 4 5 12
## Insecticida3 10.5 4.932883 4 5 16
## Insecticida4 10.0 2.581989 4 7 13
##
## Alpha: 0.05 ; DF Error: 3
##
## Critical Range
## 2 3 4
## 1.452581 1.907336 2.202606
##
## Means with the same letter are not significantly different.
##
## Crecimiento groups
## Insecticida3 10.5 a
## Insecticida4 10.0 a
## Insecticida2 9.5 a
## Insecticida1 7.5 bLa menor altura (7.5) la alcanza cuando se le suministra el Insecticida 1.
Ahora analizamos las diferencias significativas en el crecimiento de las plantas con las distintas semillas, y si el tipo de tierra influye en dicho crecimiento.
mod.tukey<- TukeyHSD(mod1, "Semilla", ordered = TRUE)
mod.tukey## Tukey multiple comparisons of means
## 95% family-wise confidence level
## factor levels have been ordered
##
## Fit: aov(formula = Crecimiento ~ Abono + Semilla + Insecticida + Terreno, data = tomate)
##
## $Semilla
## diff lwr upr p adj
## beta-gamma 0.25 -1.9526064 2.452606 0.9412719
## delta-gamma 1.75 -0.4526064 3.952606 0.0901983
## alfa-gamma 3.50 1.2973936 5.702606 0.0139151
## delta-beta 1.50 -0.7026064 3.702606 0.1305850
## alfa-beta 3.25 1.0473936 5.452606 0.0171739
## alfa-delta 1.75 -0.4526064 3.952606 0.0901983Comprobamos que únicamente hay diferencias significativas entre los tipos de semillas alfa-gamma (p-valor = 0.0139) y entre alfa-beta (p-valor =0.0171).
Para determinar si el tipo de terreno influye en el crecimiento aplicamos el método de comparaciones múltiples LSD.
LSD.test(mod1,"Terreno", p.adj="bonferroni", console=TRUE)##
## Study: mod1 ~ "Terreno"
##
## LSD t Test for Crecimiento
## P value adjustment method: bonferroni
##
## Mean Square Error: 0.4166667
##
## Terreno, means and individual ( 95 %) CI
##
## Crecimiento std r LCL UCL Min Max
## T.1 11.00 3.366502 4 9.97287 12.02713 6 13
## T.2 10.75 4.112988 4 9.72287 11.77713 6 16
## T.3 6.00 1.154701 4 4.97287 7.02713 5 7
## T.4 9.75 1.500000 4 8.72287 10.77713 8 11
##
## Alpha: 0.05 ; DF Error: 3
## Critical Value of t: 6.231543
##
## Minimum Significant Difference: 2.844297
##
## Treatments with the same letter are not significantly different.
##
## Crecimiento groups
## T.1 11.00 a
## T.2 10.75 a
## T.4 9.75 a
## T.3 6.00 bHay dos grupos de terrenos que difieren significativamente, por un lado el grupo formado por los tipos de terrenos 1,2 y 4 y por otra parte el grupo formado un solo tipo de terreno, el tipo de terreno 3. El mayor crecimiento de las plantas se produce en el tipo de terreno 1 con un crecimiento de 11 u.c.
Para analiar con qué tipo de semilla se produce el mayor crecimiento de las plantas
Newman_Keuls1 <- SNK.test(mod1,"Semilla", console=TRUE, main=" Contraste de Newman-Keuls para el factor tipo de semilla ")##
## Study: Contraste de Newman-Keuls para el factor tipo de semilla
##
## Student Newman Keuls Test
## for Crecimiento
##
## Mean Square Error: 0.4166667
##
## Semilla, means
##
## Crecimiento std r Min Max
## alfa 11.50 3.696846 4 7 16
## beta 8.25 3.201562 4 5 11
## delta 9.75 2.500000 4 7 13
## gamma 8.00 3.559026 4 5 13
##
## Alpha: 0.05 ; DF Error: 3
##
## Critical Range
## 2 3 4
## 1.452581 1.907336 2.202606
##
## Means with the same letter are not significantly different.
##
## Crecimiento groups
## alfa 11.50 a
## delta 9.75 b
## beta 8.25 c
## gamma 8.00 cEl mayor crecimiento de la planta se produce con el tipo de semilla alfa (11.50 u.c.)
Como ultimo paso se analiza si El crecimiento de las plantas es el mismo utilizando al mismo tiempo los abonos A y B que utilizando los abonos C y D
TukeyHSD(mod1, "Abono", ordered = TRUE)## Tukey multiple comparisons of means
## 95% family-wise confidence level
## factor levels have been ordered
##
## Fit: aov(formula = Crecimiento ~ Abono + Semilla + Insecticida + Terreno, data = tomate)
##
## $Abono
## diff lwr upr p adj
## B-C 0.25 -1.95260644 2.452606 0.9412719
## A-C 2.25 0.04739356 4.452606 0.0472540
## D-C 4.00 1.79739356 6.202606 0.0094879
## A-B 2.00 -0.20260644 4.202606 0.0643491
## D-B 3.75 1.54739356 5.952606 0.0114235
## D-A 1.75 -0.45260644 3.952606 0.0901983Hay diferencias significativas en el crecimiento utilizando los abonos C y D (P-valor = 0.0094), pero no las hay con los abonos A y B (P-valor= 0.064).
Analisis de varianza no parametrica
Se analiza el peso de fruto de tres variedades de arandanos, de fecha de cosecha temprana. En este diseño se quiere ver si existe alguna diferencia entre el peso de las tres variedades de arandanos, para esto se va a emplear un modelo anova donde verificaremos si lo podemos usar, y si es el caso contrario hacemos uso del test no paramétrico de Kruskal-Wallis.
library(readr)
kruskal <- read_delim("C:/Users/felip/Downloads/kruskal.csv",
";", escape_double = FALSE, trim_ws = TRUE)##
## -- Column specification --------------------------------------------------------
## cols(
## variedad = col_double(),
## peso = col_double()
## )View(kruskal)
attach(kruskal)
names(kruskal)## [1] "variedad" "peso"str(kruskal)## tibble [64 x 2] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ variedad: num [1:64] 1 1 1 1 1 1 1 1 1 1 ...
## $ peso : num [1:64] 19 18.6 18.9 17.8 18.5 18.2 18.9 19.2 18 16.9 ...
## - attr(*, "spec")=
## .. cols(
## .. variedad = col_double(),
## .. peso = col_double()
## .. )kruskal$variedad <- as.factor(kruskal$variedad)
kruskal$variedad <- as.factor(kruskal$variedad)
kruskal$variedad = factor(kruskal$variedad,labels = c(" Blue One", "Blue ribon", "Blue moon"))
class(kruskal$variedad)## [1] "factor"###Vemos la estructura de nuestros datos:
kruskal$variedad <- as.factor(kruskal$variedad)
kruskal$variedad <- as.factor(kruskal$variedad)
kruskal$variedad = factor(kruskal$variedad,labels = c(" Blue One", "Blue ribon", "Blue moon"))
class(kruskal$variedad)## [1] "factor"###Representamos nuestros datos con un boxplot:
par(mfrow=c(1,1))
boxplot(peso ~ variedad, data=kruskal, col=c("green","red","blue"),
cex.axis=0.7,las = 2, ylab="Valor", xlab="variedades", cex.lab=0.75)
### Realizamos un modelo anova pero recordamos que uno de los requisitos del anova es que su distribución sea normal, lo cual notenemos. Es por esto que estos datos se consideran erroneos y proseguimos a un modelo kruskal el cual solo tiene como condición la homogeniedad de varianazas.
nuestraanova <- lm(peso ~ variedad, data=kruskal)
nuestromodelo <- anova(nuestraanova)
nuestromodelo## Analysis of Variance Table
##
## Response: peso
## Df Sum Sq Mean Sq F value Pr(>F)
## variedad 2 7.622 3.8109 5.3291 0.007361 **
## Residuals 61 43.622 0.7151
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1###Comprobamos si nuestros residuos son normales, requisito de un ANOVA, haciendoun test de normalidad y finalmente un test de homogeniedad de varianzas.
par(mfrow=c(2,2))
plot(nuestraanova)
shapiro.test(residuals(nuestraanova))##
## Shapiro-Wilk normality test
##
## data: residuals(nuestraanova)
## W = 0.91025, p-value = 0.0001993fligner.test(peso ~ variedad, data = kruskal)##
## Fligner-Killeen test of homogeneity of variances
##
## data: peso by variedad
## Fligner-Killeen:med chi-squared = 1.593, df = 2, p-value = 0.4509library(car)
leveneTest(peso ~ variedad, data = kruskal)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.8461 0.4341
## 61Para seguir con el modelo kruskal verificamos la homogeniedad de varianzas.
library(car)
leveneTest(peso ~ variedad, data = kruskal)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.8461 0.4341
## 61Podemos ver que tiene un valor mayor al de significancia, lo cual verifica que las varianzas son iguales
Ahora realizamos un test no paramétrico de Kruskal-Wallis, y hacemos el siguiente test de comparaciones. Finalmente visualizamos los datos en un ggplot
kruskal.test(peso ~ variedad, data = kruskal)##
## Kruskal-Wallis rank sum test
##
## data: peso by variedad
## Kruskal-Wallis chi-squared = 11.549, df = 2, p-value = 0.003105Como podemos notar el p valor es mucho menor que el número de significancia, donde podemos decir que se rechaza la hipotesis nula, y se dice que hay una diferencia de peso en al menos de una de las variedades de arandanos de nuestro diseño.
Ya que no sabemos que variedad es, realizamos un Test de comparaciones múltiples para saber el desempeño de cada variedad según su peso.
library(FSA)## ## FSA v0.8.31. See citation('FSA') if used in publication.
## ## Run fishR() for related website and fishR('IFAR') for related book.##
## Attaching package: 'FSA'## The following object is masked from 'package:car':
##
## bootCasedunnTest(peso ~ variedad, data= kruskal,method="bonferroni")## Dunn (1964) Kruskal-Wallis multiple comparison## p-values adjusted with the Bonferroni method.## Comparison Z P.unadj P.adj
## 1 Blue One - Blue moon 1.686196 0.0917579808 0.275273942
## 2 Blue One - Blue ribon 3.388499 0.0007027629 0.002108289
## 3 Blue moon - Blue ribon 1.816676 0.0692667647 0.207800294Si notamos la comparación entre blue one y blue ribon, vemos el valor de 0.002108289 que indica una diferencia entre estas dos variedades, esto es lo unico que podemos deducir con estos datos, que existe una diferencia entre estas dos variedades.
library("ggplot2")
ggplot(kruskal, aes(x = variedad, y = peso)) +
geom_boxplot(fill = "grey80", colour = "black") +
scale_x_discrete() + xlab("variedad") +
ylab("peso (gr)")
Ahora que lo visualizamos en el anterior ggplot, notamos que la variedad blue one tienre un mayor peso a comparación de las otras variedades, y sobre todo de la variedad blue ribon de la cual hablamos anteriormente.
Resumen 1
Un experimento de parcela dividida es una opción a considerar, especialmente si tiene una situación en la que algunos factores tienen niveles más difíciles de cambiar que otros. Se debe tener en cuenta que 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. Algunas de las razones para cuantificar los factores difíciles de cambiar podría basarse en el presupuesto o el tiempo.
El análisis del experimento es un punto importante en el cual se debe profundizar, en la lectura propuesta afirman que un experimento de parcela dividida puede tratarse erróneamente como un CRD. Pero esto podría llevar a considerar un factor estadísticamente significativo o no significativo cuando en realidad es importante.
Una recomendación dada, es que en el experimento de parcela dividida su orden se considere cuidadosamente. Esto no significa que omita el paso de asignación al azar, sino que la forma de la asignación al azar sería más estructurada. Es importante analizar los datos como un experimento de parcela dividida. Esto conducirá a conclusiones válidas que brinden información precisa sobre los efectos. Las nuevas versiones de paquetes analíticos de JMP y Design Expert facilitan la selección de un buen experimento de parcela dividida y la realización del análisis adecuado.
Ellos afirmaban que en un experimento de parcela dividida puede producir estimaciones de parámetros de modelo más precisas que el experimento del mismo tamaño ejecutado como un CRD. Y que si se tiene en cuenta el costo del experimento, podría ser aún más fácil justificar este tipo de experimento.
También se deben considerar diferentes diseños dependiendo de los costos relativos de cambiar los factores. En una relación de costo y calidad del diseño lo cual no siempre es apropiado, con frecuencia puede permitir comparaciones más realistas y ayudar a facilitar la toma de decisiones.
Para realizar un diseño de parcela dividida se debe tener en cuenta los siguientes aspectos:
características cualitativas, como el equilibrio del número de observaciones por parcela completa y el número de niveles de cada factor; medidas cuantitativas como una buena estimación de toda la parcela y los términos de error de la subparcela, y rendimiento basado en criterios de optimización. También resaltan lo importante que es considerar que es el óptimo en el experimento y así poder hacer un análisis correcto sobre los objetivos del proyecto y sus criterios. Con esta información y diferentes estrategias podemos aumentar la cantidad de información recopilada más que con otros diseños.
Resumen 2
El diseño sólido de los experimentos y la implementación adecuada de métodos estadísticos apropiados para el análisis de datos son fundamentales para producir resultados científicos significativos que sean tanto replicables como reproducibles . Primero hay que considerar el concepto de unidad estadística, propuesto por Robinson, un término que es vago y que carece de una definición universal en la literatura general sobre el diseño de experimentos, particularmente para aplicaciones agrícolas . En cambio, definamos la unidad experimental y el observa unidad nacional , tanto formalmente como en el contexto específico de las ciencias. En este caso, es claramente la unidad experimental.
Una unidad de observación se define como la entidad física sobre la que se mide un resultado de interés en un experimento . Es la entidad que se asigna de forma independiente al tratamiento y se mide en función de los resultados, entonces sirve como unidad experimental y como unidad de observación. Por otro lado, si el resultado de interés en el ejemplo de dieta se midió en vacas individuales en cada corral, digamos producción de leche, uno encuentra una brecha natural o un desajuste entre la entidad asignada independientemente al tratamiento y la entidad medida . Una estructura de datos jerárquica se refiere a una configuración de los datos donde las observaciones no son independientes entre sí, sino que tienen una estructura de correlación impuesta por el diseño experimental.
Cabe señalar que el caso presentado en la carta de Robinson al editor carece de una descripción clara de cómo se llevó a cabo el experimento, ni tampoco se explicó claramente el proceso de recopilación de datos. Para ambos escenarios, trabajemos en el ejercicio llamado ¿Qué haría Fisher? , que fue introducido por Stroup para traducir la descripción de un diseño experimental a un shell ANOVA a un modelo estadístico real. El ejercicio WWFD es una estrategia general que, en el contexto de datos asumidos como normales, puede demostrarse que es equivalente al ejercicio tradicional ANOVA de cuadrados medios esencial para identificar el error experimental adecuado y así distinguir entre una unidad experimental y una unidad de observación. Para seguir el enfoque de WWFD, es importante tener en cuenta las posiciones relativas de las filas correspondientes a la estructura del tratamiento y las filas correspondientes a los elementos del diseño experimental , así como su combinación para caracterizar adecuadamente el proceso de recolección de datos.
En este caso, la asignación aleatoria del tratamiento se realiza dentro de un período determinado, con la reasignación de tratamientos al comienzo de cada nuevo período. Al ser la unidad experimental para el tratamiento, en un período dado define el nivel de replicación independiente para el tratamiento en la estructura jerárquica de datos y, por lo tanto, identifica el término de error experimental. Para completar, observamos que el término de error de muestreo representa la variación entre las unidades de observación, distinta del error experimental o la variación entre las unidades experimentales. así como la interacción entre el tratamiento y el período, se supone que no existe para permitir la estimación de una medida de error y por lo tanto no se considera
En este caso, se puede, y se debe, investigar la interacción entre el tratamiento y el período, siempre que se consideradas dentro de cada cuadro, particularmente si el período refleja un evento fisiológico . Además, se observa que hasta ahora se trata el período como un elemento de diseño experimental , como es consistente con la literatura general sobre diseño de experimentos. Sin embargo, específicamente para algunas aplicaciones lácteas, el período puede considerarse legítimamente como un elemento del diseño experimental o como un elemento de la estructura del tratamiento, dependiendo de la variabilidad de los días de leche dentro de un período y, por lo tanto, la etapa de lactancia . Las decisiones sobre cómo tratar el período , así como la incorporación de la interacción tratamiento período en cuadrados latinos replicados, no son triviales y deben hacerse caso por caso.
De hecho, los estudios basados en diseños de cuadrados latinos replicados mostraron que los efectos dietéticos pueden depender de la etapa de lactancia , en cuyo caso la inferencia debe centrarse en las diferencias de tratamiento dentro de los períodos como opuesto a los efectos generales del tratamiento. Además, tratar el período como un elemento del diseño experimental o de la estructura del tratamiento tiene implicaciones para la inferencia posterior porque determina cómo se define el error experimental para el tratamiento de interés . Una vez que se ha completado el ejercicio del WWFD de modo que los procesos de aleatorización y recopilación de datos en los escenarios A y B estén completamente caracterizados, se pueden transferir los elementos de fila de la sección combinada de las tablas del WWFD en un predictor lineal para especificar el modelo lineal correspondiente para cada escenario . La especificación adecuada de un modelo estadístico para el análisis de datos es un proceso riguroso para el que es fundamental tener una sólida comprensión de la estructura jerárquica de los datos mediante una comprensión profunda del proceso de recopilación de datos.
Una descripción clara y detallada de cómo se recopilaron los datos y cómo se implementó el diseño es imperativa para la especificación del modelo. Las implicaciones inferenciales de ignorar el diseño experimental y especificar un modelo estadístico que no coincide con el proceso de recolección de datos no son menores. En otras palabras, especificar incorrectamente el modelo 1 para el escenario A llevaría, en promedio, al investigador a concluir que más diferencias de tratamiento son estadísticamente significativas de lo que debería ser. Sin duda, los cambios sutiles en la forma en que se ejecuta un experimento pueden tener efectos profundos sobre cómo se especifica el modelo estadístico para el análisis de datos .
La conclusión es que el modelo estadístico debe describir un proceso plausible que da lugar a las observaciones capturando las variables independientes importantes que afectan la variable de resultado y especificando cualquier restricción en la aleatorización o cualquier otra estructura de datos que induce la correlación entre observaciones. Además, en el contexto de modelos mixtos, no hay necesidad de colapsar los datos a nivel de animal en resúmenes a nivel de corral, ni para «eliminar» los datos a nivel de animal ni para «destruir» la variación a nivel de vaca, como sugiere Robinson . Como tal, los modelos mixtos pueden reconocer simultáneamente múltiples fuentes de variabilidad aleatoria en un conjunto de datos, evaluando así la variabilidad a nivel de vaca y la variabilidad a nivel de corral al mismo tiempo, y usando uno u otro como error experimental para un tratamiento dado de interés, según corresponda. Así, Se pueden utilizar modelos mixtos para garantizar que la unidad experimental para cada tratamiento de interés se reconozca adecuadamente dentro de un diseño de estudio, lo que lleva a un reconocimiento adecuado del nivel de error experimental y, por lo tanto, a la prueba de hipótesis adecuada.
Producción y rendimiento de maíz en cuatro tipos de labranza bajo condiciones de temporal
a.La mención de la estructura factorial
En el artículo no se menciona la estructura factorial del montaje, pero se puede intuir que se trata de un diseño en parcelas divididas factorial incompleto (DPDFI) en donde se toman dos factores para evaluar el A con 2 y el B con 4 para un total de 8 niveles.
knitr::include_graphics("C:/Users/felip/Downloads/maiz.png") 
b.La razón de colocar cada factor en la parcela principal o en la subparcela
La principal razón de bloqueo y división de la parcela es si se realiza o no labranza; método de subsoleo y sin subsoleo en donde la preparación del suelo para eliminar las capas endurecidas y buscar una mejora de la humedad es la variable principal a evaluar en función al rendimiento de la semilla, luego se realiza unas subparcelas con diferentes métodos de labranza, estos bloqueos y parcelamiento creemos no fueron los más acertados para buscar el rendimiento de la semilla en función del método de labranza que se utiliza ya que estas parcelas no se pusieron en iguales condiciones climáticas, edafológicas ni nutricionales, y esto claramente afecta los resultados de rendimiento según la zona en la que se haya realizado el montaje, generando resultados muy diversos y poco precisos con lo que se busca evaluar, sumado al hecho de las diferencias significativas que marcó la labranza mínima en subsoleo y sin subsoleo se debe directamente a que la labranza mínima busca que no se haga prácticas de labranza, por lo cual en las parcelas donde se realizó subsoleo los resultados no iban a ser óptimos para esta práctica
c.La revisión de supuestos para el análisis de varianza
Creemos que los supuestos que dieron basado fueron erróneos y los corrobora el mismo artículo en donde surgen varias contradicciones.
d.La tabla del análisis de varianza
knitr::include_graphics("C:/Users/felip/Downloads/tablavar1.png") 
knitr::include_graphics("C:/Users/felip/Downloads/tablavar2.png") 
knitr::include_graphics("C:/Users/felip/Downloads/tablavar7.png") 
knitr::include_graphics("C:/Users/felip/Downloads/tablavar4.png") 
#### e.El uso de muchos análisis de varianzas en lugar de uno solo multivariante Este punto es crucial para la elaboración de los resultados, en este vemos contradicciones y resultados no certeros en los datos obtenidos en campo y lo arrojado por el análisis de la varianza, primero vemos que no se realiza un análisis multivariable de las condiciones que afectan el rendimiento de la semilla directamente, sino que basados en medias que arrojan el rendimiento de 10 mazorcas bajo 4 modos diferentes y estos a su vez se con sus pares en sitios distintos, en donde los resultados iniciales donde me informan si la labranza es crucial puede que sean no lo más indicados ya que toma medias muy diferentes que arrojaran resultados que pueden mal interpretarse según lo explicado en clase, una vez empieza el análisis de los 4 tipos diferente en cada parcela ya obtienen unos resultados más precisos, sin embargo siguen sin tomar en cuenta todas las variables que tienen que evaluar para conocer el rendimiento de una semilla de maíz.
f. El método de comparaciones de medias después del Anova
Esta comparación no fue la más indicada según nosotros dado a los supuestos y la información suministrada en clase, para esta comparaciones se usaron las medias globales de las parcelas que, 1) me tomaban las medias de rendimiento de diversas condiciones ambientales, climáticas, edafológicas, nutricionales y entre otras, 2) tomaba una media global de rendimiento de 10 mazorcas al azar en cada subparcela y no la media de rendimiento de cada mazorca, esto afecta directamente el análisis según lo explicado en clase hasta el momento, por ello para una comparación de medias más adecuada se debió tomar la media de rendimiento de cada una de las 10 mazorcas, y esta compararla solo con su contraparte en la misma zona y misma parcela, así se tendría una diferencia de medias de subparcelas en diversas condiciones.
g. La interacción de factores
El montaje experimental evalúa la interacción que surge entre arar el suelo y no ararlo y una vez esto volver a arar en diversos métodos, si bien el artículo si se centra en los resultados que se obtiene cuando estas labranzas se usan en función del rendimiento de una mazorca se puede decir que no son los resultados no son concluyentes ya que excluyen y toman diversas condiciones que me dan el rendimiento, además se generan varias contradicciones que no dejan claro si existe una interacción directa o no entre los factores.
h. La presencia de bloques
Se ve la presencia de un bloqueo principal y es la realización o no de labranza por ello aparecen 2 bloques que se desarrollan a lo largo del trabajo, estos van a estar divididos por los 4 tipos de labranza que se realizó en las parcelas labradas o no.
i. El balanceo o desbalanceo
En este se evidencia un montaje balanceado con los factores que ellos evalúan en función del rendimiento del maíz, en nuestra opinión de un modo correcto
k. Software utilizado y librería específica (en caso de ser R)
En el artículo no fue especificado el software con el que fue procesado los datos.