Diseño Factorial Simple en Areglo completamente al Azar, en Bloques completos y Bloques completos generalizados
Erik David Aldana
24/4/2021
Diseño Factorial Simple en Areglo completamente al Azar (Teniendo en cuenta el Efecto Borde)
library(readxl)
library(ggplot2)
df <- read_excel("C:\\Users\\ERICK\\OneDrive\\Documentos\\Trabajos en R\\data_DE.xlsx")dim(df)## [1] 96 5head(df)## # A tibble: 6 x 5
## x y blk trt respuesta
## <dbl> <dbl> <chr> <chr> <dbl>
## 1 0 0 B2 T3 2.11
## 2 50 0 B1 T2 2.11
## 3 100 0 B1 T2 2.11
## 4 150 0 B2 T3 2.14
## 5 200 0 B1 T2 2.13
## 6 250 0 B1 T2 2.19plot(df$x, df$y, col = df$respuesta, pch=19, )
ggplot(df)+
aes(x,y,fill=respuesta)+
geom_tile()
#Asumiendo los datos como un diseño FSCA
df$blk = as.factor(df$blk)
df$trt = as.factor(df$trt)mod1 = aov(respuesta ~ trt, df)
res_mod1 = summary(mod1)
res_mod1## Df Sum Sq Mean Sq F value Pr(>F)
## trt 2 0.0273 0.01366 0.402 0.67
## Residuals 93 3.1622 0.03400# Prueba homocedasticidad
bartlett.test(df$respuesta, df$trt)##
## Bartlett test of homogeneity of variances
##
## data: df$respuesta and df$trt
## Bartlett's K-squared = 0.44126, df = 2, p-value = 0.802# Prueba de normalidad a los residuales
shapiro.test(mod1$residuals)##
## Shapiro-Wilk normality test
##
## data: mod1$residuals
## W = 0.84972, p-value = 1.912e-08# Prueba de simetria
library(normtest)
skewness.norm.test(mod1$residuals)##
## Skewness test for normality
##
## data: mod1$residuals
## T = -0.22351, p-value = 0.323# library(symmetry)
# symmetry_test(mod1$residuals, stat = 'MOI',
# bootstrap = TRUE, k = 3, mu = 0)
# Histograma de residuales
hist(mod1$residuals)
ggplot(df)+
aes(x,y,fill=mod1$residuals)+
geom_tile()
plot(mod1$residuals, pch = 16)
#Indice de Moran para saber si los datos tienen dependencia espacial
d = as.matrix(dist(df[, 1:2]))
dim(d)## [1] 96 96d_inv = 1/d
diag(d_inv) = 0
W = d_inv
W[1:5, 1:5]## 1 2 3 4 5
## 1 0.000000000 0.020000000 0.01 0.006666667 0.005000000
## 2 0.020000000 0.000000000 0.02 0.010000000 0.006666667
## 3 0.010000000 0.020000000 0.00 0.020000000 0.010000000
## 4 0.006666667 0.010000000 0.02 0.000000000 0.020000000
## 5 0.005000000 0.006666667 0.01 0.020000000 0.000000000library(ape)
Moran.I(mod1$residuals, W)## $observed
## [1] 0.06393269
##
## $expected
## [1] -0.01052632
##
## $sd
## [1] 0.01181971
##
## $p.value
## [1] 2.984915e-10Diseño Factorial Simple en Areglo en Bloques Completos al Azar (Sin Interacción)
library(readxl)
df1<- read_excel("C:\\Users\\ERICK\\OneDrive\\Documentos\\Trabajos en R\\Ejemplos bloques.xlsx", sheet = "DATA")
head(df1)## # A tibble: 6 x 5
## x y trt bloq clorofila
## <dbl> <dbl> <chr> <dbl> <dbl>
## 1 1 10 A 10 213.
## 2 1 9 B 9 258.
## 3 1 8 A 8 214.
## 4 1 7 C 7 296.
## 5 1 6 A 6 238.
## 6 1 5 B 5 271.mod2 = aov(clorofila ~ bloq+trt, df1)
summary(mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## trt 2 0.0273 0.01366 0.402 0.67
## Residuals 93 3.1622 0.03400# omitiendo los bloques
mod2a = aov(clorofila ~ trt, df1)
summary(mod2a)## Df Sum Sq Mean Sq F value Pr(>F)
## trt 2 20220 10110 8.465 0.0014 **
## Residuals 27 32248 1194
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1St_mod2= shapiro.test(mod2$residuals); St_mod2##
## Shapiro-Wilk normality test
##
## data: mod2$residuals
## W = 0.9505, p-value = 0.1744St_mod2a= shapiro.test(mod2a$residuals); St_mod2a##
## Shapiro-Wilk normality test
##
## data: mod2a$residuals
## W = 0.95185, p-value = 0.1895ifelse(St_mod2$p.value<0.05,"No son normales","Son normales")## [1] "Son normales"ifelse(St_mod2a$p.value<0.05,"No son normales","Son normales")## [1] "Son normales"Bt_mod2= bartlett.test(df1$clorofila, df1$trt); Bt_mod2##
## Bartlett test of homogeneity of variances
##
## data: df1$clorofila and df1$trt
## Bartlett's K-squared = 1.9975, df = 2, p-value = 0.3683ifelse(Bt_mod2$p.value<0.05,"Heterocedasticidad","Homocedasticidad")## [1] "Homocedasticidad"Bt_mod2a= bartlett.test(df1$clorofila, df1$trt); Bt_mod2a##
## Bartlett test of homogeneity of variances
##
## data: df1$clorofila and df1$trt
## Bartlett's K-squared = 1.9975, df = 2, p-value = 0.3683ifelse(Bt_mod2a$p.value<0.05,"Heterocedasticidad","Homocedasticidad")## [1] "Homocedasticidad"library(ape)
d = as.matrix(dist(df1[,1:2]))
di = 1/d
diag(di) = 0
im_mod2 = Moran.I(mod2$residuals, di); im_mod2## $observed
## [1] -0.03336448
##
## $expected
## [1] -0.03448276
##
## $sd
## [1] 0.03222058
##
## $p.value
## [1] 0.9723135ifelse(im_mod2$p.value<0.05,"Hay Dependencia espacial","No hay dependencia espacial")## [1] "No hay dependencia espacial"im_mod2a = Moran.I(mod2a$residuals, di); im_mod2a## $observed
## [1] 0.3964863
##
## $expected
## [1] -0.03448276
##
## $sd
## [1] 0.03222764
##
## $p.value
## [1] 0ifelse(im_mod2a$p.value<0.05,"Hay Dependencia espacial","No hay dependencia espacial")## [1] "Hay Dependencia espacial"library(ggplot2)
ggplot(df1)+aes(x,y,fill=clorofila,col="yellow")+geom_tile()
DISEÑO FACTORIAL SIMPLE EN BLOQUES COMPLETOS GENERALIZADOS AL AZAR
library(readxl)
df2 <- read_excel("C:\\Users\\ERICK\\OneDrive\\Documentos\\Trabajos en R\\Ejemplos bloques.xlsx", sheet = "DATA 2")
head(df2)## # A tibble: 6 x 6
## trt clorofila x y bloq MO
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 A 161. 1 10 10 2.43
## 2 C 235. 1 9 9 3.19
## 3 C 283. 1 8 8 3.05
## 4 B 214. 1 7 7 3.29
## 5 C 365. 1 6 6 3.34
## 6 A 182. 1 5 5 2.28# Corriendo el DFSBCGA como la gente creeria que esta bien
mod3 = aov(clorofila ~ bloq + trt, df2)
summary(mod3)## Df Sum Sq Mean Sq F value Pr(>F)
## bloq 1 412 412 0.268 0.607
## trt 2 57186 28593 18.626 6.29e-07 ***
## Residuals 56 85966 1535
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(lattice)
bwplot(clorofila ~ trt | bloq, df2)
# Corriendo el DFSBCGA como se deberia hacer
mod3a = aov(clorofila ~ bloq * trt, df2)
mod3b = aov(clorofila ~ trt, df2)
summary(mod3a)## Df Sum Sq Mean Sq F value Pr(>F)
## bloq 1 412 412 0.274 0.603
## trt 2 57186 28593 19.001 5.65e-07 ***
## bloq:trt 2 4705 2353 1.563 0.219
## Residuals 54 81261 1505
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(mod3b)## Df Sum Sq Mean Sq F value Pr(>F)
## trt 2 57186 28593 18.87 5.15e-07 ***
## Residuals 57 86377 1515
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tapply(df2$clorofila, list(df2$bloq,df2$trt), var)## A B C
## 1 1309.1287 895.86335 7066.4767
## 2 299.0431 514.32739 2437.6010
## 3 462.3402 1191.92509 649.2597
## 4 1251.8590 242.97881 853.2411
## 5 1523.5430 844.54228 1294.7446
## 6 131.3348 578.75916 6553.4374
## 7 439.0367 4381.57283 103.9868
## 8 1267.8544 725.17645 1744.5907
## 9 208.3461 698.58351 9204.5790
## 10 4512.2318 9.43876 3669.1749res_mod3a = mod3a$residuals
res_mod3b = mod3b$residualsshapiro.test(res_mod3a)##
## Shapiro-Wilk normality test
##
## data: res_mod3a
## W = 0.96219, p-value = 0.06026shapiro.test(res_mod3b)##
## Shapiro-Wilk normality test
##
## data: res_mod3b
## W = 0.97232, p-value = 0.1893bartlett.test(res_mod3a, df2$trt)##
## Bartlett test of homogeneity of variances
##
## data: res_mod3a and df2$trt
## Bartlett's K-squared = 5.7576, df = 2, p-value = 0.0562bartlett.test(res_mod3b, df2$trt)##
## Bartlett test of homogeneity of variances
##
## data: res_mod3b and df2$trt
## Bartlett's K-squared = 5.1021, df = 2, p-value = 0.078d = as.matrix(dist(df2[,c('x','y')]))
di = 1/d
diag(di) = 0
im_mod3a = Moran.I(res_mod3a, di)
im_mod3a## $observed
## [1] 0.03988513
##
## $expected
## [1] -0.01694915
##
## $sd
## [1] 0.0161582
##
## $p.value
## [1] 0.0004358559im_mod3b = Moran.I(res_mod3b, di)
im_mod3b## $observed
## [1] 0.03930886
##
## $expected
## [1] -0.01694915
##
## $sd
## [1] 0.01616128
##
## $p.value
## [1] 0.0004994784