library(readxl)
library(outliers)
#factor: presion
#VARIABLE RESPUESTA: TIEMPO DE FRACTURA (tf)
datos = read_excel("C://Users//DELL//Documents//diseño//parcialdiseño.xlsx", sheet = 1)

View(datos)
datos$presion = as.factor(datos$presion)
datos$tiempo = as.numeric(datos$tiempo)

Analisis descriptivo

arbol de decision

library(collapsibleTree)
collapsibleTreeSummary(datos,
                      c('presion',
                         'tiempo'),
                       collapsed=FALSE)
library(ggplot2)
ggplot(data=datos,aes(x=presion, y=tiempo, fill=presion))+
  geom_col()

class(datos$tiempo)
## [1] "numeric"
ggplot(data=datos,aes(x=presion, y=tiempo, fill=presion))+
  geom_boxplot()

### Hipotesis H0:μp1=μp2=μp3=μp4 Ha:μp1≠μp2≠μp3≠μp4

mod1 = aov(tiempo ~ presion,
          datos)
sum1<-summary(mod1)
sum1 <- unlist(sum1)
sum1 = sum1[9]
summary(mod1)
##             Df   Sum Sq  Mean Sq F value  Pr(>F)    
## presion      3 0.016567 0.005522   82.83 2.3e-06 ***
## Residuals    8 0.000533 0.000067                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Conclusion se rechaza la hipotesis nula, por lo tanto por lo menos un par de medias si es diferente y tiene efecto sobre el tiempo de ruptuta del bloque de madera

#Comparacion de medias posterior al analisis de varianza

#Prueba de maxima direfencia de Tukey
par(mar=c(5,6,3,1))
tt = TukeyHSD(mod1, 'presion') 
plot(tt, las=1)
abline(v=0, lty=2, col='red', lwd=2)

tt
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = tiempo ~ presion, data = datos)
## 
## $presion
##                 diff          lwr        upr     p adj
## 0.10-0.05 0.05000000  0.028650987 0.07134901 0.0003178
## 0.20-0.05 0.06666667  0.045317653 0.08801568 0.0000396
## 0.25-0.05 0.10333333  0.081984320 0.12468235 0.0000014
## 0.20-0.10 0.01666667 -0.004682347 0.03801568 0.1344163
## 0.25-0.10 0.05333333  0.031984320 0.07468235 0.0002012
## 0.25-0.20 0.03666667  0.015317653 0.05801568 0.0025538
library(TukeyC)
tt = TukeyC(mod1, 'presion')
plot(tt)

conclusion Hay diferencias entre todos los tratamientos excepto entre presion 0.1 y presion 0.2

#Presiones que difieren y no difieren

library(agricolae)
dt = duncan.test(mod1, 'presion', console = T)
## 
## Study: mod1 ~ "presion"
## 
## Duncan's new multiple range test
## for tiempo 
## 
## Mean Square Error:  6.666667e-05 
## 
## presion,  means
## 
##         tiempo         std r  Min  Max
## 0.05 0.8800000 0.010000000 3 0.87 0.89
## 0.10 0.9300000 0.010000000 3 0.92 0.94
## 0.20 0.9466667 0.005773503 3 0.94 0.95
## 0.25 0.9833333 0.005773503 3 0.98 0.99
## 
## Alpha: 0.05 ; DF Error: 8 
## 
## Critical Range
##          2          3          4 
## 0.01537336 0.01602049 0.01638221 
## 
## Means with the same letter are not significantly different.
## 
##         tiempo groups
## 0.25 0.9833333      a
## 0.20 0.9466667      b
## 0.10 0.9300000      c
## 0.05 0.8800000      d

Revision de suspuestos

var_res = tapply(mod1$residuals, datos$presion, var)

#Igualdad de varianzas
bartlett.test(mod1$residuals, datos$presion)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  mod1$residuals and datos$presion
## Bartlett's K-squared = 0.95233, df = 3, p-value = 0.8128

Segun evidencia estadistica se concluye que con alfa mayor a 0.05 el supuesto de homoceasticidad si se cumple

#Normalidad de residuos
shapiro.test(mod1$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  mod1$residuals
## W = 0.94102, p-value = 0.5114

Segun evidencia estadistica se concluye que con alfa mayor a 0.05 el supuesto de normalidad en los datos si se cumple

Modelo matematico del diseño

Yij=μ+φi+ϵij μ:Media global del tiempo φ0.5:Efecto de la presion a 0.5 atm φ0.10:Efecto de la presion a 0.10 atm φ0.20:Efecto de la presion a 0.20 atm φ0.25:Efecto de la presion a 0.25 atm ϵij:Error

#INTERPRETACION BIOLOGICA La presión influye considerablemente en el tiempo de ruptura, teniendo en cuenta que al realizar el analisis de varianza podemos concluir que las media del tiempo para cada presión son diferentes confirmandolo con la prueba de tukey donde se obtuvieron tres agrupaciones (a,b,c) las cuales indican que las presiones de 0.25(a) y 0.05(c) difieren en gran medida en comparacion con las presiones de 0.10 (b) y 0.20 (b).

PRUEBA DE GRUBBS Mide la ausencia de atipicos #HIPOTESIS H0:NO hay datos atipicos Ha:Hay datos atipicos

library(outliers)
grub1<-grubbs.test(mod1$residuals)
ifelse(grub1[3]>0.05, "Ho es verdad","Ha es verdad")
##        p.value 
## "Ho es verdad"
grub1
## 
##  Grubbs test for one outlier
## 
## data:  mod1$residuals
## G.5 = 1.43614, U = 0.79545, p-value = 0.8393
## alternative hypothesis: lowest value -0.0100000000000004 is an outlier

Se acepta la H0, debido a que tiene un p-valor >5%, es decir, no se detectaron valores atipicos en los datos.

library(outliers)
grubbs.test(mod1$residuals)
## 
##  Grubbs test for one outlier
## 
## data:  mod1$residuals
## G.5 = 1.43614, U = 0.79545, p-value = 0.8393
## alternative hypothesis: lowest value -0.0100000000000004 is an outlier

#CUATRO MODALIDADES DEL ANALISIS DE VARIANZA - Al revisar 1000 permutaciones se confirma que los datos y el modelo mantiene una diferencia de varianzas

#HIPOTESIS ##Comparar las medias de la presion H0:μp1=μp2=μp3=μp4 Ha:μp1≠μp2≠μp3≠μp4

library(RVAideMemoire)
## *** Package RVAideMemoire v 0.9-82-2 ***
## 
## Attaching package: 'RVAideMemoire'
## The following object is masked from 'package:TukeyC':
## 
##     cv
perm1<-perm.anova(tiempo ~ presion, data = datos, nperm = 1000, progress = F)
ifelse(perm1[5]>0.05, "Ho es verdad","Ha es verdad")
##           Pr(>F)        
## presion   "Ha es verdad"
## Residuals NA

#Análisis de varianza no parametricos

krus1<- kruskal.test(datos$tiempo, datos$presion);krus1
## 
##  Kruskal-Wallis rank sum test
## 
## data:  datos$tiempo and datos$presion
## Kruskal-Wallis chi-squared = 10.274, df = 3, p-value = 0.01637
ifelse(krus1[3]>0.05, "Ho es verdad","Ha es verdad")
##        p.value 
## "Ha es verdad"