Parcial Diseno de experimentos

Ejercicio 1

Cedula: 1073243378

Para llenar el ANOVA tenemos los siguientes datos:

\[ n = 72 \\ Tratamientos = 9\\ Repeticiones =8\\ S^2 = 883\\ Sc_{tratamiento}=6000 \] Por lo que es necesario encontrar La suma de cuadrados del error.

Tenemos dos formulas que pueden ayudarnos

\[ _{(1)} S^2=\frac {SC_{Total}}{dl} \\ _{(2)}SC_{total} = SC_{trat}+ SC_{error} \]

reemplazando (2) en (1) obtenemos:

\[ S^2=\frac {SC_{trat}+ SC_{error}}{dl} \] despejamos La suma de cuadrados del error

\[ SC_{error}= S^2\cdot dl-SC_{trat} \] Con esto obtenemos la suma total de cuadrados, y también las variables faltantes para poder completar la tabla.

SC_wit = 883*71-6000; SC_wit

## [1] 56693

SC_Tot= SC_wit+6000 ; SC_Tot

## [1] 62693

df_bet= 9-1 ; df_bet

## [1] 8

df_wit= 9*(8-1) ; df_wit

## [1] 63

Msquare_Bet= 6000/df_bet ; Msquare_Bet 

## [1] 750

Msquare_wit= SC_wit/df_wit ; Msquare_wit

## [1] 899.8889

f= Msquare_Bet/Msquare_wit ; f

## [1] 0.8334362

	Suma de cuadrados	df	Cuadrados medios	F
Entre Trt	6000	8	750	0.83
———-	—————–	—-	——————	—
Inrtra Trt	56693	63	899.89
———-	—————–	—-	——————	—
Total	62693	71
———-	—————–	—-	——————	—

Si el F tabulado es 2.8. ¿qué puede decirse acerca de la Hipótesis nula de igualdad de los promedios del índice en todas las condiciones de tratamiento (use el p valor así como el cociente F calculado de la tabla para concluir)?

Hipotesis nula:

\[H_o= \mu_{Trt\ 1}= \mu_{Trt\ 2}=...= \mu_{Trt\ 9} \] Si la pobabilidad Pr(>F) es > 0.05 se acepta Ho

p_valor= pf(q=f, df1= 8, df2 = 63 , lower.tail = F ); p_valor

## [1] 0.5766744

Con un F tabulado de 2.8,la hipotesis nula no se rechaza, ya que el valor del F calculado (0.8) se encuentra por debajo de este valor. y el p valor (0.58) respalda estos resultados al ser mayor a 0.05.

Sin embargo para fiabilidad de la respuesta, se evaluara una simulación de los datos para evaluar las condiciones de normalidad, de varianza y de independencia.

b.genere unos datos para cada tratamiento que tengan una varianza de 29.21 y una media que desde el primer al último tratamiento cambie a razón de 10 unidades del índice, empezando con 408

set.seed(1995) #Semilla 

#datos de medición de clorofila
clorofila =c(rnorm(8,408,sqrt(853)),rnorm(8,418,sqrt(853)),rnorm(8,428,sqrt(853)),rnorm(8,438,sqrt(853)),rnorm(8,448,sqrt(853)),rnorm(8,458,sqrt(853)),
             rnorm(8,468,sqrt(853)),rnorm(8,478,sqrt(853)),
             rnorm(8,488,sqrt(853)))
medida= gl(9,8,72,labels=c("clor_1","clor_2","clor_3","clor_4","clor_5","clor_6","clor_7","clor_8","clor_9"))

#media de clorofila tomada en cada lote. 
med_trt = tapply(clorofila, medida, mean); med_trt 

##   clor_1   clor_2   clor_3   clor_4   clor_5   clor_6   clor_7   clor_8 
## 417.3451 421.8985 400.3088 424.5807 438.7553 455.1485 474.6046 486.6401 
##   clor_9 
## 489.8674

#media global
med_global= mean(med_trt); med_global

## [1] 445.461

#boxplot
boxplot=boxplot(clorofila~medida)
points(med_trt, col= 'blue3', pch=10)
abline(h= med_global, col='red1')

mod_1= aov(clorofila~medida);
smod=summary(mod_1); smod

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## medida       8  67810    8476   12.12 2.66e-10 ***
## Residuals   63  44064     699                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#pruebas de comparaciones multiples

Observamos que las variaciones pueden deberse a las respuestas de los tratamientos y que las medias no son iguales para los datos.

#normalidad resifuales
resid= mod_1$residuals; 
shapiro.test(resid)

## 
##  Shapiro-Wilk normality test
## 
## data:  resid
## W = 0.97502, p-value = 0.1619

se acepta que los residuales siguen una distribucion normal

#igualdad de varianza 
bartlett.test(resid~medida)

## 
##  Bartlett test of homogeneity of variances
## 
## data:  resid by medida
## Bartlett's K-squared = 10.336, df = 8, p-value = 0.2422

Se acepta la hipotesis de homogeneidad de las varianzas

#independencia
plot(resid,pch=16)

como se observa en el grafico, no hay patrones que se puedan atribuir a coorelación de los datos de manera espacial ( u otras variaciones que puedan causarla)

Aunque en el box plot observado los datos de cada tratamiento pueden estar por encima o por debajo de la media global y parecer algunas muy diferentes con respecto a otras, al evaluarlas en conjunto y al realizar los test de normalidad, varianza e independencia de los datos, no se rechaza la hipotesis nula es decir el promedio de del indice de clorofila de los 9 tratamientos es igual.

Ejercicio 2

library(readxl) 
datos=Copia_de_Datos_algodn <- read_excel("Copia de Datos algodn.xlsx")
View(datos)
#se asignan factores, niveles y bloques para los datos
perdidas= (datos$Perdida)
metodo = gl(5,6,30, labels = c("M1","M2","M3","M4","M5"))
granja= factor(rep(1:6,5))

Una primera observación que debe hacerse, consiste en que en la descripción del experimento apunta a que el investigador bloqueó las granjas de donde provino el algodón ya que quería tener en cuenta el efecto de los cultivadores en los análisis, sin embargo, luego se describe que los fardos de algodón provenientes de cada granja fueron mezclados con lo cual el bloquear por granja no tendría sentido.

#ANOVA desbalanceado 

dataframe=data.frame(perdidas,metodo,granja)
mod2_1=aov(lm(perdidas~metodo+granja)); summary(mod2_1)

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## metodo       4  47.76  11.941   8.649 0.000375 ***
## granja       5 139.66  27.932  20.233 5.16e-07 ***
## Residuals   19  26.23   1.381                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 1 observation deleted due to missingness

#cambiando el orden 
mod2_2=aov(lm(perdidas~granja+metodo)); summary(mod2_2)

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## granja       5 138.30  27.661  20.037 5.57e-07 ***
## metodo       4  49.12  12.280   8.895 0.000319 ***
## Residuals   19  26.23   1.381                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 1 observation deleted due to missingness

#Grafico para comparar los datos
plot(perdidas~metodo+granja)

Podemos observar que al cambiar el orden del bloque y los tratamientos en el ANOVA, los datos no cambian drasticamente. Observando el F value del factor de bloqueo (granja), observamos que es un numero relativamente mayor que el F value, por lo que gran parte de la variacion de los datos pueden ser atribuidos al factor de bloqueo, pudiendo existir diferencia en la calidad de la fibra algodon producida por los granjeros.

Segun el p valor, existe diferencia entre al menos dos metodos de limpieza.

b.Estimar el valor del dato faltan de del primer granjero Usando el promedio de los datos de los grandejos 2 a 6 y realice analisis de varianza para probar las diferencias en las pérdidas medias de peso para los cinco métodos de limpiado de las fibras de algodón. Compare este resultado con el caso desbalanceado (de ser posible).

#media de granja 2:6
fna=(6.75+13.05+10.26+8.01+8.42)/5 ; fna

## [1] 9.298

# Reemplazando el dato
library(readxl)
algodon_corr <- read_excel("C:/Users/Usuario/Desktop/algodon corr.xlsx")
View(algodon_corr)
metodo_c = gl(5,6,30, labels = c("M1","M2","M3","M4","M5"))
perdidas_c = algodon_corr$Perdida
granja= factor(rep(1:6,5))

df_c=data.frame(metodo_c,perdidas_c)
head(df_c)

##   metodo_c perdidas_c
## 1       M1      9.298
## 2       M1       6.75
## 3       M1      13.05
## 4       M1      10.26
## 5       M1       8.01
## 6       M1       8.42

# Modelo sin el dato faltante
d_22=aov(perdidas_c~metodo_c+granja) ; summary(d_22) 

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## metodo_c     4  51.15  12.787   9.641 0.000164 ***
## granja       5 139.36  27.873  21.015 2.45e-07 ***
## Residuals   20  26.53   1.326                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

En el modelo corregido, de igual manera se rechaza la hipotesis nula, es decir, existen diferencias entre los métodos de procesamiento del algodón.

Ejercicio 3

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.U+0.1 y de la capa inferior osciló entre 2 y 2.T+0.2. 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. Genere 50 datos en cada capa. Use la función sort.int de R para ordenar los datos de cada capa con la opción partial=25+U dentro de la propia función sort.int. Una vez cree los datos realice algún diagrama de color (preferiblemente 3D) que permita visualizar las medidas de carbono en cada capa generadas por computadora. Compare si se encuentran diferencias en la media de carbono entre capas utilizando un nivel de confianza del 95%.

cc: 1073243378

#Se fija la semilla y se crean datos con distribución uniforme, los cuales se ordenan parcialmente con la funcion sort.in

set.seed(3000)
CO_5=sort.int(runif(n = 50 ,min = 3.0, max = 3.9),partial = 33)
CO_10= sort.int(runif(n = 50, min = 2, max = 3.0),partial = 33)
CO= c(CO_5,CO_10)
#Se genera la ventana de observación del suelo 
observ= expand.grid(Longitud= seq(0,100,25),
                    Latitud=seq(0,200,length.out = 10),CO )

COR<- data.frame(Longitud = rep(observ$Longitud,2),Latitud = rep(observ$Latitud,2),profundidad = rep(c(5,10),each = 50),C.O= CO)
#grafica 3D
library(plotly)

## Warning: package 'plotly' was built under R version 4.0.3

## 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

fig <- plot_ly(COR, x = COR$Longitud, y = COR$Latitud, z = COR$profundidad,
               marker = list(color = COR$C.O, colorscale = c('#FFE1A1', '#683531'), showscale = TRUE))
fig <- fig %>% add_markers()
fig <- fig %>% layout(scene = list(xaxis = list(title = 'Longitud'),
                                   yaxis = list(title = 'Latitud'),
                                   zaxis = list(title = 'Profundidad')),
                      annotations = list(
                        x = 1.13,
                        y = 1.05,
                        text = 'CO',
                        xref = 'paper',
                        yref = 'paper',
                        showarrow = FALSE
                        ))
fig

## 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.

Comparación de las medias de CO de las dos capas a un nivel de cofianza del 95%

Hipotesis: \[H_o= \mu_{CO\ 5cm}=\mu_{CO\ 10cm}\\ H_a= \mu_{CO\ 5cm}\neq\mu_{CO\ 10cm}\]

# Se realiza la prueba T.student para datos pareados ya que son datos de dos profundidades de una misma muestra

pareado=t.test(COR$C.O~COR$profundidad, alternative= 't', paired = T )
h=ifelse(pareado$p.value<0.05, 'rechazo Ho', 'No rechazo Ho'); h

## [1] "rechazo Ho"

boxplot(COR$C.O~COR$profundidad)

La hipotesis nula se rechaza, por lo que hay diferencias en el carbono organico medido a profundidades diferentes de 5 y 10 cm. Siendo mayor el valor de la capa superior, al ser una capa con una posible mayor actividad organica y expuesta al reciclaje de material vegetal en la capa de suelo superficial.

Ejercicio 4

El siguiente diseño se corresponde con un factorial completo (3^2) en arreglo completamente al azar.

set.seed(2077)

D= expand.grid(F1= c(3.25,3.75,4.25), F2=c(4,5,6))
D=rbind(D,D)
set.seed(2075)
D= D[order(sample(1:18)),]
class(D)

## [1] "data.frame"

D$biomasa=sort.int(rnorm(18,3,0.3), partial = 9)
head(D)

##      F1 F2  biomasa
## 10 3.25  4 2.871097
## 13 3.25  5 2.623555
## 11 3.75  4 2.858383
## 16 3.25  6 2.267441
## 12 4.25  4 2.722977
## 6  4.25  5 2.763510

Modelo del diseño

\[ Y_{ijk}= \mu+\gamma_i+\delta_j+(\gamma\delta)_{ij}+\varepsilon_{ijk}\\i:1..3\\j:1..2\\k:1..2 \] \(\mu\) es la media global \(\gamma_i\) es el efecto del \(i\)-ésimo nivel del factor A \(\delta_i\) es el efecto del \(i\)-ésimo nivel del factor B \((\gamma\delta)_{ij}\) es el efecto causado por la interaccion del nivel \(i\) del factor A y del nivel \(j\) del factor B \(\varepsilon_{ijk}\) es el termino del error

2. Realice el Anova para este diseño y de ser necesario realice la prueba de comparaciones de medias para los efectos principales de F1: dosis de un insecticida que se sospecha tiene un efecto de disminución del crecimiento (biomasa) y F2: número de aplicaciones durante el desarrollo del cultivo.

dosis =as.factor(D$F1)
aplicaciones = as.factor(D$F2)
biomasa=as.vector(D$biomasa)

framebiom=data.frame(dosis,aplicaciones,biomasa); head(framebiom)

##   dosis aplicaciones  biomasa
## 1  3.25            4 2.871097
## 2  3.25            5 2.623555
## 3  3.75            4 2.858383
## 4  3.25            6 2.267441
## 5  4.25            4 2.722977
## 6  4.25            5 2.763510

aovbiom= aov(biomasa~dosis*aplicaciones); summary(aovbiom)

##                    Df Sum Sq Mean Sq F value Pr(>F)
## dosis               2 0.0359 0.01796   0.078  0.925
## aplicaciones        2 0.2385 0.11924   0.519  0.612
## dosis:aplicaciones  4 0.1952 0.04880   0.212  0.925
## Residuals           9 2.0675 0.22972

No se evidencia interaccion entre las dosis de insecticida con las aplicaciones ya que el p valor es mayor a 0.05, por lo tanto las hipotesis nulas para la interaccion y para el efecto de los dos factores no se rechazan ya que el efecto de los factores no es observable. Sin embargo el F value muestra un valor menor a 1, lo que puede significar que la variación en los tratamientos puede estar siendo alta intra tratamientos, lo cual puede ser señal de que el experimento no se está realizando de la manera adecuada, bloqueando y manejando los factores externos.

por lo tanto procedemos a realizar la prueba de Tukey de comparación de medias y entender quien es el causante de la variacion.

TukeyHSD(aovbiom)

##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = biomasa ~ dosis * aplicaciones)
## 
## $dosis
##                  diff        lwr       upr     p adj
## 3.75-3.25  0.10931116 -0.6632919 0.8819142 0.9183410
## 4.25-3.25  0.05066391 -0.7219392 0.8232670 0.9817273
## 4.25-3.75 -0.05864725 -0.8312503 0.7139558 0.9756088
## 
## $aplicaciones
##            diff        lwr       upr     p adj
## 5-4 -0.27210356 -1.0447066 0.5004995 0.6048618
## 6-4 -0.20000237 -0.9726054 0.5726007 0.7566730
## 6-5  0.07210119 -0.7005019 0.8447043 0.9634176
## 
## $`dosis:aplicaciones`
##                       diff       lwr      upr     p adj
## 3.75:4-3.25:4 -0.224872300 -2.120982 1.671238 0.9998467
## 4.25:4-3.25:4 -0.032385058 -1.928495 1.863725 1.0000000
## 3.25:5-3.25:4 -0.523345339 -2.419455 1.372765 0.9614880
## 3.75:5-3.25:4 -0.201238124 -2.097348 1.694872 0.9999328
## 4.25:5-3.25:4 -0.348984574 -2.245095 1.547125 0.9966334
## 3.25:6-3.25:4 -0.365993022 -2.262103 1.530117 0.9954037
## 3.75:6-3.25:4 -0.135294445 -2.031404 1.760816 0.9999968
## 4.25:6-3.25:4 -0.355976997 -2.252087 1.540133 0.9961641
## 4.25:4-3.75:4  0.192487242 -1.703623 2.088597 0.9999519
## 3.25:5-3.75:4 -0.298473039 -2.194583 1.597637 0.9988328
## 3.75:5-3.75:4  0.023634176 -1.872476 1.919744 1.0000000
## 4.25:5-3.75:4 -0.124112274 -2.020222 1.771998 0.9999983
## 3.25:6-3.75:4 -0.141120722 -2.037231 1.754989 0.9999955
## 3.75:6-3.75:4  0.089577855 -1.806532 1.985688 0.9999999
## 4.25:6-3.75:4 -0.131104697 -2.027215 1.765005 0.9999975
## 3.25:5-4.25:4 -0.490960280 -2.387070 1.405150 0.9728042
## 3.75:5-4.25:4 -0.168853066 -2.064963 1.727257 0.9999822
## 4.25:5-4.25:4 -0.316599515 -2.212709 1.579510 0.9982491
## 3.25:6-4.25:4 -0.333607964 -2.229718 1.562502 0.9975055
## 3.75:6-4.25:4 -0.102909387 -1.999019 1.793201 0.9999996
## 4.25:6-4.25:4 -0.323591939 -2.219702 1.572518 0.9979688
## 3.75:5-3.25:5  0.322107214 -1.574003 2.218217 0.9980311
## 4.25:5-3.25:5  0.174360765 -1.721749 2.070471 0.9999772
## 3.25:6-3.25:5  0.157352316 -1.738758 2.053462 0.9999896
## 3.75:6-3.25:5  0.388050893 -1.508059 2.284161 0.9933132
## 4.25:6-3.25:5  0.167368341 -1.728742 2.063478 0.9999833
## 4.25:5-3.75:5 -0.147746449 -2.043856 1.748364 0.9999936
## 3.25:6-3.75:5 -0.164754898 -2.060865 1.731355 0.9999852
## 3.75:6-3.75:5  0.065943679 -1.830166 1.962054 1.0000000
## 4.25:6-3.75:5 -0.154738873 -2.050849 1.741371 0.9999909
## 3.25:6-4.25:5 -0.017008448 -1.913118 1.879102 1.0000000
## 3.75:6-4.25:5  0.213690128 -1.682420 2.109800 0.9998948
## 4.25:6-4.25:5 -0.006992424 -1.903102 1.889118 1.0000000
## 3.75:6-3.25:6  0.230698577 -1.665411 2.126809 0.9998149
## 4.25:6-3.25:6  0.010016025 -1.886094 1.906126 1.0000000
## 4.25:6-3.75:6 -0.220682552 -2.116793 1.675427 0.9998665

Ya que todos los “p-valor ajustado” son mayores a 0.05%, se muestra que no ocurre ninguna interaccion entre los factores

#se ordenan los graficos que seran usados en el grafico
graf_biom = aggregate(x=list(meanbiom=framebiom$biomasa), by=list(Dosis = framebiom$dosis, Aplicacion=framebiom$aplicaciones),
  FUN=mean, na.rm=TRUE)
#grafico para observar interaccion

library(ggplot2)
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

graf_biom %>% 
  ggplot(aes(Dosis,meanbiom,color= Aplicacion))+
  geom_line(aes(group = Aplicacion))

Segun el grafico, si existe interaccion entre los factores y muestra al tratamiento 4 con la dosis 3.25 como el mejor tratamiento. sin embargo esta no es la conclusion mas fiable ya que como se mostro antes, existe una alta variabilidad entre los datos intra tratamiento, lo que afecta la respuesta del experimento. Sin embargo, seria posible llegar a un acuerdo al discutir que la dosis mas recomendada puede ser la de 3.75 para cualquier tratamiento y realizar el experimento de nuevo corrigiendo posibles efectos externos.

El investigador quiso colocar como covariable el contenido de arcilla(expansible) en el suelo utilizado en cada unidad experimental. Genere unos datos con la distribución uniforme cuya medida oscile entre 0.20 y 0.40 , ordene estas medidas en forma decreciente y meta dentro del análisis esta variable.

# Crea los datos de la covariable arcilla, con distribución uniforme.
Ar <- sort.default(runif(18, 0.20, 0.40), decreasing = TRUE) 

dataar = data.frame(D, Ar);head(dataar)

##      F1 F2  biomasa        Ar
## 10 3.25  4 2.871097 0.3987843
## 13 3.25  5 2.623555 0.3833869
## 11 3.75  4 2.858383 0.3803615
## 16 3.25  6 2.267441 0.3771692
## 12 4.25  4 2.722977 0.3758448
## 6  4.25  5 2.763510 0.3654590

Modelo 2

Modelo del diseño

\[ Y_{ijk}= \mu+\gamma_i+\delta_j+(\gamma\delta)_{ij}+\beta(x_{ijk}−\bar{x})+\varepsilon_{ijk}\\i:1..3\\j:1..2\\k:1..2 \]

Especifique nuevamente el modelo y realice el análisis de covarianza respectivo ¿se justifica el uso de la covariable? Construya nuevamente el gráfico de interacción y compare con el caso sin covariable (discuta el resultado). Revise en internet los supuestos que deben tener las covariables para ser utilizadas en el modelo. ¿Se está incumpliendo en nuestros datos alguno de los supuestos necesarios? Revise los supuestos sobre los residuales tanto del ANOVA como del ANCOVA ¿qué puede percibir? ¿recomendaría el uso de arcillas para minimizar el efecto sobre el contenido de biomasa que puede ocasionar el uso del insecticida?

Epecificacion del modelo

\(\mu\) es la media global \(\gamma_i\) es el efecto del \(i\)-ésimo nivel del factor A \(\delta_i\) es el efecto del \(i\)-ésimo nivel del factor B \((\gamma\delta)_{ij}\) es el efecto causado por la interaccion del nivel \(i\) del factor A y del nivel \(j\) del factor B \(\varepsilon_{ijk}\) es el termino del error \(\bar{x}\) es la media de los valores \(x_{ijk}\) $$ es el coeficiente de regresion, la relacion de \(Y_{ijk}\) con la covariable \(x_{ijk}\)

Ar <- sort.default(runif(18, 0.20, 0.40), decreasing = TRUE) 

dataar = data.frame(D, Ar);head(dataar)

##      F1 F2  biomasa        Ar
## 10 3.25  4 2.871097 0.3969166
## 13 3.25  5 2.623555 0.3919383
## 11 3.75  4 2.858383 0.3812244
## 16 3.25  6 2.267441 0.3722197
## 12 4.25  4 2.722977 0.3592836
## 6  4.25  5 2.763510 0.3331260

arc = as.factor(dataar$F1)
Aplicaciones_Arc = as.factor(dataar$F2)
biomasa_Arc = as.vector(dataar$biomasa)

model_ar= aov(data=dataar ,biomasa_Arc~Ar+arc+Aplicaciones_Arc); summary(model_ar)

##                  Df Sum Sq Mean Sq F value   Pr(>F)    
## Ar                1 1.7425  1.7425  64.021 3.75e-06 ***
## arc               2 0.0408  0.0204   0.749  0.49358    
## Aplicaciones_Arc  2 0.4271  0.2136   7.846  0.00662 ** 
## Residuals        12 0.3266  0.0272                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

El analisis de covarianza muestra que la arcilla puede influir en la respuesta de la biomasa, ya que se muestran diferencias en las medidas de biomasa.

\

Ejercicio 5

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.

Datos5 <- read_excel("C:/Users/Usuario/Desktop/Datos5.xlsx")
View(Datos5)

df2 = data.frame(Finca=as.factor(Datos5$Finca),Parcelas=Datos5$Parcelas,Variedad=Datos5$Variedad,Test=Datos5$Test,Rendimiento=Datos5$Rendimiento)
df2

##    Finca Parcelas Variedad Test Rendimiento
## 1      1        1        1    1        9.76
## 2      1        1        1    2        9.24
## 3      1        1        2    1       11.91
## 4      2        1        1    1        9.02
## 5      1        2        1    1       10.65
## 6      1        2        1    2        7.77
## 7      1        2        2    1       10.00
## 8      2        2        1    1       13.69
## 9      1        3        1    1        6.50
## 10     1        3        1    2        6.26
## 11     1        3        2    1        8.02
## 12     2        3        1    1        7.95
## 13     1        4        1    1        8.08
## 14     1        4        1    2        5.28
## 15     1        4        2    1        9.15
## 16     2        4        1    1        7.46
## 17     1        5        1    1        7.84
## 18     1        5        1    2        5.91
## 19     1        5        2    1        7.43
## 20     2        5        1    1        6.11
## 21     1        6        1    1        9.00
## 22     1        6        1    2        8.38
## 23     1        6        2    1        7.01
## 24     2        6        1    1        8.58
## 25     1        7        1    1       12.81
## 26     1        7        1    2       13.58
## 27     1        7        2    1       11.13
## 28     2        7        1    1       10.00
## 29     1        8        1    1       10.62
## 30     1        8        1    2       11.71
## 31     1        8        2    1       14.07
## 32     2        8        1    1       14.56
## 33     1        9        1    1        4.88
## 34     1        9        1    2        4.96
## 35     1        9        2    1        4.08
## 36     2        9        1    1        4.76
## 37     1       10        1    1        9.38
## 38     1       10        1    2        8.02
## 39     1       10        2    1        6.73
## 40     2       10        1    1        6.99
## 41     1       11        1    1        5.91
## 42     1       11        1    2        5.79
## 43     1       11        2    1        6.59
## 44     2       11        1    1        6.55
## 45     1       12        1    1        7.19
## 46     1       12        1    2        7.22
## 47     1       12        2    1        5.77
## 48     2       12        1    1        8.33
## 49     1       13        1    1        7.93
## 50     1       13        1    2        6.48
## 51     1       13        2    1        8.12
## 52     2       13        1    1        7.43
## 53     1       14        1    1        3.70
## 54     1       14        1    2        2.86
## 55     1       14        2    1        3.95
## 56     2       14        1    1        5.92
## 57     1       15        1    1        4.64
## 58     1       15        1    2        5.70
## 59     1       15        2    1        5.96
## 60     2       15        1    1        5.88
## 61     1       16        1    1        5.94
## 62     1       16        1    2        6.28
## 63     1       16        2    1        4.18
## 64     2       16        1    1        5.24
## 65     1       17        1    1        9.50
## 66     1       17        1    2        8.00
## 67     1       17        2    1       11.25
## 68     2       17        1    1       11.14
## 69     1       18        1    1       10.93
## 70     1       18        1    2       12.16
## 71     1       18        2    1        9.51
## 72     2       18        1    1       12.71
## 73     1       19        1    1       11.95
## 74     1       19        1    2       10.58
## 75     1       19        2    1       16.79
## 76     2       19        1    1       13.08
## 77     1       20        1    1        4.34
## 78     1       20        1    2        5.45
## 79     1       20        2    1        7.51
## 80     2       20        1    1        5.21

La tabla en la que se relacionan las Fincas, Parcelas, Variedads, Test y rendimientos fue creada en Excel y fue importada para posteriormente, crear el respectivo data frame.

A continuación se creó un diagrama de árbol para ilustrar el caso de estudio, correspondiente a un diseño anidad, factorial incompleto.

library(collapsibleTree)   

## Warning: package 'collapsibleTree' was built under R version 4.0.3

collapsibleTree(df2,hierarchy=c("Finca","Variedad","Test"))

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

mod1 = aov(Rendimiento ~ Parcelas + Parcelas:Finca + Parcelas:Finca:Variedad, data = df2)
summary(mod1)

##                         Df Sum Sq Mean Sq F value Pr(>F)
## Parcelas                 1    5.3   5.282   0.608  0.438
## Parcelas:Finca           1    3.1   3.077   0.354  0.554
## Parcelas:Finca:Variedad  1    5.2   5.179   0.596  0.443
## Residuals               76  660.7   8.693

Revisando los resultados del anova, se puede ver como todos los pvalores son mayores al 0.05 lo cual indicaría que ninguna variable tiene una varianza atribuíble bien sea a la parcela, la finca o la variedad.

Uso de la libreria lme4

Al cargar los datos desde Excel se presentó un inconveniente con el uso de la función lme4, ya que las variables fueron reconocidas como caracteres y la función trabaja con valores numéricos, por lo que fue necesario el uso de la función as.numeric para poder ejecutar el código.

df2$Rendimiento = as.numeric(df2$Rendimiento)
df2$Parcelas = as.numeric(df2$Parcelas)
df2$Finca = as.numeric(df2$Finca)
df2$Variedad = as.numeric(df2$Variedad)
df2$Test = as.numeric(df2$Test)

library(lme4)

## Warning: package 'lme4' was built under R version 4.0.3

## Loading required package: Matrix

## 
## Attaching package: 'lme4'

## The following object is masked from 'package:daewr':
## 
##     cake

mod2 <- lmer( Rendimiento ~ 1 + (1|Parcelas) + (1|Parcelas:Finca)+
+ (1|Parcelas:Finca:Variedad), data = df2)

## boundary (singular) fit: see ?isSingular

summary(mod2)

## Linear mixed model fit by REML ['lmerMod']
## Formula: Rendimiento ~ 1 + (1 | Parcelas) + (1 | Parcelas:Finca) + +(1 |  
##     Parcelas:Finca:Variedad)
##    Data: df2
## 
## REML criterion at convergence: 326
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.92753 -0.39932  0.00922  0.43797  1.65397 
## 
## Random effects:
##  Groups                  Name        Variance Std.Dev.
##  Parcelas:Finca:Variedad (Intercept) 1.2305   1.1093  
##  Parcelas:Finca          (Intercept) 0.0000   0.0000  
##  Parcelas                (Intercept) 7.0127   2.6481  
##  Residual                            0.8795   0.9378  
## Number of obs: 80, groups:  
## Parcelas:Finca:Variedad, 60; Parcelas:Finca, 40; Parcelas, 20
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)   8.2369     0.6188   13.31
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular

Calculo de porcentaje de variación en parcelas

Este calculo se obtiene mediante la siguente fórmula:

variacionp = (100*7.0127)/(7.0127+1.2305+0.8795)
variacionp

## [1] 76.87088

Este resultado demuestra que el 76.8 % de la variación se presenta entre parcelas.

Ejercicio 8

¿Cuándo debería considerar un diseño en parcelas divididas? (Christine M. Anderson – Cook) (http://207.67.83.164/quality-progress/2007/10/laboratory/when-should-you-consider-a- split-plot-design.html) En este artículo la autora explica a grandes rasgos las consideraciones a tener en cuenta y las ventajas que representa el uso de diseños en parcelas divididas en las investigaciones. En un primer vistazo al origen de este diseño, la investigadora plantea la procedencia del nombre del diseño e infiere que al ser un diseño originario del ámbito investigativo agrícola, mucha de su terminología está asociada con parcelas de tierra, además de esto esboza los principales conceptos básicos asociados como factores de parcela completa (aquellos factores que son difíciles de cambiar o manejar, como por ejemplo el riego) y factores de subparcela (los factores más fáciles de cambiar o controlar, como por ejemplo la fertilización), cabe aclarar que estos factores podrían variar de acuerdo con la situación del investigador pero que es necesario tenerlos en cuenta a la hora de usar un diseño en parcelas divididas. Posteriormente, el texto trae a colación el hecho de que es posible confundirse e interpretar este tipo de diseños como diseños completamente aleatorizados lo que podría resultar en análisis equívocos por parte del experimentador y reconoce que si se hace un uso adecuado teniendo todas las consideraciones del diseño, este puede representar una gran ventaja frente a otros métodos, debido a que al tener los factores separados pueden obtenerse dos resultados que pueden ser interpretados de manera dividida.

Sobre el reconocimiento de la unidad experimental adecuada en los estudios con animales en las ciencias lácteas. (https://www.sciencedirect.com/science/article/pii/S002203021630621X) (https://online.stat.psu.edu/stat502/lesson/6/6.1-0) Esta publicación evidencia los conflictos que se dan dentro de un principio inherente a la investigación científica, el hecho de que sea reproducible. A grandes rasgos y yéndose al campo del diseño experimental los autores expresan una inconformidad respecto a la forma en la que pueden llegarse a usar de manera indiscriminada los conceptos de unidad experimental (también llamada unidad de replicación) y la unidad de observación (conocida también como unidad de muestreo) y el problema que representa para la inferencia estadística este uso confuso de conceptos en la industria láctea. A lo largo de la lectura, se hace necesario un replanteamiento de lo que se viene haciendo en torno a la planificación y ejecución de los diseños experimentales, desde la jerarquización de los factores hasta el posterior análisis, haciendo uso del concepto de jerarquización en el diseño; se nos cuenta cómo este concepto hace referencia a la configuración de datos que presentan una estructura correlacionada impuesta por el mismo diseño y además se enfatiza en la importancia de la replicabilidad como requisito fundamental para la inferencia experimental. Todo lo anterior con la preocupación y el afán de realizar diseños experimentales dentro de la industria láctea que puedan ser replicados de acuerdo a las necesidades de los investigadores y los usuarios planteando la siguiente pregunta ¿Cuál es el efecto de un tratamiento en un resultado específico de interés?

Fundamentos del diseño experimental: pautas para diseñar experimentos exitosos (Michael D. Casler) (https://acsess.onlinelibrary.wiley.com/doi/full/10.2134/agronj2013.0114) En esta publicación el autor expone una problemática que se presenta en las áreas biológicas en cuanto a la planeación y la realización de diseños experimentales, que consiste principalmente en el “amaño” que tienen algunos investigadores a la hora de elegir un tipo de diseño, ya que estos se encuentran casados con los mismos modelos que han utilizados siempre, y como el mismo autor lo expresa, usan estos diseños como si fueran una receta, lo cual suscita una creciente preocupación en su área debido a que en las investigaciones biológicas pueden encontrarse múltiples fuentes de variación que inciden en los resultados finales. Para solucionar este inconveniente, el artículo realiza una descripción detallada sobre los que se consideran los cuatro pilares del diseño experimental: 1. Replicación (importante para las estimaciones del error experimental y la precisión de los experimentos), 2. Aleatorización (clave para no obtener estimaciones sesgadas en las medias de los tratamientos y en los errores experimentales), 3. Bloqueo (oportuno para contrarrestar posibles fuentes de variación en los resultados) y 4. Tamaño de las unidades experimentales (uso de la ley de Smith y uno de los 4 pilares menos estudiados); poniendo de manifiesto algunas de sus virtudes y consideraciones a tener en cuenta en un ámbito en el que la planeación de los diseños muchas veces se ve paleada por los recursos necesarios para su realización y por ende, el poco lugar que dejan estás condiciones económicas altas para fallar en los ensayos. Además de esto, el autor invita a los nuevos investigadores que empiezan a sumergirse en el mundo de la ciencia y el diseño experimental a probar nuevos diseños y retroalimentarse de sus hallazgos para enriquecer todo lo que se conoce sobre el diseño experimental hoy día.