\[ H_0: \mu_{CBD(A1)} = \mu_{CBD(A1)} = \mu_{CBD(A2)} = \mu_{CBD(A3)} = \mu_{CBD(A4)} \\ H_a: At~least~one~is~different\]

set.seed(1001280374)
CBD = c(rnorm(n = 30, mean = 60, sd = 2),
        rnorm(n = 30, mean = 50, sd = 1.5),
        runif(n = 30, min = 30, max = 40),
        rnorm(n = 30, mean = 40, sd = 3))

Accesion = gl(n = 4, k = 30, length = 120, c("A1", "A2", "A3", "A4"))


xy = expand.grid(x = 1:10, y = 1:12)


df1 = data.frame( CBD, Accesion)


aleat = sample(1:120, 120, FALSE)
df1$cord_x = xy$x[aleat]
df1$cord_y = xy$y[aleat]
head(df1)
##        CBD Accesion cord_x cord_y
## 1 58.90225       A1      7      5
## 2 61.73631       A1      8      8
## 3 59.67711       A1      9      9
## 4 62.05145       A1      4      4
## 5 61.45322       A1      9     12
## 6 63.64721       A1      6      6

##Estadística descriptiva

library(dplyr)
df1 %>% 
  group_by(Accesion) %>% 
  summarise(media = mean(CBD),
            desv = sd(CBD),
            cv = 100*desv/media)
## # A tibble: 4 x 4
##   Accesion media  desv    cv
##   <fct>    <dbl> <dbl> <dbl>
## 1 A1        60.0  2.22  3.70
## 2 A2        49.7  1.61  3.24
## 3 A3        34.7  2.71  7.81
## 4 A4        39.9  2.11  5.28

Estadística inferencial

mod1 = aov(CBD ~ Accesion, data = df1)
s_mod1 = summary(mod1)
p_valor0= s_mod1[[1]][1,5]
ifelse(p_valor0<0.05, "Rechazo H0", "No rechazo H0")
## [1] "Rechazo H0"
# Revisando suspuestos del modelo

# Normalidad de residuales
shapiro = shapiro.test(mod1$residuals)
p_valor1= shapiro$p.value
ifelse(p_valor1<0.05, "Rechazo H0", "No rechazo H0")
## [1] "No rechazo H0"
# Igualadad de varianzas
barlett = bartlett.test(mod1$residuals, df1$Accesion)
p_valor2= barlett$p.value
ifelse(p_valor2<0.05, "Rechazo H0", "No rechazo H0")
## [1] "No rechazo H0"

Matriz de distacia

dist_matrix <- as.matrix(dist(cbind(df1$cord_x, df1$cord_y)))#;dist_matrix
which.max(dist_matrix) # Es la posición de la máxima
## [1] 5228
max(dist_matrix) # Mayor valor 
## [1] 14.21267
min(dist_matrix) # Menor valor
## [1] 0
dim(dist_matrix)
## [1] 120 120

Inversa de la distancia

dist_matrix_inv <- 1 / dist_matrix # Element wise
diag(dist_matrix_inv) <- 0
which.max(dist_matrix_inv) # Es la posición de la máxima
## [1] 45
max(dist_matrix_inv) # Mayor valor 
## [1] 1
min(dist_matrix_inv) # Menor valor
## [1] 0

Indice de Moran

Moran = Moran.I(mod1$residuals, dist_matrix_inv)# p.value < 0.05 <- Dependencia espacial
p_valor3= Moran$p.value
ifelse(p_valor3<0.05, "Rechazo H0", "No rechazo H0")
## [1] "No rechazo H0"

\[Conclusión: Se~ cumplen~ todos~ los~ supuestos~ menos~ el~ de~ ANOVA\]

library(writexl)
## Warning: package 'writexl' was built under R version 4.1.3
write_xlsx(df1, 'datos_CE_090521.xlsx')