1. Prueba t-Dos muestras independientes:

Se midió la conductancia estomática (gs: mol/m2s) en dos cultivares de papa diploide (Colombia y Ocarina) bajo una condición de déficit de riego. Parte de los datos se muestran en la siguiente tabla:

Nota: En mi caso la U=7

Colombia (0.45, 0.41, 0.4U, 0.46, 0.39, 0.44, 0.48, 0.42, 0.44, 0.48, 0.50, 0.47, 0.44, 0.52)

Ocarina (0.28, 0.25, 0.32, 0.34, 0.36, 0.40, 0.3U, 0.36, 0.39, 0.41, 0.37, 0.42, 0.41)

Determinar al 95% de nivel de confianza si las dos medias obtenidas para los cultivares son estadísticamente iguales. Utilice la información del artículo mostrado en clase para decidir si las varianzas pueden considerarse iguales o no. (más información del tema en: https://revistas.unal.edu.co/index.php/agrocol/article/view/78982/74726)

Desarrollo actividad:

Según el artículo los resultados obtenidos en la conductancia estomatica las diferencias que se encontraron fueron entre plantas con estres hidrico respecto a los que no tenian estres hidrico, por dichas razones las varianzas serian iguales, esto porque nuestros datos solo representan ha plantas con estres hidrico y ambos cultivares fueron medidos con equipos similares y en los mismos tiempos, en conclusión el artículo define que las varianzas son iguales; igualmente para tener una mayor seguridad de que las varianzas sean iguales, adelante se calcula el test de Shapiro-Wilk para igualdad de varianzas, ya que igualmetne para aplicar la prueba T, se requiere que nuestros datos tengan una varianza igual.

Tambien antes de continuar para este ejercicio como se debe realizar una prueb de T, las hipotesis serian las siguentes:

\[H_0:\mu_1=\mu_2\] \[H_1:\mu_1\neq\mu_2\]

Colombia=c(0.45,0.41,0.47,0.46,0.39,0.44,0.48,0.42,0.44,0.48,0.50,0.47,0.44,0.52)  
Ocarina=c(0.28,0.25,0.32,0.34,0.36,0.40,0.37,0.36,0.39,0.41,0.37,0.42,0.41)
shapiro.test(Colombia)
## 
##  Shapiro-Wilk normality test
## 
## data:  Colombia
## W = 0.98666, p-value = 0.9971
shapiro.test(Ocarina)
## 
##  Shapiro-Wilk normality test
## 
## data:  Ocarina
## W = 0.90993, p-value = 0.183

Como podemos observar, el p-valor de ambos cultivares (Colombia=0.9971 y Ocarina= 0.183) se sitúan por encima de 0.05. Esto significa que aceptamos la hipótesis nula del test de shapiro-wilk y consideramos que ambas se distribuyen siguiendo una distribución normal.

colombia=c(0.45,0.41,0.47,0.46,0.39,0.44,0.48,0.42,0.44,0.48,0.50,0.47,0.44,0.52)  
ocarina=c(0.28,0.25,0.32,0.34,0.36,0.40,0.37,0.36,0.39,0.41,0.37,0.42,0.41)
grup = gl(2,14,27,c("Colombia","Ocarina"))
df = data.frame(c(Colombia,Ocarina), grup)
colnames(df) = c("rto", "variedad")
df
##     rto variedad
## 1  0.45 Colombia
## 2  0.41 Colombia
## 3  0.47 Colombia
## 4  0.46 Colombia
## 5  0.39 Colombia
## 6  0.44 Colombia
## 7  0.48 Colombia
## 8  0.42 Colombia
## 9  0.44 Colombia
## 10 0.48 Colombia
## 11 0.50 Colombia
## 12 0.47 Colombia
## 13 0.44 Colombia
## 14 0.52 Colombia
## 15 0.28  Ocarina
## 16 0.25  Ocarina
## 17 0.32  Ocarina
## 18 0.34  Ocarina
## 19 0.36  Ocarina
## 20 0.40  Ocarina
## 21 0.37  Ocarina
## 22 0.36  Ocarina
## 23 0.39  Ocarina
## 24 0.41  Ocarina
## 25 0.37  Ocarina
## 26 0.42  Ocarina
## 27 0.41  Ocarina
library(ggplot2)
ggplot(df, aes(x = rto, fill = variedad))+
  geom_density()

Igualmente al realizar un grafico de densidad encontramos que los dos cultivares presentan una forma de campana de gauss, donde confirma que presenta normalidad.

t.test(colombia,ocarina, alternative = "t", mu = 0, 
       conf.level = 0.95)
## 
##  Welch Two Sample t-test
## 
## data:  colombia and ocarina
## t = 5.5539, df = 21.041, p-value = 1.63e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.05943233 0.13056767
## sample estimates:
## mean of x mean of y 
##     0.455     0.360

Finalmente al realizar la prueba T par dos muestras independientes, se rechaza la hipótesis; como podemos observar el P-valor es menor a nuestros nivel de significancia del 0.05, por lo cual las dos medias no son iguales para los dos cultivares, Para este ejercicio el valor T encontrado es 5.5539 y que nuestro intervalo de confianza al 95% va desde 0.0594 a 0.1305, donde no se incluye el cero y esto también demuestra que las medias no son iguales.

  1. Prueba t-Dos muestras dependientes/pareadas

Se propuso un plan de fertilización en papa criolla tal como se muestra a continuación:

N(100 kg/ha); P2O2(50 kg/ha); K2O(100 kg/ha); CaO(30 kg/ha); MgO(24 kg/ha); S(20 kg/ha); Fe(0 kg/ha); Mn(0 kg/ha); Cu(0.30 kg/ha); Zn(2.0 kg/ha); B(1.2 Kg/ha)

y se midió a los 45 y 77 días después de la siembra el peso de tubérculos (Kg/ha) más raíces encontrando los siguientes datos:

45dds (69, 66, 72, 68, 65, 66, 67, 68, 69, 6U, 66, 68, 64, 67, 60, 68)

77dds (873, 850, 832, 834, 843, 840, 8U5, 790, 905, 910, 920, 840, 832, 800, 759, 812)

Determinar al 95% de nivel de confianza si se incrementó la medida de rendimiento en las dos evaluaciones registradas. Haga una representación gráfica para ilustrar el comportamiento de ambas medidas. Calcule el cambio relativo porcentual promedio entre ambos tiempos de evaluación. Calcule el coeficiente de correlación de Pearson entre ambas medidas. Explique sus resultados. (más información del tema en: http://www.scielo.org.co/pdf/agc/v32n1/v32n1a09.pdf )

PTR=data.frame(med1=c(69,66,72,68,65,66,67,68,69,67,66,68,64,67,60,68), 
               med2=c(873,850,832,834,843,840,875,790,905,910,920,840,832,800,759,812))
colnames(PTR) = c("45dds", "77dds")
PTR
##    45dds 77dds
## 1     69   873
## 2     66   850
## 3     72   832
## 4     68   834
## 5     65   843
## 6     66   840
## 7     67   875
## 8     68   790
## 9     69   905
## 10    67   910
## 11    66   920
## 12    68   840
## 13    64   832
## 14    67   800
## 15    60   759
## 16    68   812
shapiro.test(PTR$`45dds`)
## 
##  Shapiro-Wilk normality test
## 
## data:  PTR$`45dds`
## W = 0.91678, p-value = 0.1497
shapiro.test(PTR$`77dds`)
## 
##  Shapiro-Wilk normality test
## 
## data:  PTR$`77dds`
## W = 0.9645, p-value = 0.7435

El p-valor de ambos días donde se midio (45dds=0.1497 y 77dds= 0.7435) se sitúa por encima de 0.05. Esto significa que aceptamos la hipótesis nula del test de shapiro-wilk y consideramos que ambas se distribuyen siguiendo una distribución normal.

t.test(PTR$`45dds`,PTR$`77dds`, paired = TRUE, alternative = "t", mu=0, conf.level = 0.95) 
## 
##  Paired t-test
## 
## data:  PTR$`45dds` and PTR$`77dds`
## t = -71.814, df = 15, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -800.8982 -754.7268
## sample estimates:
## mean of the differences 
##               -777.8125

Con un p-value = 2.2e-16 menor de 0.05 podemos rechazar la hipótesis nula \(H_0\) de igualdad de medias. Podemos concluir que existen diferencias entre la media de los valores en la primera medición y la de los valores en la segunda medición.

med1=c(69,66,72,68,65,66,67,68,69,67,66,68,64,67,60,68)
med2=c(873,850,832,834,843,840,875,790,905,910,920,840,832,800,759,812)
Renbox <- data.frame(med1, med2)
summary(Renbox)
##       med1            med2      
##  Min.   :60.00   Min.   :759.0  
##  1st Qu.:66.00   1st Qu.:827.0  
##  Median :67.00   Median :840.0  
##  Mean   :66.88   Mean   :844.7  
##  3rd Qu.:68.00   3rd Qu.:873.5  
##  Max.   :72.00   Max.   :920.0
boxplot(Renbox, border=3, col=8, horizontal=T, xlab="(kg/ha)", ylab="Dias")

El grafico de boxplot presenta las dos mediciones realizados, observandose obiamente un crecimiento en el rendimiento bastante alto alos 77 dds respecto a los 45dds

med1=c(69,66,72,68,65,66,67,68,69,67,66,68,64,67,60,68)
med2=c(873,850,832,834,843,840,875,790,905,910,920,840,832,800,759,812)
grup2=gl(2,16,32, c("45dds", "75dds"))
PTR2=data.frame(c(med1,med2), grup2) 
colnames(PTR2) = c("datos", "tomadatos")
PTR2
##    datos tomadatos
## 1     69     45dds
## 2     66     45dds
## 3     72     45dds
## 4     68     45dds
## 5     65     45dds
## 6     66     45dds
## 7     67     45dds
## 8     68     45dds
## 9     69     45dds
## 10    67     45dds
## 11    66     45dds
## 12    68     45dds
## 13    64     45dds
## 14    67     45dds
## 15    60     45dds
## 16    68     45dds
## 17   873     75dds
## 18   850     75dds
## 19   832     75dds
## 20   834     75dds
## 21   843     75dds
## 22   840     75dds
## 23   875     75dds
## 24   790     75dds
## 25   905     75dds
## 26   910     75dds
## 27   920     75dds
## 28   840     75dds
## 29   832     75dds
## 30   800     75dds
## 31   759     75dds
## 32   812     75dds
(medias = tapply(PTR2$datos,PTR2$tomadatos, mean))
##    45dds    75dds 
##  66.8750 844.6875
Cambiorelativo = (844.6875-66.8750)/(844.6875)
Cambiorelativo
## [1] 0.9208287

Al calcular el cambio relativo se observa que hubo un cambio alrededor del 92%, entonces desde la toma de datos a los 45dds hasta los 77dds, se encuentra un incremento de rendimiento bastante aceptable con el plan de fertilización que se tiene para la papa criolla.

library(PerformanceAnalytics)
## Loading required package: xts
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
## 
##     legend
med1=c(69,66,72,68,65,66,67,68,69,67,66,68,64,67,60,68)
med2=c(873,850,832,834,843,840,875,790,905,910,920,840,832,800,759,812)
dat1 <- data.frame(med1, med2)
chart.Correlation(dat1)

En este caso el coficiente de correlación de Pearson es o.333 lo que indica que es mayor a cero (r > 0) Entonces el 45dds y 77dds tienen correlación positiva: Ambas variables se correlacionan en un sentido directo.

  1. Prueba de Wilcoxon de la suma de rangos-Dos muestras independientes

Se está evaluando la calidad de frito mediante la textura de las hojuelas de papa criolla en dos tipos de aceite (palma y maíz) utilizado para freír en condiciones controladas de tiempo y temperatura. Al final se recolectaron las hojuelas y se evaluó en una escala diagramática la calidad de frito (escala de 1 a 5, desde (1) no crujiente hasta (5) bastante crujientes). Los datos se muestran a continuación:

Palma (3, 4, 3, 4, 4, 3, 3, 4, 4, 3, 4, 4, 2, 4, 3, 4, 3, 3, 3, 4, 4)

Maíz (3, 4, 4, 4, 4, 4, 3, 4, 3, 4, 4, 4, 4, 3, 4, 4, 4, 3, 3, 4, 3)

Determinar al 95% de nivel de confianza si existen diferencias estadísticas en las medianas de la textura para los dos tipos de aceite. Haga una representación gráfica para ilustrar el comportamiento de ambas medidas. Explique sus resultados.
más información relacionada con el tema: https://core.ac.uk/download/pdf/61543501.pdf

\[H_0:La mediana de las diferencias de cada par de datos es cero. Mediana(diferencias)=0.\] \[H_a:La mediana de las diferencias entre cada par de datos es diferente de cero. Mediana(diferencias)≠0.\]

Palma=c(3,  4,  3,  4,  4,  3,  3,  4,  4,  3,  4,  4,  2,  4,  3,  4,  3,  3,  3,  4, 4)
Maiz=c(3,   4,  4,  4,  4,  4,  3,  4,  3,  4,  4,  4,  4,  3,  4,  4,  4,  3,  3,  4,  3)
diferencias <-c(Palma-Maiz)
rbind(Palma, Maiz,diferencias)
##             [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
## Palma          3    4    3    4    4    3    3    4    4     3     4     4
## Maiz           3    4    4    4    4    4    3    4    3     4     4     4
## diferencias    0    0   -1    0    0   -1    0    0    1    -1     0     0
##             [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
## Palma           2     4     3     4     3     3     3     4     4
## Maiz            4     3     4     4     4     3     3     4     3
## diferencias    -2     1    -1     0    -1     0     0     0     1
Evaluacion<-data.frame(Palma=Palma, Maiz=Maiz)
Evaluacion
##    Palma Maiz
## 1      3    3
## 2      4    4
## 3      3    4
## 4      4    4
## 5      4    4
## 6      3    4
## 7      3    3
## 8      4    4
## 9      4    3
## 10     3    4
## 11     4    4
## 12     4    4
## 13     2    4
## 14     4    3
## 15     3    4
## 16     4    4
## 17     3    4
## 18     3    3
## 19     3    3
## 20     4    4
## 21     4    3
par(mfrow=c(1,1))
hist(Palma,col = "green", xlab = "Calidad del frito", main = "Aceite de Palma")

hist(Palma,col = "blue", xlab = "Calidad del frito", main = "Aceite de Maíz")

wilcox.test(Palma, Maiz, mu=0, alternative = "two.sided", paired = FALSE, conf.level = 0.95)
## Warning in wilcox.test.default(Palma, Maiz, mu = 0, alternative = "two.sided", :
## cannot compute exact p-value with ties
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  Palma and Maiz
## W = 185.5, p-value = 0.3111
## alternative hypothesis: true location shift is not equal to 0

Como el P-valor es mayor al 0.05, no existen diferencias estadisticas en la fritura con los dos tipos de aceites, esto también se puede observar de igual manera en las graficas de histogramas donde se observa un comportamiento igual, entonces la fritura de la papa con los dos aceites se obtendra el mismo resultado.

  1. Prueba F- AOV-FSCA-B

Un laboratorio emplea 3 métodos (Bray,Olsen; Mehlich-3) para determinar el contenido de fósforo en suelos. Surge la pregunta: “¿Difieren las medias determinaciones entre los métodos de análisis?” Para responder a esta pregunta, se analizaron 10 muestras de una misma unidad de manejo agronómico. Los resultados se muestran (en mg de fósforo/Kg de suelo) a continuación:

Bray Olsen Mehlich-3 7.1 7.4 6.0 6.6 9.1 8.7 6.8 6.6 6.2 6.6 7.1 7.6 6.6 6.7 6.4 6.7 7.8 7.8 6.7 6.8 6.6 6.4 7.3 7.8 6.8 7.0 6.9 6.2 7.6 7.2 6.7 6.U 6.6 6.6 7.8 7.2 6.9 6.7 6.4 6.3 7.4 7.U 6.8 7.1 6.4 5.8 7.3 7.3 6.7 6.6 6.5 6.5 7.3 7.5 6.6 6.6 6.3 6.2 7.1 7.1 6.9 6.9 6.4 6.4 8.0 7.1

Realice un análisis descriptivo de los datos (tablas de medias y coeficientes de variación). Construya gráficos apropiados para visualizar todos los métodos a la vez.

Bray=c(7.1,6.8,6.6,6.7,6.8,6.7,6.9,6.8,6.7,6.6, 6.9,7.4,6.6,6.7,6.8,7.0,6.7,6.7,7.1,6.6,6.6,6.9)
Olsen=c(6.0,6.2,6.4,6.6,6.9,6.6,6.4,6.4,6.5,6.3,6.4,6.6,6.6,6.7,6.4,6.2,6.6,6.3,5.8,6.5,6.2,6.4)
Mehlich=c(9.1,7.1,7.8,7.3,7.6,7.8,7.4,7.3,7.3,7.1,8.0,8.7,7.6,7.8,7.8,7.2,7.2,7.7,7.3,7.5,7.1,7.1)

grup = gl(3,22,66,c("Bray","Olsen","Mehlich"))
df = data.frame(c(Bray,Olsen,Mehlich), grup)
colnames(df) = c("dato", "metodo")
df
##    dato  metodo
## 1   7.1    Bray
## 2   6.8    Bray
## 3   6.6    Bray
## 4   6.7    Bray
## 5   6.8    Bray
## 6   6.7    Bray
## 7   6.9    Bray
## 8   6.8    Bray
## 9   6.7    Bray
## 10  6.6    Bray
## 11  6.9    Bray
## 12  7.4    Bray
## 13  6.6    Bray
## 14  6.7    Bray
## 15  6.8    Bray
## 16  7.0    Bray
## 17  6.7    Bray
## 18  6.7    Bray
## 19  7.1    Bray
## 20  6.6    Bray
## 21  6.6    Bray
## 22  6.9    Bray
## 23  6.0   Olsen
## 24  6.2   Olsen
## 25  6.4   Olsen
## 26  6.6   Olsen
## 27  6.9   Olsen
## 28  6.6   Olsen
## 29  6.4   Olsen
## 30  6.4   Olsen
## 31  6.5   Olsen
## 32  6.3   Olsen
## 33  6.4   Olsen
## 34  6.6   Olsen
## 35  6.6   Olsen
## 36  6.7   Olsen
## 37  6.4   Olsen
## 38  6.2   Olsen
## 39  6.6   Olsen
## 40  6.3   Olsen
## 41  5.8   Olsen
## 42  6.5   Olsen
## 43  6.2   Olsen
## 44  6.4   Olsen
## 45  9.1 Mehlich
## 46  7.1 Mehlich
## 47  7.8 Mehlich
## 48  7.3 Mehlich
## 49  7.6 Mehlich
## 50  7.8 Mehlich
## 51  7.4 Mehlich
## 52  7.3 Mehlich
## 53  7.3 Mehlich
## 54  7.1 Mehlich
## 55  8.0 Mehlich
## 56  8.7 Mehlich
## 57  7.6 Mehlich
## 58  7.8 Mehlich
## 59  7.8 Mehlich
## 60  7.2 Mehlich
## 61  7.2 Mehlich
## 62  7.7 Mehlich
## 63  7.3 Mehlich
## 64  7.5 Mehlich
## 65  7.1 Mehlich
## 66  7.1 Mehlich
(medias = tapply(df$dato,df$metodo, mean))
##     Bray    Olsen  Mehlich 
## 6.804545 6.409091 7.581818
(varianzas = tapply(df$dato,df$metodo, var))
##       Bray      Olsen    Mehlich 
## 0.04140693 0.05800866 0.26251082
(desv = tapply(df$dato,df$metodo, sd))
##      Bray     Olsen   Mehlich 
## 0.2034869 0.2408499 0.5123581
(Coenfvariació= tapply(df$dato,df$metodo,sd)/tapply(df$dato,df$metodo, mean))*100
##     Bray    Olsen  Mehlich 
## 2.990456 3.757941 6.757721
boxplot(df$dato~df$metodo)
points(medias, pch = 19, col = 'red')

Las medias entre los métodos varian o son diferentes, con esto se puede intuir que un método para medir el nivel de fosforo pueder ser mejor, aunque se debe tener en cuenta que en nuestro ejercicio las medias no son tan diferentes. En cuanto al coeficiente de varaición, en los tres metodos que se usan para analizar el contenido de fosforo presentan una varavilidad muy baja, esto porque presentan un CV entre 0 y 0.1.

Realice un análisis de Varianza para contrastar la hipótesis:

\[Ho:\tau_B=\tau_O=\tau_M=0\] \[H_a:H_0 es falsa\]

Presente el valor p tal como se presenta en la tabla de salida de R. Interprete el resultado.

mod1 = aov(df$dato~df$metodo)
summary(mod1)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## df$metodo    2  15.66   7.831   64.91 4.97e-16 ***
## Residuals   63   7.60   0.121                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res1 = mod1$residuals

El p-valor es significativo, al ser significativo se rechaza la hipotesis nula, entonces al menos un par de medias son diferentes

norm = shapiro.test(res1)
norm
## 
##  Shapiro-Wilk normality test
## 
## data:  res1
## W = 0.86946, p-value = 5.098e-06

No hay normalidad en lo datos, esto porque el p-valor es muy pequeña, por lo tanto no tienen una distribución normal

hom = bartlett.test(res1, df$metodo)
hom
## 
##  Bartlett test of homogeneity of variances
## 
## data:  res1 and df$metodo
## Bartlett's K-squared = 21.061, df = 2, p-value = 2.67e-05

En este caso tambien no se comple la homocedasticidad de los datos debido a un p-valor muy cerca de cero.

l_P = log(x = df$dato,base = 10)
mod1l = aov(l_P ~ df$metodo)
mod1l
## Call:
##    aov(formula = l_P ~ df$metodo)
## 
## Terms:
##                  df$metodo  Residuals
## Sum of Squares  0.05912181 0.02565021
## Deg. of Freedom          2         63
## 
## Residual standard error: 0.02017787
## Estimated effects may be unbalanced
norml = shapiro.test(mod1l$residuals)
norml
## 
##  Shapiro-Wilk normality test
## 
## data:  mod1l$residuals
## W = 0.91247, p-value = 0.000195

En esta trasformación se rechaza normalidad

boxplot(mod1$residuals, main = 'Boxplor Residuales')

En el anterior gráfico hay dos datos atipicos Hipotesis: primero que no es un dato atipico y alternativa es atipica y buscamos el dato atipico en los residuales

\[H_0: Dato_{\ max}\] \[H_a: Dato_{\ max}\]

library(outliers)
dat = grubbs.test(x = mod1$residuals,two.sided = TRUE)
dat
## 
##  Grubbs test for one outlier
## 
## data:  mod1$residuals
## G.45 = 4.43977, U = 0.69208, p-value = 8.785e-05
## alternative hypothesis: highest value 1.51818181818182 is an outlier

P valor < 0.05, rechazo \(H_0\), el valor 1.52 es un dato atipico.

which.max(mod1$residuals)
## 45 
## 45
(medias=tapply(df$dato,df$metodo, mean))
##     Bray    Olsen  Mehlich 
## 6.804545 6.409091 7.581818
round(res1[45:57],1)
##   45   46   47   48   49   50   51   52   53   54   55   56   57 
##  1.5 -0.5  0.2 -0.3  0.0  0.2 -0.2 -0.3 -0.3 -0.5  0.4  1.1  0.0
df$dato[45]=medias[3]
df$dato[56]=medias[3]
pander::pander(df)
dato metodo
7.1 Bray
6.8 Bray
6.6 Bray
6.7 Bray
6.8 Bray
6.7 Bray
6.9 Bray
6.8 Bray
6.7 Bray
6.6 Bray
6.9 Bray
7.4 Bray
6.6 Bray
6.7 Bray
6.8 Bray
7 Bray
6.7 Bray
6.7 Bray
7.1 Bray
6.6 Bray
6.6 Bray
6.9 Bray
6 Olsen
6.2 Olsen
6.4 Olsen
6.6 Olsen
6.9 Olsen
6.6 Olsen
6.4 Olsen
6.4 Olsen
6.5 Olsen
6.3 Olsen
6.4 Olsen
6.6 Olsen
6.6 Olsen
6.7 Olsen
6.4 Olsen
6.2 Olsen
6.6 Olsen
6.3 Olsen
5.8 Olsen
6.5 Olsen
6.2 Olsen
6.4 Olsen
7.582 Mehlich
7.1 Mehlich
7.8 Mehlich
7.3 Mehlich
7.6 Mehlich
7.8 Mehlich
7.4 Mehlich
7.3 Mehlich
7.3 Mehlich
7.1 Mehlich
8 Mehlich
7.582 Mehlich
7.6 Mehlich
7.8 Mehlich
7.8 Mehlich
7.2 Mehlich
7.2 Mehlich
7.7 Mehlich
7.3 Mehlich
7.5 Mehlich
7.1 Mehlich
7.1 Mehlich 
mod1c = aov(df$dato~df$metodo)
summary(mod1c)
##             Df Sum Sq Mean Sq F value Pr(>F)    
## df$metodo    2 12.446   6.223   105.1 <2e-16 ***
## Residuals   63  3.729   0.059                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res2 = mod1c$residuals
normc = shapiro.test(mod1c$residuals)
normc
## 
##  Shapiro-Wilk normality test
## 
## data:  mod1c$residuals
## W = 0.98471, p-value = 0.593

Al realizar los cambios de los datos atipicos finalmente conseguimos que tengan normalidad, entonces podemos concluir que un P.valor < 0.05, existen diferencias estadisticas en el metodo usado para determinar el contenido de fosforo. Se rechaza la hipótesis nula de que las medias de los 5 tratamientos no son iguales, esto significa que al menos una media es diferente.

tk = TukeyHSD(mod1c)
tk
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = df$dato ~ df$metodo)
## 
## $`df$metodo`
##                     diff        lwr        upr   p adj
## Olsen-Bray    -0.3954545 -0.5715379 -0.2193712 3.3e-06
## Mehlich-Bray   0.6574380  0.4813546  0.8335214 0.0e+00
## Mehlich-Olsen  1.0528926  0.8768092  1.2289759 0.0e+00

Podemos observar que los P-valor son <0.05, lo que indica que todos los tratamientos comparados son diferentes

(varianza=tapply(df$dato,df$metodo, var))
##       Bray      Olsen    Mehlich 
## 0.04140693 0.05800866 0.07817109
data<-varianza
plot(medias, varianza, col="blue", pch=16, cex=3)

data.aov<-aov(dato~metodo, df)
plot(data.aov)

El gráfico Q-Q ayuda a comparar gráficamente dos distribuciones, es decir permite comparar la distribución de los residuos con la distribución normal teórica. Para este los residuos tienen una distribución normal porque se observa que siguen aproximadamente la línea recta diagonal en el gráfico Q-Q normal.

El siguente gráfico que es constant leverage es util para ayudar a identificar puntos de datos influyentes en un modelo, es decir que los valores atípicos pueden ser influyentes, aunque no necesariamente tienen que hacerlo y algunos puntos dentro de un rango normal en un modelo podrían ser muy influyentes.

anova(data.aov)
## Analysis of Variance Table
## 
## Response: dato
##           Df  Sum Sq Mean Sq F value    Pr(>F)    
## metodo     2 12.4461  6.2230  105.13 < 2.2e-16 ***
## Residuals 63  3.7293  0.0592                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

El cuadrado medio del error de nuestro análisis de varianza es 3.6977, mientras que el del ajustado es 3.698 y como el cuadrado medio del error mide la variabilidad dentro de cada nivel, es decir, la variabilidad que no es debida a variable factor, entonces en la nueva análisis de varianza y en la ajsutando respecto a la variabilidad dentro de cada nivel no se observan varaibilidad.

boxplot(df$dato~df$metodo)
points(medias, pch = 19, col = 'red')

oneway.test(df$dato~df$metodo, var.equal = F)
## 
##  One-way analysis of means (not assuming equal variances)
## 
## data:  df$dato and df$metodo
## F = 88.532, num df = 2.000, denom df = 41.309, p-value = 1.158e-15

El test oneway muestra un mejor estadistico respecto con el modelo original.

matrix(sample(df$metodo), ncol = 6)
##       [,1]      [,2]      [,3]      [,4]      [,5]      [,6]     
##  [1,] "Bray"    "Bray"    "Mehlich" "Mehlich" "Bray"    "Mehlich"
##  [2,] "Bray"    "Bray"    "Bray"    "Bray"    "Bray"    "Olsen"  
##  [3,] "Bray"    "Bray"    "Olsen"   "Olsen"   "Mehlich" "Mehlich"
##  [4,] "Mehlich" "Bray"    "Bray"    "Olsen"   "Bray"    "Olsen"  
##  [5,] "Mehlich" "Mehlich" "Mehlich" "Bray"    "Bray"    "Mehlich"
##  [6,] "Mehlich" "Olsen"   "Bray"    "Olsen"   "Olsen"   "Olsen"  
##  [7,] "Olsen"   "Mehlich" "Olsen"   "Olsen"   "Mehlich" "Olsen"  
##  [8,] "Olsen"   "Bray"    "Mehlich" "Olsen"   "Bray"    "Olsen"  
##  [9,] "Bray"    "Mehlich" "Olsen"   "Olsen"   "Mehlich" "Mehlich"
## [10,] "Mehlich" "Mehlich" "Olsen"   "Olsen"   "Bray"    "Mehlich"
## [11,] "Olsen"   "Olsen"   "Mehlich" "Mehlich" "Bray"    "Bray"