Autores:
Jeffer Guerrero Velasco
Laura Sofia Hernandez Sanchez
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 813. Con esta información complete la tabla del ANOVA que se muestra a continuación
library(readxl)
ANOVA <- read_excel("C:/Users/Sofia Hernandez/Documents/2020-2/Disennio de experimentos/parcial 1.xlsx")
head(ANOVA)
## # A tibble: 3 x 5
## ...1 `Sum of squares` df `Mean square (Var)` F
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 Between 6000 8 750 0.914
## 2 Within 51723 63 821 NA
## 3 Total 57723 71 NA NA
Por medio de la tabla de ANOVA se realizo el F calculado, siendo de 0.9135, lo cual indica que la variabilidad de las repeticiones es mayor que la correspondiente a los tratamientos, sin embargo, es necesario observar el comportamiento de la clorofila en las plantas ya que al ser cercano a 1 el investigador puede concluir de acuerdo a sus conocimientos bajo toda la unidad experimental. Asimismo, el F tabulado suministrado de 2.8 muestra que se acepta la hipotesis nula ya que el F calculado es menor que el F tabulado, lo que se confirma mediante el p-valor de 51.13% (distancia entre el Fcal y el Ftab) lo que apoya la no diferencia entre los tratamientos.
Fcal=0.91352
Ftab=2.8
pvalor=pf(q = 0.9135,df1 = 8,df2 = 63,lower.tail = F); pvalor # P-value
## [1] 0.5113565
d<-df(seq(0,4,0.1),8,72)
plot(seq(0,4,0.1),d,type="l",main="Curva F")
abline(v=Ftab)
abline(v=Fcal)
text(1,0.2, "No rechazo Ho")
text(3.2,0.2, "Rechazo Ho")
text(3.5,0.02, "5%")
text(0.913,0.45, "Fc=0.913")
text(1.5,0.8, "pvalor=51%")
text(2.8,0.4,"Ftab=2.8")
Con base en lo anterior, no vale la pena usar la prueba de Tukey ya que no hay una diferencia considerable entre las medias de los tratamientos.
Para aplicar la prueba de ANOVA se deben cumplir con tres supuestos escenciales.
-Las poblaciones de donde proceden las muestras son normales.
-Las muestras en las que se aplican los tratamientos deben ser independientes o elegidas al azar.
-Las poblaciones tienen todas igual varianza (homoscedasticidad).
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.
Respuesta: Pérdidas en peso (en kg) después de la limpieza.
Factor: Metodo de limpieza (M2 y M3 son mecánicos; M4 y M5 combinación mecánica y química.)
Bloque: Granjero
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, las mediciones de la granja 1 procesada por el limpiador M1 se perdieron.
library(readxl)
Algodon <- read_excel("C:/Users/Sofia Hernandez/Documents/2020-2/Disennio de experimentos/parcial 1.xlsx", sheet = 2)
head(Algodon)
## # A tibble: 5 x 7
## M G1 G2 G3 G4 G5 G6
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 M1 NA 6.75 13.0 10.3 8.01 8.42
## 2 M2 5.54 3.53 11.2 7.21 3.24 6.45
## 3 M3 7.67 4.15 9.79 8.27 6.75 5.5
## 4 M4 7.89 1.97 8.97 6.12 4.22 7.84
## 5 M5 9.27 7.39 13.4 9.13 9.2 7.13
M <- factor(c(rep("M1",6), rep("M2",6),rep("M3",6), rep("M4",6), rep("M5",6)))
G <- c(rep("G1",1),rep("G2",1),rep("G3",1),rep("G4",1),rep("G5",1),rep("G6",1))
pe_peso<-c(NA, 6.75,13.05,10.26,8.01,8.42,5.54,3.53,11.20,7.21,3.24,6.45, 7.67,4.15,9.79,8.27,6.75,5.50,7.89,1.97,8.97,6.12,4.22,7.84,9.27,7.39,13.44,9.13,9.20,7.13)
dfa<-data.frame(M,G,pe_peso)
View(dfa)
#res<-list(Algodon[,2:6])
#res_block_fact
anov1<-aov(lm(dfa$pe_peso ~ dfa$G +dfa$M))
summary(anov1)
## Df Sum Sq Mean Sq F value Pr(>F)
## dfa$G 5 120.18 24.037 17.49 1.58e-06 ***
## dfa$M 4 57.83 14.457 10.52 0.000115 ***
## Residuals 19 26.11 1.374
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 1 observation deleted due to missingness
#res_fact_block
anov2<-aov(lm(dfa$pe_peso ~ dfa$M +dfa$G))
summary(anov2)
## Df Sum Sq Mean Sq F value Pr(>F)
## dfa$M 4 56.96 14.239 10.36 0.000126 ***
## dfa$G 5 121.06 24.212 17.62 1.5e-06 ***
## Residuals 19 26.11 1.374
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 1 observation deleted due to missingness
De acuerdo con las tablas ANOVA al cambiar el orden de colocación de los efectos del modelo hace que no varien los resultados, sin embargo, se coloca primero el bloque “Granjero”, ya que este modifica las variables, posteriormente se añade el factor que en este caso fue el “Metodo”.
#Evaluavion de supuestos
plot(anov1,1:3)
b
Donde “a” representa la media y se añade en G1 y M1, haciendo el caso balanceado.
a<-c(6.75,13.05,10.26,8.01,8.42)
mean(a)
## [1] 9.298
pe_peso2<-c(9.298,6.75,13.05,10.26,8.01,8.42,5.54,3.53,11.20,7.21,3.24,6.45, 7.67, 4.15,9.79,8.27,6.75,5.50,7.89,1.97,8.97,6.12,4.22,7.84,9.27,7.39,13.44,9.13,9.20,7.13)
dfa2<-data.frame(M,G,pe_peso2)
#res_block_fact
anov3<-aov(dfa2$pe_peso2 ~ dfa$G + dfa$M)
summary(anov3)
## Df Sum Sq Mean Sq F value Pr(>F)
## dfa$G 5 120.88 24.176 18.39 7.2e-07 ***
## dfa$M 4 59.98 14.994 11.41 5.5e-05 ***
## Residuals 20 26.29 1.315
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En el analisis de varianza desbalanceado y balanceado se observo que tanto el bloqueo “granjero” como el modelo presentan un p.valor significativo, por lo que ambos afectan el peso final del algodon. Además, en el caso balanceado se obtuvo valores F mayores que el desbalanceado.
Use la función de R para generar de la distribución uniforme unos datos de carbono orgánico del suelo medida a 5 cm y 10 cm de profundidad.Suponga que la medida de la capa superior osciló entre 3.0 y 3.7 y de la capa inferior osciló entre 2 y 2.3.
set.seed(2014)
prof_5cm<-sort.int(runif(50, 3, 3.7), partial = 31)
prof_10cm<-sort.int(runif(50, 2, 2.3), partial = 31)
df1<-data.frame(prof_5cm, prof_10cm);df1
## prof_5cm prof_10cm
## 1 3.200064 2.029445
## 2 3.118236 2.010974
## 3 3.026865 2.042876
## 4 3.216780 2.041207
## 5 3.384892 2.035195
## 6 3.059382 2.030054
## 7 3.372033 2.031547
## 8 3.304012 2.022191
## 9 3.069244 2.058200
## 10 3.109104 2.057141
## 11 3.258024 2.058321
## 12 3.036855 2.056546
## 13 3.389134 2.114462
## 14 3.379940 2.090976
## 15 3.164862 2.153751
## 16 3.283415 2.140008
## 17 3.346621 2.068577
## 18 3.041619 2.072300
## 19 3.265430 2.132302
## 20 3.145261 2.062642
## 21 3.309820 2.125228
## 22 3.151100 2.120795
## 23 3.360300 2.139609
## 24 3.150527 2.100963
## 25 3.420730 2.115136
## 26 3.423666 2.067597
## 27 3.411438 2.116963
## 28 3.438139 2.140012
## 29 3.438539 2.170407
## 30 3.452974 2.172817
## 31 3.455680 2.175629
## 32 3.458720 2.176831
## 33 3.467841 2.176965
## 34 3.472839 2.211328
## 35 3.473779 2.188921
## 36 3.522733 2.202171
## 37 3.496175 2.207665
## 38 3.515967 2.222072
## 39 3.676899 2.184555
## 40 3.574713 2.221320
## 41 3.695735 2.207510
## 42 3.682269 2.250930
## 43 3.601742 2.285423
## 44 3.575732 2.272785
## 45 3.646556 2.250908
## 46 3.602918 2.236188
## 47 3.688961 2.241371
## 48 3.582438 2.282666
## 49 3.638406 2.294294
## 50 3.580440 2.298066
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.Dado que se tienen 50 datos por capa de profundidad, es necesario que la ventana de observación tenga 50 bloques, de este modo cada bloque tendrá dos datos, uno por capa.
posicion<-expand.grid(latitud =seq(0,200,length.out = 10),longitud=seq(0,100,length.out = 5));posicion
## latitud longitud
## 1 0.00000 0
## 2 22.22222 0
## 3 44.44444 0
## 4 66.66667 0
## 5 88.88889 0
## 6 111.11111 0
## 7 133.33333 0
## 8 155.55556 0
## 9 177.77778 0
## 10 200.00000 0
## 11 0.00000 25
## 12 22.22222 25
## 13 44.44444 25
## 14 66.66667 25
## 15 88.88889 25
## 16 111.11111 25
## 17 133.33333 25
## 18 155.55556 25
## 19 177.77778 25
## 20 200.00000 25
## 21 0.00000 50
## 22 22.22222 50
## 23 44.44444 50
## 24 66.66667 50
## 25 88.88889 50
## 26 111.11111 50
## 27 133.33333 50
## 28 155.55556 50
## 29 177.77778 50
## 30 200.00000 50
## 31 0.00000 75
## 32 22.22222 75
## 33 44.44444 75
## 34 66.66667 75
## 35 88.88889 75
## 36 111.11111 75
## 37 133.33333 75
## 38 155.55556 75
## 39 177.77778 75
## 40 200.00000 75
## 41 0.00000 100
## 42 22.22222 100
## 43 44.44444 100
## 44 66.66667 100
## 45 88.88889 100
## 46 111.11111 100
## 47 133.33333 100
## 48 155.55556 100
## 49 177.77778 100
## 50 200.00000 100
plot(posicion, lty = NULL, col = "white")
grid(10,5,lty = 1, col = "mediumturquoise")
df2<-data.frame(longitud_m = rep(posicion$longitud), latitud_m = rep(posicion$latitud), profundidad_cm = rep(c(5,10),each = 50), co = c(prof_5cm,prof_10cm));df2
## longitud_m latitud_m profundidad_cm co
## 1 0 0.00000 5 3.200064
## 2 0 22.22222 5 3.118236
## 3 0 44.44444 5 3.026865
## 4 0 66.66667 5 3.216780
## 5 0 88.88889 5 3.384892
## 6 0 111.11111 5 3.059382
## 7 0 133.33333 5 3.372033
## 8 0 155.55556 5 3.304012
## 9 0 177.77778 5 3.069244
## 10 0 200.00000 5 3.109104
## 11 25 0.00000 5 3.258024
## 12 25 22.22222 5 3.036855
## 13 25 44.44444 5 3.389134
## 14 25 66.66667 5 3.379940
## 15 25 88.88889 5 3.164862
## 16 25 111.11111 5 3.283415
## 17 25 133.33333 5 3.346621
## 18 25 155.55556 5 3.041619
## 19 25 177.77778 5 3.265430
## 20 25 200.00000 5 3.145261
## 21 50 0.00000 5 3.309820
## 22 50 22.22222 5 3.151100
## 23 50 44.44444 5 3.360300
## 24 50 66.66667 5 3.150527
## 25 50 88.88889 5 3.420730
## 26 50 111.11111 5 3.423666
## 27 50 133.33333 5 3.411438
## 28 50 155.55556 5 3.438139
## 29 50 177.77778 5 3.438539
## 30 50 200.00000 5 3.452974
## 31 75 0.00000 5 3.455680
## 32 75 22.22222 5 3.458720
## 33 75 44.44444 5 3.467841
## 34 75 66.66667 5 3.472839
## 35 75 88.88889 5 3.473779
## 36 75 111.11111 5 3.522733
## 37 75 133.33333 5 3.496175
## 38 75 155.55556 5 3.515967
## 39 75 177.77778 5 3.676899
## 40 75 200.00000 5 3.574713
## 41 100 0.00000 5 3.695735
## 42 100 22.22222 5 3.682269
## 43 100 44.44444 5 3.601742
## 44 100 66.66667 5 3.575732
## 45 100 88.88889 5 3.646556
## 46 100 111.11111 5 3.602918
## 47 100 133.33333 5 3.688961
## 48 100 155.55556 5 3.582438
## 49 100 177.77778 5 3.638406
## 50 100 200.00000 5 3.580440
## 51 0 0.00000 10 2.029445
## 52 0 22.22222 10 2.010974
## 53 0 44.44444 10 2.042876
## 54 0 66.66667 10 2.041207
## 55 0 88.88889 10 2.035195
## 56 0 111.11111 10 2.030054
## 57 0 133.33333 10 2.031547
## 58 0 155.55556 10 2.022191
## 59 0 177.77778 10 2.058200
## 60 0 200.00000 10 2.057141
## 61 25 0.00000 10 2.058321
## 62 25 22.22222 10 2.056546
## 63 25 44.44444 10 2.114462
## 64 25 66.66667 10 2.090976
## 65 25 88.88889 10 2.153751
## 66 25 111.11111 10 2.140008
## 67 25 133.33333 10 2.068577
## 68 25 155.55556 10 2.072300
## 69 25 177.77778 10 2.132302
## 70 25 200.00000 10 2.062642
## 71 50 0.00000 10 2.125228
## 72 50 22.22222 10 2.120795
## 73 50 44.44444 10 2.139609
## 74 50 66.66667 10 2.100963
## 75 50 88.88889 10 2.115136
## 76 50 111.11111 10 2.067597
## 77 50 133.33333 10 2.116963
## 78 50 155.55556 10 2.140012
## 79 50 177.77778 10 2.170407
## 80 50 200.00000 10 2.172817
## 81 75 0.00000 10 2.175629
## 82 75 22.22222 10 2.176831
## 83 75 44.44444 10 2.176965
## 84 75 66.66667 10 2.211328
## 85 75 88.88889 10 2.188921
## 86 75 111.11111 10 2.202171
## 87 75 133.33333 10 2.207665
## 88 75 155.55556 10 2.222072
## 89 75 177.77778 10 2.184555
## 90 75 200.00000 10 2.221320
## 91 100 0.00000 10 2.207510
## 92 100 22.22222 10 2.250930
## 93 100 44.44444 10 2.285423
## 94 100 66.66667 10 2.272785
## 95 100 88.88889 10 2.250908
## 96 100 111.11111 10 2.236188
## 97 100 133.33333 10 2.241371
## 98 100 155.55556 10 2.282666
## 99 100 177.77778 10 2.294294
## 100 100 200.00000 10 2.298066
Es necesario primero hacer el Data Frame, se replica por dos los datos de la longitud y la latitud, ya que cada uno de estos debe tener dos datos, uno de la capa de 5 cm y otro de la capa de 10 cm.
library(plotly)
## Loading required package: ggplot2
## Warning: package 'ggplot2' was built under R version 4.0.3
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
library(ggplot2)
graf3d<-plot_ly(df2,x = df2$longitud_m,y = df2$latitud_m,z = df2$profundidad_cm, color = df2$co, colors = c("salmon","salmon1","salmon2","salmon3","salmon4"), width = 8, type = "scatter3d") %>% layout( scene = list(xaxis = list(title = 'Longitud'), yaxis = list(title = 'Latitud'), zaxis = list(title = 'Profundidad'), textfont.color = 'carbono'));graf3d
## No scatter3d mode specifed:
## Setting the mode to markers
## Read more about this attribute -> https://plot.ly/r/reference/#scatter-mode
## Warning: `arrange_()` is deprecated as of dplyr 0.7.0.
## Please use `arrange()` instead.
## See vignette('programming') for more help
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.
Compare si se encuentran diferencias en la media de carbono entre capas utilizando un nivel de confianza del 95%.Dado que son dos factores independientes, se puede realizar la prueba T pareada.
\[H_o:\mu_{prof_5cm}=\mu_{prof10cm}\\H_a:\mu_{prof_5cm}\neq\mu_{prof10cm}\]
t.test(prof_5cm,prof_10cm,"two.sided",paired = T)
##
## Paired t-test
##
## data: prof_5cm and prof_10cm
## t = 64.386, df = 49, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 1.200788 1.278160
## sample estimates:
## mean of the differences
## 1.239474
Una vez realizada la prueba T pareada,con un nivel de confianza de 95%, tenemos la información de que las medias de ambas capas son estadísticamente distintas, siendo mayor la de 5 cm, resultados que no contrastan con diferentes reportes en la literatura, que informan esa tendencia general en los suelo.
El diseño más simple del sistema 3^k es el diseño 3^2, el cual tiene dos factores, cada uno con tres niveles, obteniendo un total de 9 tratamientos diferentes.
El modelo estadístico para el diseño 3^2 se puede escribir considerando el efecto individual de cada factor y de la interacción entre ambos, como se presenta:
\[y_{ijk}=\mu+\alpha_i+\tau_j+(\alpha\tau)_{ij}+\epsilon_{ijk}\\i,j=0,1,2\\k=1...,r.\] \(y_{ijk}\) :Variable respuesta \(\alpha_i\) : Es el fecto del factor \(A\) \(\tau_j\) : Es el efecto del factor \(B\) \((\alpha\tau)_{ij}\) : representa el efecto de la interacción entre los dos factores.
En contraposición la Hipotesis
\[H_0:(\alpha\beta)_{ij}=0\] Muestra que no hay efecto de interacción de los factores A y B sobre la variable respuesta. Si esto no se rechaza, se prosigue a las hipotesis \(H_0: \alpha_i=0\) No hay efecto del factor A, sobre la variable respuesta \(H_0:\beta_j=0\) No hay efecto del factor B, sobre la variable respuesta
D<-expand.grid(F1=c(3.25, 3.75,4.25), F2= c(4,5,6))
D<-rbind(D,D)
set.seed(2020)
D<-D[order(sample(1:18)),]
class(D)
## [1] "data.frame"
D$biomasa=sort.int(rnorm(18,3,0.3), partial=9)
names(D)
## [1] "F1" "F2" "biomasa"
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
D1<-rename(D, Insecticida=F1, Aplicacion=F2, biomasa=biomasa); D1
## Insecticida Aplicacion 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
Se realiza un diagrama de arbol teniendo encuenta:
-Los efectos principales de F1: dosis de un insecticida que se sospecha tiene un efecto de disminución del crecimiento (biomasa).
-F2: número de aplicaciones durante el desarrollo del cultivo.
-Respuesta: Crecimiento (biomasa).
library(collapsibleTree)
collapsibleTreeSummary(D1,(c("Insecticida", "Aplicacion")),nodeSize = "leafCount", maxPercent = 50,attribute = "biomasa",linkLength =150, fillFun = colorspace::terrain_hcl,tooltip = T)
FC_CA<-aov(biomasa~Insecticida*Aplicacion, D1)
summary(FC_CA)
## Df Sum Sq Mean Sq F value Pr(>F)
## Insecticida 1 0.3803 0.3803 4.548 0.0511 .
## Aplicacion 1 0.3160 0.3160 3.780 0.0722 .
## Insecticida:Aplicacion 1 0.0013 0.0013 0.015 0.9042
## Residuals 14 1.1705 0.0836
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Con el ANOVA realizado se concluye que no existen evidencias estadisticas para decir que la dosis de insecticida y sus aplicaciones, junto con su interacción afecta la variable respuesta de crecimiento de las plantas (biomasa). Así pues, al no tener un p.valor<0.05, no se requiere una prueba de comparación de medias.
#Supuestos del ANOVA
shapiro.test(FC_CA$residuals) # Cumple con norvalidad
##
## Shapiro-Wilk normality test
##
## data: FC_CA$residuals
## W = 0.96739, p-value = 0.7473
#leveneTest(FC_CA$residuals~Insecticida*Aplicacion)
library(sjPlot)
## Warning: package 'sjPlot' was built under R version 4.0.3
library(sjmisc)
## Warning: package 'sjmisc' was built under R version 4.0.3
library(ggplot2)
theme_set(theme_sjplot())
D1$Aplicacion <- to_factor(D1$Aplicacion)
# fit model with interaction
fit <- lm(biomasa ~ Insecticida + Aplicacion, data = D1)
plot_model(fit, type = "pred", terms = c("Insecticida", "Aplicacion"))
library(lattice)
bwplot(D1$biomasa~D1$Insecticida|D1$Aplicacion)
tabla=tapply(D1$biomasa, list(D1$Insecticida,D1$Aplicacion), mean)
tabla
## 4 5 6
## 3.25 2.965294 2.634592 2.592757
## 3.75 3.080186 3.128949 2.801410
## 4.25 3.344110 2.894908 3.021714
tabla1<-tapply(D1$biomasa,D1$Insecticida, mean);tabla1
## 3.25 3.75 4.25
## 2.730881 3.003515 3.086911
tabla2<-tapply(D1$biomasa,D1$Aplicacion, mean);tabla2
## 4 5 6
## 3.129864 2.886150 2.805294
Por medio de los graficos y del calculo de las medias en ambos factores, se observa que no existe interacción entre las dosis del insecticida y sus aplicaciones.
Adicionalmente, se busca minimizar el efecto de las variables no controlables denominados “Covariables”, las cuales estabilizan y minimizan la variabilidad de las respuestas, identifIcando además los factores que contribuyen a las mayores causas de variabilidad. Al incluir una covariable en el modelo se puede reducir el error en el modelo para incrementar la potencia de las pruebas de los factores.
set.seed(123)
D1$arcilla<-sort(runif(18,0.2,0.4),decreasing = TRUE)
D1
## Insecticida Aplicacion biomasa arcilla
## 2 3.75 4 2.708826 0.3913667
## 10 3.25 4 2.772692 0.3880935
## 16 3.25 6 2.143359 0.3799650
## 4 3.25 5 2.560519 0.3784838
## 13 3.25 5 2.708666 0.3766035
## 6 4.25 5 2.773705 0.3576610
## 17 3.75 6 2.770350 0.3355141
## 8 3.75 6 2.832470 0.3145267
## 14 3.75 5 2.898280 0.3102870
## 5 3.75 5 3.359619 0.3056211
## 9 4.25 6 3.054099 0.2913229
## 1 3.25 4 3.157896 0.2906668
## 12 4.25 4 3.487669 0.2817954
## 11 3.75 4 3.451547 0.2575155
## 15 4.25 5 3.016111 0.2492175
## 7 3.25 6 3.042156 0.2205849
## 3 4.25 4 3.200552 0.2091113
## 18 4.25 6 2.989329 0.2084119
\[ Y_{ijk}= \mu + I_{i} + A_{j} +\beta(X_{ijk} - \overline{X_{...}}) + \epsilon_{ijk} \\ i=1..3 \\ j=1..3 \\ k=1..2\]
ancova1<-aov(biomasa~arcilla+(Insecticida*Aplicacion), D1)
summary(ancova1)
## Df Sum Sq Mean Sq F value Pr(>F)
## arcilla 1 0.8518 0.8518 17.103 0.00166 **
## Insecticida 1 0.0354 0.0354 0.711 0.41713
## Aplicacion 2 0.3970 0.1985 3.986 0.04989 *
## Insecticida:Aplicacion 2 0.0361 0.0181 0.362 0.70398
## Residuals 11 0.5478 0.0498
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Mediante la covarianza se hallo que la arcilla presenta evidencia estadistica para influir en la biomasa final de la planta al tener un p.valor menor de 0.05, esto puede ocurrir. Así tambien, esta ayuda a la dismincion de error; se observa en los siguientes graficos.
fit <- lm(biomasa ~ Insecticida * Aplicacion * arcilla, data = D1)
plot_model(fit, type = "pred", terms = c("arcilla"))
library(ggplot2)
D1 %>%
ggplot(aes(x = arcilla, y = biomasa)) +
geom_point()+geom_abline()
Linealidad entre la covariable y la variable de resultado en cada nivel de la variable de agrupación (diagrama de dispersión agrupado de la covariable y la variable de resultado)
Homogeneidad de las pendientes de regresión . Las pendientes de las líneas de regresión, formadas por la covariable y la variable de resultado, deben ser las mismas para cada grupo. Esta suposición evalúa que no hay interacción entre el resultado y la covariable. Las líneas de regresión trazadas por grupos deben ser paralelas.
La variable de resultado debe tener una distribución aproximadamente normal . Esto se puede verificar usando la prueba de normalidad de Shapiro-Wilk en los residuos del modelo.
Homoscedasticidad u homogeneidad de la varianza de residuos para todos los grupos. Se supone que los residuos tienen una varianza constante (homocedasticidad).
Tomado de: DATANOVIA
library(tidyverse)
library(ggpubr)
library(rstatix)
library(broom)
# Supuesto de linealidad
ggscatter( D1, x = "arcilla", y = "biomasa", facet.by = c("Insecticida", "Aplicacion"), short.panel.labs = FALSE
)+ stat_smooth(method = "loess", span = 0.9)
# Independencia
Mod_b2<-aov(D1$arcilla~ D1$Insecticida*D1$Aplicacion)
summary(Mod_b2)
## Df Sum Sq Mean Sq F value Pr(>F)
## D1$Insecticida 1 0.01591 0.015905 4.270 0.0611 .
## D1$Aplicacion 2 0.00455 0.002273 0.610 0.5593
## D1$Insecticida:D1$Aplicacion 2 0.00095 0.000475 0.127 0.8815
## Residuals 12 0.04470 0.003725
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
No es significativo 0.969 porque es mayor a 0.05, por lo que cumple el supuesto de independencia.
# Homogeneidad de las pendientes de regresión
D1 %>%
anova_test(biomasa ~ arcilla + Insecticida + Aplicacion +
Insecticida*Aplicacion + arcilla*Insecticida +
arcilla*Aplicacion + arcilla*Aplicacion*Insecticida
)
## Coefficient covariances computed by hccm()
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 arcilla 1 6 15.032 0.008 * 0.715
## 2 Insecticida 1 6 0.650 0.451 0.098
## 3 Aplicacion 2 6 5.726 0.041 * 0.656
## 4 Insecticida:Aplicacion 2 6 2.586 0.155 0.463
## 5 arcilla:Insecticida 1 6 7.697 0.032 * 0.562
## 6 arcilla:Aplicacion 2 6 0.785 0.498 0.207
## 7 arcilla:Insecticida:Aplicacion 2 6 0.189 0.832 0.059
Hubo homogeneidad de las pendientes de regresión ya que los términos de interacción, entre la covariable arcilla y las variables de agrupamiento dosis de insecticida y numero de aplicacion no fue estadísticamente significativa, p> 0.05.
# Normalidad de los residuos
model <- lm(biomasa ~ arcilla + Insecticida*Aplicacion, data = D1)
# Inspect the model diagnostic metrics
model.metrics <- augment(model)
head(model.metrics, 3)
## # A tibble: 3 x 11
## .rownames biomasa arcilla Insecticida Aplicacion .fitted .resid .std.resid
## <chr> <dbl> <dbl> <dbl> <fct> <dbl> <dbl> <dbl>
## 1 2 2.71 0.391 3.75 4 2.81 -0.102 -0.565
## 2 10 2.77 0.388 3.25 4 2.80 -0.0305 -0.184
## 3 16 2.14 0.380 3.25 6 2.36 -0.220 -1.40
## # ... with 3 more variables: .hat <dbl>, .sigma <dbl>, .cooksd <dbl>
shapiro.test(model$residuals)
##
## Shapiro-Wilk normality test
##
## data: model$residuals
## W = 0.97409, p-value = 0.8701
La prueba de Shapiro no fue significativa (p> 0.05), entonces, se asume normalidad de los residuos.
# Homogeneidad de variaciones
library(car)
## Warning: package 'car' was built under R version 4.0.3
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:purrr':
##
## some
## The following object is masked from 'package:dplyr':
##
## recode
Se cumple la homogeneidad de las pendientes de regresión arcilla:Insecticida:Aplicacion es no significativo, por lo que se cumple las pendientes de regresion.
Se concluye que el ANCOVA es un modelo de análisis de regresión lineal que implica que la relación de variables (dependientes e independientes) tiene que ser lineal. En contraposicion el propósito de ANOVA es verificar si los datos de varios grupos tienen una media común o no. (mort-sure, 2020).
En un diseño completamente anidado, cada nivel del factor superior conduce a dos (o más) niveles de cada factor o etapa sucesiva. En un diseño anidado escalonado, por el contrario, solo uno de los dos niveles del factor sucesivo conduce a la siguiente etapa de dos niveles (Lawson, 2015).
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.
La tabla en la que se recogen los datos se muestra acontinuación:
polymer <- read_excel("C:/Users/Sofia Hernandez/Documents/2020-2/Disennio de experimentos/parcial 1.xlsx", sheet = 5)
View(polymer)
El ANOVA del diseño es:
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.design
##
## Attaching package: 'daewr'
## The following object is masked _by_ '.GlobalEnv':
##
## polymer
mod2<-lm(respuesta ~ ue + ue:finca + ue:finca:variedad, data = polymer)
summary(mod2)
##
## Call:
## lm(formula = respuesta ~ ue + ue:finca + ue:finca:variedad, data = polymer)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.7708 -2.0801 -0.5685 1.6336 8.5275
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.62978 0.68482 12.601 <2e-16 ***
## ue -0.12338 0.10769 -1.146 0.256
## ue:fincaB 0.10718 0.11007 0.974 0.333
## ue:fincaA:variedad 0.05203 0.06741 0.772 0.443
## ue:fincaB:variedad NA NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.948 on 76 degrees of freedom
## Multiple R-squared: 0.02008, Adjusted R-squared: -0.0186
## F-statistic: 0.5191 on 3 and 76 DF, p-value: 0.6704
collapsibleTreeSummary(polymer,(c("finca", "variedad", "test")),nodeSize = "leafCount", maxPercent = 50,attribute = "respuesta",linkLength =150, fillFun = colorspace::terrain_hcl,tooltip = T)
c = (100*(7.0127)/(7.0127+1.2305+0.8795))
c
## [1] 76.87088
Por medio de la siguiente formula se obtuvo que la mayor variación esta dada por las parcelas.
\[(100 × (7.2427) / (7.2427 + 1.0296 + 0.6568)) = 76\% \]
La variación total se debe a la variabilidad entre parcelas ya que estas causan mayor impresición, por otro lado, la variabilidad dentro de la parcela o de finca a finca es minima.
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).
Se sacan los datos de la libreria
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':
##
## bootstrap
data(K)
media<-tapply(K$K,K$lab,mean)
desv<-tapply(K$K,K$lab,sd)
varia<-tapply(K$K,K$lab,var)
media_global<-mean(K$K)
df6<-data.frame(Media = (round(media, 2)), Desv_esta = (round(desv, 2)), varianza = (round(varia, 2)))
library(DT)
tabla8<-datatable(df6,class = 'cell-border stripe');tabla8
En la tabla es posible observar la media, la desviación estándar y la varianza de los datos obtenidos por los laboratorios, de entrada la media más alta es del laboratorio J, la más baja la del laboratorio H, los resultados menos dispersos son los del laboratorio G, y los más dispersos los del laboratorio E.
graf_1<-barplot(df6$Media, axes = F, axisnames = F, ylim = c(1,450), col = c('green3', 'darkcyan','tomato2','navajowhite1','khaki','thistle3', 'tan1','mediumspringgreen'), main = "Cantidad de Potasio", xlab = "Laboratorio", ylab = "Media")
axis(1,labels = c("B","D","E","F","G","H","I","J"), at=graf_1)
axis(2,at =seq(0,450,by = 25))
segments(graf_1-0.1,df6$Media-df6$Desv_esta,graf_1+0.1,df6$Media-df6$Desv_esta,lwd = 2)
segments(graf_1-0.1,df6$Media+df6$Desv_esta,graf_1+0.1,df6$Media+df6$Desv_esta,lwd = 2)
segments(graf_1,df6$Media-df6$Desv_esta,graf_1,df6$Media+df6$Desv_esta,lwd = 2)
segments(0,media_global,10,media_global, lty = 10, lwd = 3)
text(6.7,340,"Media Global")
A Continuación se presenta un gráfico de barras, en la cual se grafica el promedio de las resultados por laboratorio, además se añade también la media global, esto para poder tener una idea de que tan cerca se encuentra cada laboratorio, se ve una fuerte diferencia en los resultados del laboratorio H y la media global, lo cual ratifica lo comentado anteriormente.
violin<- K %>%plot_ly(x = K$lab,y = K$K,split = K$lab,type = 'violin',box = list(visible = T),meanline = list(visible = T))
violin <- violin %>%layout(xaxis = list(title = "Laboratorio"),yaxis = list(title = "K_potasio",zeroline = F));violin
En este diagrama de violines se pueden observar datos de suma importancia, además de las medias de cada laboratorio, se puede observar la dispersión de los datos, podemos validar gráficamente que el laboratorio con los datos menos dispersos es el G, y el más disperso el laboratorio E.
Una vez analizados estos gráficos y tablas, se puede empezar a pensar en ciertas cosas, la primera es la clara diferencia que hay entre los datos del laboratorio H y el resto, y que dada la poca dispersión de los datos y lo cercano que está a la media global el laboratorio podria ser el mas confiable de los 8.
Para validar si hay alguna diferencia estadística es importante validar estadísticamente esto, lo primero será validar si los datos cumplen con una distribución normal, para lo cual se realizará la prueba Shapiro
\[H_o: Datos \ con \ distribucion \ normal\\H_a: Datos\ con\ distribucion\ no\ normal \]
shapiro.test(K$K)
##
## Shapiro-Wilk normality test
##
## data: K$K
## W = 0.9616, p-value = 0.02723
Dado el resultado de la prueba es posible deducir que la distribución de los datos no es normal, dado que el p valor de la prueba fue de 0.027.
Dado esto tenemos unos datos que no son paramétricos, por lo cual no podemos realizar la prueba ANOVA para validar si los resultados son estadísticamente iguales o si estos difieren, por este motivo la prueba que se va usar en este caso es la prueba Kruskal-Wallis, la cual es precisamente para conocer si hay diferencia entre los datos, ya que estos no son paramétricos, además de esto se usa dado que los grupos a comparar son 8, por lo cual otras pruebas quedan descartadas.
\[H_o:\mu_{labB}=\mu_{labD}=\mu_{labE}\cdots=\mu_{labJ}\\H_a:\mu_{labB}\neq\mu_{labD}\neq\mu_{labE}\cdots\neq\mu_{labJ}\]
kruskal.test(K$K,K$lab)
##
## Kruskal-Wallis rank sum test
##
## data: K$K and K$lab
## Kruskal-Wallis chi-squared = 24.482, df = 7, p-value = 0.000937
Los resultados de la prueba indican que al menos hay diferencia en unos de los resultados de las pruebas de los laboratorios, una vez conocido esto lo que hay que saber es cual de los 8 presenta esta diferencia, para esto se realiza la prueba nemenyi, que es un aprueba diseñada para que después de saber que existe una diferencia entre datos no paramétricos, se pueda calcular cual muestra o muestras presentan esta diferencia
library(PMCMR)
## PMCMR is superseded by PMCMRplus and will be no longer maintained. You may wish to install PMCMRplus instead.
posthoc.kruskal.nemenyi.test(K$K,K$lab, method = "none")
## Warning in posthoc.kruskal.nemenyi.test.default(K$K, K$lab, method = "none"):
## 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 and 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: none
Los resultados obtenidos por esta prueba indican que el dato que es diferente es del laboratorio H, lo cual confirma la idea que se tenía después de ver las gráficas anteriores.
Para finalizar dada la información recolectada en el ejercicio, el laboratorio el cual tiene las pruebas más confiables es el G, si bien cierto que en las medias los laboratorios, excluyendo al H, no varían mucho, este laboratorio es el más preciso, la dispersión de sus datos es muy pequeña, item que le da un valor agregado, convirtiéndolo en el laboratorio más confiable de los 8.
Se quiere estudiar el rendimiento de 4 variedades distintas de lechuga batavia, a tres planes de fertilización distinta, sin embargo también se quiere conocer en qué condiciones geomorfológicas responde mejor, es por eso que se siembran en dos partes, una en llanura, la zona 1, y la otra en la falda de una montaña, la zona de 2, además de esto, dadas las características de la especie, se mide la porosidad del suelo a la hora de la siembra de 240 plantas, 120 por cada zona, 60 por cada variedad, 80 por plan de fertilización.
set.seed(2006)
library(DT)
rto_lech<-sort.int(rnorm(240,700,120), partial = 20) # Datos de rendimiento
bloque<-gl(2,120,240,paste0('zona',1:2)) # Las Zonas que seran los bloques
vari<-gl(4,30,240,paste0('var',1:4)) # El factor de la parcela principal, las variedades
Fert<-gl(3,10,240,paste0('F',1:3)) # El factor de la subparcela
poros<-sort.int(rnorm(240,36,3), decreasing = F ) # la covariable de porosidad
df_pd<-data.frame(Rendimiento = rto_lech, Zona = bloque, Variedad = vari,Fertilizacion = Fert,Porosidad = poros);df_pd
## Rendimiento Zona Variedad Fertilizacion Porosidad
## 1 359.1996 zona1 var1 F1 27.12308
## 2 397.0892 zona1 var1 F1 29.57477
## 3 388.6930 zona1 var1 F1 30.32433
## 4 471.2730 zona1 var1 F1 30.57925
## 5 458.3187 zona1 var1 F1 30.69360
## 6 436.2487 zona1 var1 F1 30.86479
## 7 458.5413 zona1 var1 F1 31.16034
## 8 416.2528 zona1 var1 F1 31.28301
## 9 462.4292 zona1 var1 F1 31.34434
## 10 472.3965 zona1 var1 F1 31.55054
## 11 466.0759 zona1 var1 F2 31.62091
## 12 448.4592 zona1 var1 F2 31.70162
## 13 495.2088 zona1 var1 F2 31.84480
## 14 494.1111 zona1 var1 F2 31.94492
## 15 488.0672 zona1 var1 F2 31.98275
## 16 485.3443 zona1 var1 F2 32.08756
## 17 499.8269 zona1 var1 F2 32.10607
## 18 498.9484 zona1 var1 F2 32.19697
## 19 481.8313 zona1 var1 F2 32.20416
## 20 505.4609 zona1 var1 F2 32.30211
## 21 558.6078 zona1 var1 F3 32.37136
## 22 508.3664 zona1 var1 F3 32.46108
## 23 560.9864 zona1 var1 F3 32.54396
## 24 541.1602 zona1 var1 F3 32.54677
## 25 560.6774 zona1 var1 F3 32.54772
## 26 543.3227 zona1 var1 F3 32.54825
## 27 551.3675 zona1 var1 F3 32.65368
## 28 521.2263 zona1 var1 F3 32.69139
## 29 546.7770 zona1 var1 F3 32.72613
## 30 528.7961 zona1 var1 F3 32.74982
## 31 551.6053 zona1 var2 F1 32.78356
## 32 527.8929 zona1 var2 F1 32.84244
## 33 559.0164 zona1 var2 F1 32.88219
## 34 569.2245 zona1 var2 F1 32.99930
## 35 566.5412 zona1 var2 F1 33.21488
## 36 572.3302 zona1 var2 F1 33.23983
## 37 580.1103 zona1 var2 F1 33.26809
## 38 578.2475 zona1 var2 F1 33.31017
## 39 587.6821 zona1 var2 F1 33.34225
## 40 594.5340 zona1 var2 F1 33.38624
## 41 581.9141 zona1 var2 F2 33.41786
## 42 574.5588 zona1 var2 F2 33.45161
## 43 573.5604 zona1 var2 F2 33.49585
## 44 578.0379 zona1 var2 F2 33.56603
## 45 585.0095 zona1 var2 F2 33.62498
## 46 591.0305 zona1 var2 F2 33.62934
## 47 592.3624 zona1 var2 F2 33.64558
## 48 591.3959 zona1 var2 F2 33.64826
## 49 787.8322 zona1 var2 F2 33.65967
## 50 679.8890 zona1 var2 F2 33.68622
## 51 688.8563 zona1 var2 F3 33.74554
## 52 712.8978 zona1 var2 F3 33.76077
## 53 856.8427 zona1 var2 F3 33.79119
## 54 854.1386 zona1 var2 F3 33.88183
## 55 748.4855 zona1 var2 F3 33.89975
## 56 677.9226 zona1 var2 F3 33.92400
## 57 719.1168 zona1 var2 F3 33.92532
## 58 725.0253 zona1 var2 F3 33.97454
## 59 839.2867 zona1 var2 F3 34.01369
## 60 714.8916 zona1 var2 F3 34.04606
## 61 610.0416 zona1 var3 F1 34.04885
## 62 858.1947 zona1 var3 F1 34.13014
## 63 706.0661 zona1 var3 F1 34.14125
## 64 615.6746 zona1 var3 F1 34.14244
## 65 648.0567 zona1 var3 F1 34.15064
## 66 719.6349 zona1 var3 F1 34.16978
## 67 887.3120 zona1 var3 F1 34.24114
## 68 763.8637 zona1 var3 F1 34.24630
## 69 860.1709 zona1 var3 F1 34.26926
## 70 659.1101 zona1 var3 F1 34.36629
## 71 596.3549 zona1 var3 F2 34.44162
## 72 725.1510 zona1 var3 F2 34.47971
## 73 845.4188 zona1 var3 F2 34.53295
## 74 683.9488 zona1 var3 F2 34.58593
## 75 774.6607 zona1 var3 F2 34.60176
## 76 693.6911 zona1 var3 F2 34.61117
## 77 855.3419 zona1 var3 F2 34.65168
## 78 646.0588 zona1 var3 F2 34.66557
## 79 614.7456 zona1 var3 F2 34.68080
## 80 659.6315 zona1 var3 F2 34.71739
## 81 675.5881 zona1 var3 F3 34.72408
## 82 851.4293 zona1 var3 F3 34.88023
## 83 744.4961 zona1 var3 F3 34.90825
## 84 621.3244 zona1 var3 F3 34.91209
## 85 737.2088 zona1 var3 F3 34.95509
## 86 626.3964 zona1 var3 F3 34.95848
## 87 699.9759 zona1 var3 F3 35.00977
## 88 731.5481 zona1 var3 F3 35.03801
## 89 698.2007 zona1 var3 F3 35.08569
## 90 787.2269 zona1 var3 F3 35.11630
## 91 606.6997 zona1 var4 F1 35.12464
## 92 699.9355 zona1 var4 F1 35.13048
## 93 643.9997 zona1 var4 F1 35.14689
## 94 725.4493 zona1 var4 F1 35.15961
## 95 759.4513 zona1 var4 F1 35.17581
## 96 691.1683 zona1 var4 F1 35.27012
## 97 887.4207 zona1 var4 F1 35.38039
## 98 753.2683 zona1 var4 F1 35.41621
## 99 635.8285 zona1 var4 F1 35.43147
## 100 862.5129 zona1 var4 F1 35.46302
## 101 652.2302 zona1 var4 F2 35.49271
## 102 753.2454 zona1 var4 F2 35.49326
## 103 818.2147 zona1 var4 F2 35.55754
## 104 669.9743 zona1 var4 F2 35.57672
## 105 809.1687 zona1 var4 F2 35.57828
## 106 810.6600 zona1 var4 F2 35.58561
## 107 646.0039 zona1 var4 F2 35.60021
## 108 769.1871 zona1 var4 F2 35.60025
## 109 699.1226 zona1 var4 F2 35.60328
## 110 690.7941 zona1 var4 F2 35.60929
## 111 726.8520 zona1 var4 F3 35.61384
## 112 667.6053 zona1 var4 F3 35.62880
## 113 813.3228 zona1 var4 F3 35.63083
## 114 675.2839 zona1 var4 F3 35.68140
## 115 680.0754 zona1 var4 F3 35.70383
## 116 840.3919 zona1 var4 F3 35.74308
## 117 597.8687 zona1 var4 F3 35.79086
## 118 735.6058 zona1 var4 F3 35.80248
## 119 608.0425 zona1 var4 F3 35.92316
## 120 846.9249 zona1 var4 F3 35.96110
## 121 775.8710 zona2 var1 F1 35.97269
## 122 680.7381 zona2 var1 F1 35.97928
## 123 790.2587 zona2 var1 F1 36.01659
## 124 625.4226 zona2 var1 F1 36.03259
## 125 609.9509 zona2 var1 F1 36.03275
## 126 812.8573 zona2 var1 F1 36.09098
## 127 639.3629 zona2 var1 F1 36.13890
## 128 662.5485 zona2 var1 F1 36.17188
## 129 812.2802 zona2 var1 F1 36.17506
## 130 781.7712 zona2 var1 F1 36.20469
## 131 721.7405 zona2 var1 F2 36.21864
## 132 756.5226 zona2 var1 F2 36.22679
## 133 741.5489 zona2 var1 F2 36.23755
## 134 608.5795 zona2 var1 F2 36.25633
## 135 756.2224 zona2 var1 F2 36.29088
## 136 774.7485 zona2 var1 F2 36.30176
## 137 660.8723 zona2 var1 F2 36.36144
## 138 678.6423 zona2 var1 F2 36.42499
## 139 887.7329 zona2 var1 F2 36.42730
## 140 678.9000 zona2 var1 F2 36.43296
## 141 730.7997 zona2 var1 F3 36.48830
## 142 777.1893 zona2 var1 F3 36.50799
## 143 715.3083 zona2 var1 F3 36.56631
## 144 743.0494 zona2 var1 F3 36.59425
## 145 745.2280 zona2 var1 F3 36.62287
## 146 785.5601 zona2 var1 F3 36.63109
## 147 750.5512 zona2 var1 F3 36.64806
## 148 609.1007 zona2 var1 F3 36.67678
## 149 828.2083 zona2 var1 F3 36.86516
## 150 844.3507 zona2 var1 F3 36.86851
## 151 694.6310 zona2 var2 F1 36.91171
## 152 639.0627 zona2 var2 F1 36.94131
## 153 853.7050 zona2 var2 F1 37.00603
## 154 595.0875 zona2 var2 F1 37.01737
## 155 866.3271 zona2 var2 F1 37.03294
## 156 690.1224 zona2 var2 F1 37.04803
## 157 753.1908 zona2 var2 F1 37.05416
## 158 667.7656 zona2 var2 F1 37.06913
## 159 697.7283 zona2 var2 F1 37.12342
## 160 820.8346 zona2 var2 F1 37.18262
## 161 676.6315 zona2 var2 F2 37.22176
## 162 664.5169 zona2 var2 F2 37.25704
## 163 725.6989 zona2 var2 F2 37.36139
## 164 840.2466 zona2 var2 F2 37.40699
## 165 752.3770 zona2 var2 F2 37.40757
## 166 741.1694 zona2 var2 F2 37.41162
## 167 712.8859 zona2 var2 F2 37.42639
## 168 796.7892 zona2 var2 F2 37.43759
## 169 741.3119 zona2 var2 F2 37.43930
## 170 611.8249 zona2 var2 F2 37.44231
## 171 732.9357 zona2 var2 F3 37.46315
## 172 635.3194 zona2 var2 F3 37.48588
## 173 731.9325 zona2 var2 F3 37.51942
## 174 597.0427 zona2 var2 F3 37.68780
## 175 597.2748 zona2 var2 F3 37.74677
## 176 772.6769 zona2 var2 F3 37.76209
## 177 649.4134 zona2 var2 F3 37.80629
## 178 714.7854 zona2 var2 F3 37.83898
## 179 679.7268 zona2 var2 F3 37.84220
## 180 821.8485 zona2 var2 F3 37.90367
## 181 620.1480 zona2 var3 F1 37.91111
## 182 642.2723 zona2 var3 F1 37.91835
## 183 828.8958 zona2 var3 F1 37.99561
## 184 646.1875 zona2 var3 F1 38.01773
## 185 863.4576 zona2 var3 F1 38.02202
## 186 760.4752 zona2 var3 F1 38.04253
## 187 801.1727 zona2 var3 F1 38.04429
## 188 637.5760 zona2 var3 F1 38.06474
## 189 604.8123 zona2 var3 F1 38.06776
## 190 621.7117 zona2 var3 F1 38.09484
## 191 631.3023 zona2 var3 F2 38.09966
## 192 851.5299 zona2 var3 F2 38.11824
## 193 685.8981 zona2 var3 F2 38.12974
## 194 816.3736 zona2 var3 F2 38.19821
## 195 718.5891 zona2 var3 F2 38.21828
## 196 706.7845 zona2 var3 F2 38.24870
## 197 640.9183 zona2 var3 F2 38.27139
## 198 885.4539 zona2 var3 F2 38.34889
## 199 727.1640 zona2 var3 F2 38.35177
## 200 774.7659 zona2 var3 F2 38.70147
## 201 671.5392 zona2 var3 F3 38.70761
## 202 764.6918 zona2 var3 F3 38.74121
## 203 708.5020 zona2 var3 F3 38.76147
## 204 687.4206 zona2 var3 F3 38.82990
## 205 598.6932 zona2 var3 F3 38.85536
## 206 654.5246 zona2 var3 F3 39.00754
## 207 601.6366 zona2 var3 F3 39.03320
## 208 755.6607 zona2 var3 F3 39.16088
## 209 834.7377 zona2 var3 F3 39.18319
## 210 603.8741 zona2 var3 F3 39.34468
## 211 796.1826 zona2 var4 F1 39.38262
## 212 845.7794 zona2 var4 F1 39.38938
## 213 701.4651 zona2 var4 F1 39.58753
## 214 873.2401 zona2 var4 F1 39.66896
## 215 617.2646 zona2 var4 F1 39.72169
## 216 809.9268 zona2 var4 F1 39.82298
## 217 602.0658 zona2 var4 F1 39.92928
## 218 655.7099 zona2 var4 F1 39.98433
## 219 838.2384 zona2 var4 F1 39.98643
## 220 892.4937 zona2 var4 F1 39.99260
## 221 1078.8332 zona2 var4 F2 40.26763
## 222 918.8794 zona2 var4 F2 40.32615
## 223 1049.7948 zona2 var4 F2 40.37543
## 224 937.5212 zona2 var4 F2 40.38792
## 225 893.4678 zona2 var4 F2 40.49265
## 226 908.9728 zona2 var4 F2 40.57494
## 227 933.3720 zona2 var4 F2 40.64259
## 228 942.3207 zona2 var4 F2 40.88538
## 229 931.9715 zona2 var4 F2 40.92257
## 230 915.6419 zona2 var4 F2 41.29410
## 231 924.8835 zona2 var4 F3 41.63557
## 232 971.0014 zona2 var4 F3 41.86286
## 233 924.0416 zona2 var4 F3 41.98487
## 234 893.9353 zona2 var4 F3 41.99595
## 235 928.9008 zona2 var4 F3 42.43944
## 236 928.7130 zona2 var4 F3 43.22740
## 237 919.0423 zona2 var4 F3 43.29603
## 238 909.5980 zona2 var4 F3 44.99739
## 239 892.3183 zona2 var4 F3 45.14669
## 240 922.4868 zona2 var4 F3 45.32652
tablani<-datatable(df_pd,class = 'cell-border stripe');tablani
library(ggpubr)
ggscatter(df_pd, x = "Fertilizacion", y = "Rendimiento", facet.by = c("Zona", "Variedad"), short.panel.labs = FALSE
)+ stat_smooth(method = "loess", span = 0.5)
## `geom_smooth()` using formula 'y ~ x'
En este esquema se puede observar el diseño del experimento, dos zonas, que son el bloqueo espacial, 4 parcelas principales por zona, y 3 subparcelas por parcelas principales para los tratamientos.
Para realizar un análisis y poder validar si hay diferencias entre los resultados, es necesario confirmar si los datos presentan una distribución normal, por lo cual se hará la prueba Shapiro
\[H_o: Datos \ con \ distribución \ normal\\H_a: Datos\ con\ distribución\ no\ normal \]
shapiro.test(df_pd$Rendimiento)
##
## Shapiro-Wilk normality test
##
## data: df_pd$Rendimiento
## W = 0.99498, p-value = 0.6177
Dados los resultados de la prueba, se rechaza la hipótesis alternativa, es decir los datos poseen una distribución normal, despues de esto, se va correr la prueba de análisis de varianza sin incluir la covariable de porosidad, con dos funciones distintas, con la función aov, y con la función sp.plot de la librería agricole,
#modelo aov
modelo_pd<-aov(rto_lech~bloque+vari*Fert+Error(bloque/vari), data = df_pd)
summary(modelo_pd)
##
## Error: bloque
## Df Sum Sq Mean Sq
## bloque 1 796663 796663
##
## Error: bloque:vari
## Df Sum Sq Mean Sq F value Pr(>F)
## vari 3 1151189 383730 2.145 0.273
## Residuals 3 536728 178909
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Fert 2 117577 58789 9.027 0.00017 ***
## vari:Fert 6 103490 17248 2.649 0.01678 *
## Residuals 224 1458750 6512
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#modelo agricole
library(agricolae)
##
## Attaching package: 'agricolae'
## The following object is masked from 'package:PMCMR':
##
## durbin.test
modelo_pdagr<-sp.plot(bloque,vari,Fert,rto_lech)
##
## ANALYSIS SPLIT PLOT: rto_lech
## Class level information
##
## vari : var1 var2 var3 var4
## Fert : F1 F2 F3
## bloque : zona1 zona2
##
## Number of observations: 240
##
## Analysis of Variance Table
##
## Response: rto_lech
## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 1 796663 796663 4.4529 0.1253
## vari 3 1151189 383730 2.1448 0.2735
## Ea 3 536728 178909 6.4470 2.874e-06 ***
## Fert 2 117577 58789 1.7803 0.2293
## vari:Fert 6 103490 17248 0.5223 0.7778
## Eb 8 264168 33021 6.4470 2.874e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 60.2 %, cv(b) = 25.9 %, Mean = 702.0989
Con los obtenidos, pareciera ser que la principales fuentes de diferencias son las variedades,los planes de fertilización, además de esto el efecto de los bloques parece poder llegar a generar alguna diferencia.
Sin embargo aún falta analizar los resultados con la covariable, pero es necesario primero validar si existe una relación entre la covariable y la variable respuesta, a este propósito es útil el índice de correlación de Pearson, el cual tiene un rango de [-1,1], siendo 0, indiferencia entre la covariable y la variable respuesta, y siendo 1, una fuerte correlación entres estos dos, que puede ser positiva o negativa, esto dependiendo el signo.
cor(x=df_pd$Rendimiento, y=df_pd$Porosidad, method = "pearson")
## [1] 0.7050391
Como la prueba lo indica el coeficiente es de 0.85, lo cual indica una correlación marcada entre el rendimiento y la porosidad del suelo, a mayor porosidad en el suelo, mayor el rendimiento.
library(plotly)
rer<-plot_ly(data = df_pd, x = df_pd$Porosidad, y = df_pd$Rendimiento,marker = list(size = 10,color = 'cyan',line = list(color = 'cyan3',width = 2)));rer
## No trace type specified:
## Based on info supplied, a 'scatter' trace seems appropriate.
## Read more about this trace type -> https://plot.ly/r/reference/#scatter
## No scatter mode specifed:
## Setting the mode to markers
## Read more about this attribute -> https://plot.ly/r/reference/#scatter-mode
rer <- rer %>% layout(title = 'Correlacion Porosidad vs Rendimiento',yaxis = list(title = 'Rendimiento'),xaxis = list(title = 'Porosidad'));rer
## No trace type specified:
## Based on info supplied, a 'scatter' trace seems appropriate.
## Read more about this trace type -> https://plot.ly/r/reference/#scatter
## No scatter mode specifed:
## Setting the mode to markers
## Read more about this attribute -> https://plot.ly/r/reference/#scatter-mode
Nos apoyamos en esta gráfica de correlación entre la porosidad y el rendimiento, donde podemos observar más fácilmente la relación que describe el coeficiente de pearson, se puede también observar cómo al extremo de la porosidad los resultados empiezan a dispersarse
Una vez confirmada la correlación entre estas dos variables, se realiza el análisis de varianza con la covariable, modelo el cual presenta el siguiente modelo:
\[Y_{ijk}=\mu+ \alpha_i+ \gamma_k+ \eta_{ ik} + \beta_j + (\alpha\beta){ij} + \epsilon{ijk}\] \(Y_{ijk}\) = Rendimiento
\(\alpha_i\) = Efecto fijo de la variedad
\(\gamma_k\) = Efecto fijo del bloque
\(\eta_{ik}\) = Error en parcela principal
\(\beta_j\) = Efecto fijo de la fertilización \(\alpha\beta_{ij}\) = Efecto interacción entre parcela principal y subparcela \(\epsilon_{ijk}\) = Error en subparcela
#modelo covar
mod_pdcoa<-aov(rto_lech~poros+bloque+vari*Fert+Error(bloque/vari), data = df_pd)
summary(mod_pdcoa)
##
## Error: bloque
## Df Sum Sq Mean Sq
## poros 1 796663 796663
##
## Error: bloque:vari
## Df Sum Sq Mean Sq F value Pr(>F)
## poros 1 1318239 1318239 8.085 0.105
## vari 3 43570 14523 0.089 0.960
## Residuals 2 326109 163054
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## poros 1 122630 122630 19.273 1.75e-05 ***
## Fert 2 29570 14785 2.324 0.1003
## vari:Fert 6 108735 18123 2.848 0.0108 *
## Residuals 223 1418882 6363
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Con el resultado de la prueba de análisis de covarianza, vemos que el valor F y el valor p indican que muy posiblemente la causa principal de diferencias en el experimento es la porosidad, y en menor medida el plan de fertilización y que las variedades
Estos resultados no son nada alentador para nuestro experimento, dado que las variables que estudiamos resultaron no ser las principales causantes de diferencias, en este caso es muy importante incluir la variable en el análisis de los resultados. Sin embargo dado que hay que escoger alguna opción una alternativa es validar que zona tiene mejor porosidad para el cultivo.
La porosidad en este caso es muy importante ya que el cultivo que es lechuga batavia, es una especie con raíces poco desarrolladas, por lo cual no tiene la fuerza para penetrar en suelos muy compactos, sin embargo en suelos muy porosos es posible que por lavados, el fertilizante que se aplicaba se pudiera perder, esa podria ser una explicación del porque en la gráfica de relación en la parte extrema, de mayor porosidad, los datos de rendimientos fueron tan dispersos.
corre<- plot_ly(data = df_pd, x = df_pd$Rendimiento, y = df_pd$Porosidad, color = df_pd$Zona)
corre<- corre %>% layout(title = 'Correlacion Porosidad vs Rendimiento por zonas',yaxis = list(title = 'Rendimiento'),xaxis = list(title = 'Porosidad'));corre
## No trace type specified:
## Based on info supplied, a 'scatter' trace seems appropriate.
## Read more about this trace type -> https://plot.ly/r/reference/#scatter
## No scatter mode specifed:
## Setting the mode to markers
## Read more about this attribute -> https://plot.ly/r/reference/#scatter-mode
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels
zonapor<-tapply(df_pd$Porosidad,df_pd$Zona,mean)
zonarto<-tapply(df_pd$Rendimiento,df_pd$Zona,mean)
data.frame(Porosidad = zonapor, Rendimiento = zonarto)
## Porosidad Rendimiento
## zona1 33.83104 644.4844
## zona2 38.33130 759.7134
La gráfica y la tabla muestran claramente una tendencia, a que la zona 2, que es la de falda de montaña presenta más porosidad y mayor rendimiento, en su literatura se reporta que la lechuga batavia crece de buena manera con una porosidad de entre 36 y 38, vemos que la media de la zona 2 es de 38.55, estando muy cerca de ese límite, con esto podríamos aconsejar al menos en este caso, sembrar lechuga batavia en la zona 2, dado que por la porosidad que hay en este suelo, el rendimiento podria ser mayor.
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-
Los diseños en parcelas divididas aumentan la información que se puede extraer de un experimento diseñado, siendo útil y práctico. El Modelo constan básicamente de la parcela principal o completa compuestas por los factores difíciles de cambiar y los factores de subparcela que son más fáciles de cambiar. De esta manera el número de parcelas completas es la cantidad de veces que restablece los factores de la parcela completa y el total de experimentos corresponde a las subparcelas.
En este tipo de experimento se cumple con la aleatorización, sin embargo, este es más estructurado, ya que ordena la ejecución del experimento al tener dos aleatorizaciones separadas las cuales determinan el orden en el que ejecuta las parcelas completas y recopila observaciones en cada subparcela, conduciendo a dos términos de error en el modelo general, que afecta el análisis. Una de las ventajas es que se puede estimar los parámetros del modelo por separado del error de la parcela completa y de la subparcela, que comúnmente son difíciles de estimar;
Se debe tener en cuenta el tipo de estudio en general que se desea desarrollar para las características cualitativas y cuantitativas en ambos niveles para tener un mejor desarrollo de este. En muchas situaciones, un experimento de parcela dividida puede producir estimaciones de parámetros de modelo más precisas que el experimento del mismo tamaño ejecutado en un diseño completamente aleatorizado.
Sobre lo que significa unidad experimental y unidad de observación
Short communication: On recognizing the proper experimental unit in animal studies in the dairy sciences
(https://www.sciencedirect.com/science/article/pii/S002203021630621X)
La unidad experimental también llamada unidad de replicación, es la unidad más pequeña a la que se aplica un tratamiento, las unidades experimentales muchas veces no difieren de manera fundamental, de modo que se obtienen inferencias confiables independientemente del tratamiento de cada unidad. Por otro lado, la unidad de observación, también conocida como unidad de muestreo, suele estar contenida en una unidad experimental, por lo que es una parte de esta, es la unidad sobre la que se realizan las observaciones o mediciones. Si las unidades de observación se anidan dentro de una unidad experimental, las unidades de observación se denominan comúnmente submuestras, pseudorreplicadas o réplicas técnicas. Cabe resaltar que estos dos tipos de unidades pueden coincidir siendo a la vez una unidad experimental y observacional. Según el ejemplo dado, la escala del diseño va a depender propiamente de lo que se quiere investigar, ya que se pueden tener unidades experimentales grandes como un corral o una finca con cierto número de vacas o, por el contrario, que la vaca sea la unidad experimental. En agronomía se podría tener un lote en donde se evalué el crecimiento de tomate en cien plantas, teniendo como unidad observacional 30 plantas para cada tratamiento.
https://online.stat.psu.edu/stat502/lesson/6/6.1-0
Guía para diseñar experimentos exitosos
El diseño de experimento es mejor verlo como una serie de decisiones y no una serie de pasos establecidos o receta como el mismo autor lo llama, existen tres tipos de experimentos, de observación, de medición y comparativos.Para los diseños comparativos se plantean hipótesis sobre un tema puntual, luego de esto se crean los modelos estadísticos y luego se llevan a campo en diseños experimentales, los cuales brindaran un respuesta a la hipótesis inicial.
Dada la importancia de este tipo de diseños, la inversión de capital y la inversión emocional a la hora de realizar estos experimentos, siempre se busca que los resultados sean positivos, es decir que se pueda demostrar que hay alguna diferencia, no es seguro que se pueda llegar a este tipo de conclusiones pero lo que sí se puede hacer es efectuar correctamente el diseño del experimento para mitigar cualquier fuente de error y que el resultado sea lo más verídico posible, el autor plantea cuatro pilares fundamentales en este proceso la replicación, la aleatorización, el bloqueo y las unidades experimentales.
Replicación: es un pilar muy importante en el diseño de experimentos, brinda la posibilidad de estimar el error experimental que es importante a la hora de validar si una hipótesis es verdadera o no , también brinda precisión sobre los resultados, por esto mismo es común que a mayor repeticiones mejor, igualmente a mayores repeticiones mayor rango se obtiene del experimento, llegando a detectar mayores condiciones sobre este, y finalmente se obtiene mayor control sobre el error del experimento.
Aleatorización: En este aspecto toma relevancia la elección de la muestra, de la población representativa, ya que cada tipo de individuo de la población debe tener un representante en la muestra, es por esto que la toma de muestra es tan importante en un diseño comparativo. A la hora ya de realizar el experimento y aplicar los tratamientos es igualmente importante que los tratamientos se hagan de forma aleatoria para evitar cualquier sesgo y perturbación sobre el experimento.
Bloques: Se usa en principio para poder moldear el tamaño y homogeneizar las unidades experimentales, el diseño RBCD es el más común, el tamaño del bloque debe ser igual a de los tratamientos, en otros casos como parcelas divididas, se usan por practicidad del experimento, cuando hay factores más grandes y más difíciles de controlar.
Unidad experimental: Este es un ítem el cual consiste en el tamaño de la parcela, por ejemplo, que tanto conviene tener una parcela grande o pequeña y esto como afecta a la varianza de los datos, hay un punto en el que aumentar el tamaño de la parcela tendrá pequeños cambios en la varianza, en ese caso es importante definir si es realmente necesario ese aumento, o si por el contrario el cambio de la varianza es grande.
La investigación biológica, que se apoya en experimentos comparativos, son estudios caros en su mayoría , y difíciles de controlar en su totalidad, la única manera de obtener resultados que se puedan aplicar y que sean válidos es controlando lo que más se pueda con el mínimo de presupuesto, es decir optimizando los procesos y decisiones que se tomen a la hora de diseñar y llevar a cabo el experimento
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:
Articulo:Producción y rendimiento de maíz en cuatro tipos de labranza bajo condiciones de temporal.
En el artículo informan acerca de las zonas en las que van a trabajar, que son dos, indican las características físicas y climáticas de estas, y porque no son controlables, ya que la textura del suelo y el clima de la zona son bastante dispendiosas y caras de controlar, sin embargo no es explícito que estos sean los bloques del experimento, informan que van usar 4 parcelas y 4 subparcelas cada una con dos repeticiones, pero el esquema del experimento no es muy claro especialmente en la parcela de La Croix, ya que usan un a nomenclatura que no está definida en el texto y que no usan en las otras parcelas, siendo un poco confuso.
El factor de las parcelas principales es el subsoleo, el cual consiste en labrar el suelo a 40 cm de profundidad, y el factor de las subparcelas es el tipo de labranza, el cual es superficial en todos los casos, por este motivo deciden que este sea el subfactor, dado es menos profundo que el factor principal, no supera los 15 cm.
El artículo carece de la revisión de supuestos, lo cual es una grave falla, ya que no realizan pruebas para validar el tipo de distribución de los datos , ni si se presenta igualdad de varianzas, por lo cual no es seguro realizar la Anova ya que no hay seguridad que los datos se encuentren en la curva F, esto le quita algo de credibilidad de los datos, es importante verificar estos supuestos antes de cualquier análisis.
Las tablas no muestran los grados de libertad, estos son importantes para conocer cuantos datos manejaron para realizar el análisis, otro aspecto para detallar es que no aparece el error de las repeticiones explícito en las tablas.
El artículo presenta 6 tablas de análisis de varianza, dos por parcelas principales, ambas con subsoleo, la diferencia es que el alfa en una es de 0.01 y en la otra es de 0.05, en el primer caso no hay diferencias significativas, pero en el segundo si y 4 por subparcelas, dos con subsoleo y dos sin subsoleo, igualmente la diferencia es el alfa. No presentan una tabla de análisis de varianza con los 2 factores, que puedan informar acerca de la interacción de estos, esta tabla es muy útil en este caso ya que evitaría tantas tablas en ela rticulo, además de esto no se entiende porque usan dos alfas, especialmente en el caso de las subparcelas, ya que aquí la conclusión es la misma ya sea con un alfa de 0.01 o de 0.05.
No usan un método para comparar las medias, pero es necesario hacer una prueba en el caso de las parcelas principales, dado que con un alfa de 0.01 si hay diferencia entre parcelas, sin embargo no se hizo, para las subparcelas dado que no hay diferencia no era necesario hacer esta comparación de medias.
Presentan tablas pero separadas, esto no ayuda a analizar las interacciones entre factores, sería más útil y práctico si hubieran hecho la Anova con los dos factores directamente.
Estos están definidos pero no se presentan explícitamente como bloques, y en el análisis de varianza no se ven por ningún lado.
En todos los factores se encuentra que es balanceada, ya que todos tenían la misma cantidad de datos, o por lo menos eso indica el artículo
Este punto es de los mejores desarrollados por el artículo, explican a cabalidad los procedimientos que van a llevar a cabo para establecer la unidad, la distancia entre surcos, las medidas de los lotes, las áreas de muestreo, entre otros.
En ninguna parte del artículo informan que software utilizan, es más no indican si hicieron el análisis por computadora o a mano, esto no ayuda a validar si los datos fueron bien tratados, ya que no se puede hacer seguimiento del procedimiento llevado a cabo.
Es interesante observar la distribución de las parcelas, ya que están en distintas zonas, y son de diferentes formas, otro punto es que para determinar el rendimiento midieron los granos por mazorca y las mazorcas por planta, además de esto en el estudio le dan más relevancia al valor F que al valor P.