FSBGCA & diseƱo desvalanceado
Anderson Bayardo Palma Oliva
2023-04-25
#Factorial simple con bloques generalizados completamente al azar #DiseƱo desbalanceado
modelo \[y_{ijk} = \mu + \tau_i + \beta_{j} + \tau\beta_{ij} + \epsilon _ {ijk}\] \(i=1,2,\dots,a:\text{numero de tratamientos}\) \(a\) numero de nieveles del factor \(j=1,2,\dots,b\) \(k=1,2,\dots,r_i\) \(r_i\)repeticion de cada tratamiento
set.seed(123)#el numero de datos no es igual en cada tratamiento
aceite = c(
rnorm(12, 10, 0.8),
rnorm(12, 11, 0.78),
rnorm(11, 9, 0.70),
rnorm(12, 10, 0.8),
rnorm(12, 11, 0.78),
rnorm(11, 9, 0.70)
)
bloque = gl(2,35,70, c('b1','b2'))
metodo = rep(rep(c('T1','T2','T3'), c(12,12, 11)),2 )
datos = data.frame(metodo, bloque, aceite)
head(datos)## metodo bloque aceite
## 1 T1 b1 9.551619
## 2 T1 b1 9.815858
## 3 T1 b1 11.246967
## 4 T1 b1 10.056407
## 5 T1 b1 10.103430
## 6 T1 b1 11.372052table(datos$metodo, datos$bloque)##
## b1 b2
## T1 12 12
## T2 12 12
## T3 11 11#tratamiento: 3 metodos de extraer aceite escencial mod1 = aov(aceite ~ bloque * metodo, datos)
summary(mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 1 0.05 0.054 0.107 0.745
## metodo 2 41.24 20.622 40.834 3.72e-12 ***
## bloque:metodo 2 0.65 0.323 0.640 0.531
## Residuals 64 32.32 0.505
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#El unico p-valor interpretable es metodo
#como es menor a 0.5 se rechaza la hipotesis nula (las medias no son iguales de aceite escencial extraido por todos lo metodos )Corriendo como si fuera balanceado}
mod1 = aov(aceite ~ bloque * metodo, datos)
summary(mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 1 0.05 0.054 0.107 0.745
## metodo 2 41.24 20.622 40.834 3.72e-12 ***
## bloque:metodo 2 0.65 0.323 0.640 0.531
## Residuals 64 32.32 0.505
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1mod1a = aov(aceite ~ metodo * bloque, datos)
summary(mod1a) #es exactamente lo mismo, solo cambia el orden (solo cuando es un solo factor) ## Df Sum Sq Mean Sq F value Pr(>F)
## metodo 2 41.24 20.622 40.834 3.72e-12 ***
## bloque 1 0.05 0.054 0.107 0.745
## metodo:bloque 2 0.65 0.323 0.640 0.531
## Residuals 64 32.32 0.505
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Corriendo desbalanceado
mod2 = anova(lm(aceite ~ bloque * metodo, datos))
mod2## Analysis of Variance Table
##
## Response: aceite
## Df Sum Sq Mean Sq F value Pr(>F)
## bloque 1 0.054 0.0541 0.1071 0.7446
## metodo 2 41.243 20.6217 40.8342 3.716e-12 ***
## bloque:metodo 2 0.646 0.3231 0.6399 0.5307
## Residuals 64 32.321 0.5050
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1mod2a = anova(lm(aceite ~ metodo * bloque, datos))
mod2a## Analysis of Variance Table
##
## Response: aceite
## Df Sum Sq Mean Sq F value Pr(>F)
## metodo 2 41.243 20.6217 40.8342 3.716e-12 ***
## bloque 1 0.054 0.0541 0.1071 0.7446
## metodo:bloque 2 0.646 0.3231 0.6399 0.5307
## Residuals 64 32.321 0.5050
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#los bloques no son factores experimentales, son simplemente restricciones para aleatorizar lote <- c(rep("lote1",1),
rep("lote2",1),
rep("lote3",1),
rep("lote4",1),
rep("lote5",1))
genotipo <- c(rep("gtA",5),
rep("gtB",5),
rep("gtC",5),
rep("gtD",5),
rep("gtE",5))
prov <- c("A","E","C","B","D",
"C","B","A","D","E",
"B","C","D","E","A",
"D","A","E","C","B",
"E","D","B","A" ,"C")
biom <- c(42,45,41,56,47,
47,54,46,52,49,
55,52,57,49,45,
51,44,47,50,54,
44,50,48 ,43,46)
mydata <- data.frame(lote, genotipo, prov, biom);mydata## lote genotipo prov biom
## 1 lote1 gtA A 42
## 2 lote2 gtA E 45
## 3 lote3 gtA C 41
## 4 lote4 gtA B 56
## 5 lote5 gtA D 47
## 6 lote1 gtB C 47
## 7 lote2 gtB B 54
## 8 lote3 gtB A 46
## 9 lote4 gtB D 52
## 10 lote5 gtB E 49
## 11 lote1 gtC B 55
## 12 lote2 gtC C 52
## 13 lote3 gtC D 57
## 14 lote4 gtC E 49
## 15 lote5 gtC A 45
## 16 lote1 gtD D 51
## 17 lote2 gtD A 44
## 18 lote3 gtD E 47
## 19 lote4 gtD C 50
## 20 lote5 gtD B 54
## 21 lote1 gtE E 44
## 22 lote2 gtE D 50
## 23 lote3 gtE B 48
## 24 lote4 gtE A 43
## 25 lote5 gtE C 46library(collapsibleTree)
collapsibleTreeSummary(mydata,
c('lote',
'prov',
'genotipo',
'biom'),
collapsed = FALSE)library(ggplot2)
ggplot(mydata)+
aes(biom, genotipo)+
geom_point(size=5, shape=15)+
facet_grid(lote ~ prov)
ggplot(mydata)+
aes(lote, genotipo, fill=biom)+
geom_tile()+
coord_equal()+
facet_wrap( ~ prov, nrow=1)+
theme(axis.text.x = element_text(angle=90))
library(lattice)
bwplot(biom ~ genotipo | prov + lote,
mydata)
Modelo \[ y_{ijk}= \mu + \tau_i + \beta_{j} +
\delta_{k} + \epsilon_{ijk}\\
i=1\dots\]
tbl = matrix(mydata$prov, 5)
colnames(tbl) = unique(mydata$genotipo)
rownames(tbl) = unique(mydata$lote)
tbl## gtA gtB gtC gtD gtE
## lote1 "A" "C" "B" "D" "E"
## lote2 "E" "B" "C" "A" "D"
## lote3 "C" "A" "D" "E" "B"
## lote4 "B" "D" "E" "C" "A"
## lote5 "D" "E" "A" "B" "C"mod <- lm(biom ~ lote + genotipo + prov,
mydata)
anova(mod)## Analysis of Variance Table
##
## Response: biom
## Df Sum Sq Mean Sq F value Pr(>F)
## lote 4 17.76 4.440 0.7967 0.549839
## genotipo 4 109.36 27.340 4.9055 0.014105 *
## prov 4 286.16 71.540 12.8361 0.000271 ***
## Residuals 12 66.88 5.573
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\[H0 : \mu_{B_{g_1}} = \mu_{B_{g_2}} = \mu_{B_{g_3}} = \mu_{B_{g_4}} = \mu_{B_{g_5}} \]
bwplot(biom ~genotipo | prov,
mydata)
library(ggplot2)
ggplot(mydata)+
aes(genotipo,
biom,
fill=prov)+
geom_col(
position = 'dodge')
res_mod = mod$residuals
#1.Normalidad de residuales
shapiro.test(res_mod) #SE cumple el supuesto##
## Shapiro-Wilk normality test
##
## data: res_mod
## W = 0.97691, p-value = 0.8178#2. Igualdad de varainzas (se cumple supuesto, variables iguales)
bartlett.test(res_mod,
genotipo)##
## Bartlett test of homogeneity of variances
##
## data: res_mod and genotipo
## Bartlett's K-squared = 5.9223, df = 4, p-value = 0.205library(TukeyC)
tt = TukeyC(mod, 'genotipo')
plot(tt, lwd=2, cex=2)