require(faraway)
## Loading required package: faraway
data(rats)

rats
##    time poison treat
## 1  0.31      I     A
## 2  0.82      I     B
## 3  0.43      I     C
## 4  0.45      I     D
## 5  0.45      I     A
## 6  1.10      I     B
## 7  0.45      I     C
## 8  0.71      I     D
## 9  0.46      I     A
## 10 0.88      I     B
## 11 0.63      I     C
## 12 0.66      I     D
## 13 0.43      I     A
## 14 0.72      I     B
## 15 0.76      I     C
## 16 0.62      I     D
## 17 0.36     II     A
## 18 0.92     II     B
## 19 0.44     II     C
## 20 0.56     II     D
## 21 0.29     II     A
## 22 0.61     II     B
## 23 0.35     II     C
## 24 1.02     II     D
## 25 0.40     II     A
## 26 0.49     II     B
## 27 0.31     II     C
## 28 0.71     II     D
## 29 0.23     II     A
## 30 1.24     II     B
## 31 0.40     II     C
## 32 0.38     II     D
## 33 0.22    III     A
## 34 0.30    III     B
## 35 0.23    III     C
## 36 0.30    III     D
## 37 0.21    III     A
## 38 0.37    III     B
## 39 0.25    III     C
## 40 0.36    III     D
## 41 0.18    III     A
## 42 0.38    III     B
## 43 0.24    III     C
## 44 0.31    III     D
## 45 0.23    III     A
## 46 0.29    III     B
## 47 0.22    III     C
## 48 0.33    III     D

Contexto

El experimento evaluado consiste en una población de 48 ratas a las que se les suministraron una dosis de tres venenos disponibles a cada una y luego se les suministraró una dosis de 4 tratamientos disponibles a cada una.

require(ggplot2)
## Loading required package: ggplot2
ggplot(rats,aes(x=poison, y=time))+geom_point()+theme_bw()

ggplot(rats,aes(x=treat, y=time))+geom_point()+theme_bw()

Para la parte exploratorioa de este documento se realizaron los gráficos anteriores para visualizar como están distribuidos los datos y algunas de sus relaciones entre sí. En términos generales, el veneno III era mucho más agresivo con las ratas pues les daba menos de una hora de vida; mientras que el veneno II aunque podía ser bastante fuerte para algunas ratas, permitia que sobrevivieran por más tiempo con un récord máximo de aproximadamente 1 hora con 15 minutos.

En cuanto a los tratamientos, considerando a todas las ratas en general, el tratamiento A fue el menos efectivo pues ninguna rata sobrevivió por más de una hora, este es seguido por el C, el B y por último el tratamiento B siendo el más efectivo pues logró que las ratas pudieran sobrevivir por más de una hora aunque esta clasificación puede no ser muy apropiada dado que los datos de por ejemplo las ratas que recibieron el tratamiento B están bastante dispersos.

Ajuste del modelo

Indica los coeficientes de las diferentes variables del experimento-, en este caso son los venenos y los tratamientos.

(times=0.45229+(-0.07312*poisonII)+…)

mod=lm(time~poison+treat,data=rats)
summary(mod)
## 
## Call:
## lm(formula = time ~ poison + treat, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.25167 -0.09625 -0.01490  0.06177  0.49833 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.45229    0.05592   8.088 4.22e-10 ***
## poisonII    -0.07312    0.05592  -1.308  0.19813    
## poisonIII   -0.34125    0.05592  -6.102 2.83e-07 ***
## treatB       0.36250    0.06458   5.614 1.43e-06 ***
## treatC       0.07833    0.06458   1.213  0.23189    
## treatD       0.22000    0.06458   3.407  0.00146 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1582 on 42 degrees of freedom
## Multiple R-squared:  0.6503, Adjusted R-squared:  0.6087 
## F-statistic: 15.62 on 5 and 42 DF,  p-value: 1.123e-08
data("rats")
rats
##    time poison treat
## 1  0.31      I     A
## 2  0.82      I     B
## 3  0.43      I     C
## 4  0.45      I     D
## 5  0.45      I     A
## 6  1.10      I     B
## 7  0.45      I     C
## 8  0.71      I     D
## 9  0.46      I     A
## 10 0.88      I     B
## 11 0.63      I     C
## 12 0.66      I     D
## 13 0.43      I     A
## 14 0.72      I     B
## 15 0.76      I     C
## 16 0.62      I     D
## 17 0.36     II     A
## 18 0.92     II     B
## 19 0.44     II     C
## 20 0.56     II     D
## 21 0.29     II     A
## 22 0.61     II     B
## 23 0.35     II     C
## 24 1.02     II     D
## 25 0.40     II     A
## 26 0.49     II     B
## 27 0.31     II     C
## 28 0.71     II     D
## 29 0.23     II     A
## 30 1.24     II     B
## 31 0.40     II     C
## 32 0.38     II     D
## 33 0.22    III     A
## 34 0.30    III     B
## 35 0.23    III     C
## 36 0.30    III     D
## 37 0.21    III     A
## 38 0.37    III     B
## 39 0.25    III     C
## 40 0.36    III     D
## 41 0.18    III     A
## 42 0.38    III     B
## 43 0.24    III     C
## 44 0.31    III     D
## 45 0.23    III     A
## 46 0.29    III     B
## 47 0.22    III     C
## 48 0.33    III     D
table(rats$time)
## 
## 0.18 0.21 0.22 0.23 0.24 0.25 0.29  0.3 0.31 0.33 0.35 0.36 0.37 0.38  0.4 0.43 
##    1    1    2    3    1    1    2    2    3    1    1    2    1    2    2    2 
## 0.44 0.45 0.46 0.49 0.56 0.61 0.62 0.63 0.66 0.71 0.72 0.76 0.82 0.88 0.92 1.02 
##    1    3    1    1    1    1    1    1    1    2    1    1    1    1    1    1 
##  1.1 1.24 
##    1    1
table(rats$poison)
## 
##   I  II III 
##  16  16  16
table(rats$treat)
## 
##  A  B  C  D 
## 12 12 12 12
rats$time_=as.factor(rats$time)
rats$poison_=as.factor(rats$poison)
rats$treat_=as.factor(rats$treat)

rats$time_
##  [1] 0.31 0.82 0.43 0.45 0.45 1.1  0.45 0.71 0.46 0.88 0.63 0.66 0.43 0.72 0.76
## [16] 0.62 0.36 0.92 0.44 0.56 0.29 0.61 0.35 1.02 0.4  0.49 0.31 0.71 0.23 1.24
## [31] 0.4  0.38 0.22 0.3  0.23 0.3  0.21 0.37 0.25 0.36 0.18 0.38 0.24 0.31 0.23
## [46] 0.29 0.22 0.33
## 34 Levels: 0.18 0.21 0.22 0.23 0.24 0.25 0.29 0.3 0.31 0.33 0.35 0.36 ... 1.24
rats$poison_
##  [1] I   I   I   I   I   I   I   I   I   I   I   I   I   I   I   I   II  II  II 
## [20] II  II  II  II  II  II  II  II  II  II  II  II  II  III III III III III III
## [39] III III III III III III III III III III
## Levels: I II III
rats$treat_
##  [1] A B C D A B C D A B C D A B C D A B C D A B C D A B C D A B C D A B C D A B
## [39] C D A B C D A B C D
## Levels: A B C D

(##y=time) (#x1=poison) (#x2=treat)

mod_7=lm(time~poison+treat, data=rats)
summary(mod_7)
## 
## Call:
## lm(formula = time ~ poison + treat, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.25167 -0.09625 -0.01490  0.06177  0.49833 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.45229    0.05592   8.088 4.22e-10 ***
## poisonII    -0.07312    0.05592  -1.308  0.19813    
## poisonIII   -0.34125    0.05592  -6.102 2.83e-07 ***
## treatB       0.36250    0.06458   5.614 1.43e-06 ***
## treatC       0.07833    0.06458   1.213  0.23189    
## treatD       0.22000    0.06458   3.407  0.00146 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1582 on 42 degrees of freedom
## Multiple R-squared:  0.6503, Adjusted R-squared:  0.6087 
## F-statistic: 15.62 on 5 and 42 DF,  p-value: 1.123e-08
ggplot(rats,aes(x=treat, y=time, fill=poison))+geom_boxplot()+theme_bw()

ggplot(rats,aes(x=poison, y=time, fill=poison))+geom_boxplot()+theme_bw()

require(ggplot2)
ggplot(rats, aes(x=treat,y=time,fill=treat))+geom_boxplot()+theme_bw()

x=as.numeric(rats$treat=="yes")
y=rats$time

mod_0=lm(y~x)
summary(mod_0)
## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.29938 -0.17938 -0.07938  0.14312  0.76062 
## 
## Coefficients: (1 not defined because of singularities)
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   0.4794     0.0365   13.13   <2e-16 ***
## x                 NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2529 on 47 degrees of freedom

Las gráficas anteriores también cumplen la función de visualizar los datos del experimento.Podemos ver en el primer diagrama de cajas y bigotes que los diferentes venenos tienen diferentes efectos en cada rata según el tratamiento: Podemos ver que por ejemplo el veneno III es con certeza muy letal para las ratas sin importar el tratamiento pero el veneno II aunque es el menos letal para las ratas con tratamiento B, para el tratamiento C el veneno menos letal fue el I.

En la segunda gáfica de cajas y bigotes se puede ilustrat claramente que tan efectivos son los venenos en cuanto al tiempo que tardan en hacer efecto matándo al roedor.El veneno más efectivo es el III pues actúa muy rápidamente. Los otros dos venenos son muy parecidos en que las ratas pueden permanecer vivas por mucho más tiempo aunque los efectos en cada individuo sun más variados en comparación con el veneno III.

En el tercer diagrama de cajas y bigotes se considera el tratamiento recibido por cada rata. Podemos ver que el tratamiento con menor efecto fue el A pues la mortalidad de las ratas se presentaba bastante rápido mientras que el B permitia que las ratas pudieran sobrevibir por más tiempo aunque los efectos variaban en cada animal haciendo que el efecto se diera igual de rápido en algunas ratas que como en las ratas del tratamiento A.

#Modelo de diseño - ANOVA

require(faraway)

data(rats)

rats
##    time poison treat
## 1  0.31      I     A
## 2  0.82      I     B
## 3  0.43      I     C
## 4  0.45      I     D
## 5  0.45      I     A
## 6  1.10      I     B
## 7  0.45      I     C
## 8  0.71      I     D
## 9  0.46      I     A
## 10 0.88      I     B
## 11 0.63      I     C
## 12 0.66      I     D
## 13 0.43      I     A
## 14 0.72      I     B
## 15 0.76      I     C
## 16 0.62      I     D
## 17 0.36     II     A
## 18 0.92     II     B
## 19 0.44     II     C
## 20 0.56     II     D
## 21 0.29     II     A
## 22 0.61     II     B
## 23 0.35     II     C
## 24 1.02     II     D
## 25 0.40     II     A
## 26 0.49     II     B
## 27 0.31     II     C
## 28 0.71     II     D
## 29 0.23     II     A
## 30 1.24     II     B
## 31 0.40     II     C
## 32 0.38     II     D
## 33 0.22    III     A
## 34 0.30    III     B
## 35 0.23    III     C
## 36 0.30    III     D
## 37 0.21    III     A
## 38 0.37    III     B
## 39 0.25    III     C
## 40 0.36    III     D
## 41 0.18    III     A
## 42 0.38    III     B
## 43 0.24    III     C
## 44 0.31    III     D
## 45 0.23    III     A
## 46 0.29    III     B
## 47 0.22    III     C
## 48 0.33    III     D
mod=lm(time~poison,data=rats)
summary(mod)
## 
## Call:
## lm(formula = time ~ poison, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.31438 -0.15922 -0.03125  0.08594  0.69563 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.61750    0.05234  11.799 2.28e-15 ***
## poisonII    -0.07313    0.07401  -0.988    0.328    
## poisonIII   -0.34125    0.07401  -4.611 3.32e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2093 on 45 degrees of freedom
## Multiple R-squared:  0.3438, Adjusted R-squared:  0.3146 
## F-statistic: 11.79 on 2 and 45 DF,  p-value: 7.656e-05
anova(mod)
## Analysis of Variance Table
## 
## Response: time
##           Df Sum Sq Mean Sq F value    Pr(>F)    
## poison     2 1.0330 0.51651  11.786 7.656e-05 ***
## Residuals 45 1.9721 0.04382                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
mod=lm(time~treat,data=rats)
summary(mod)
## 
## Call:
## lm(formula = time ~ treat, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.38667 -0.15292 -0.01417  0.12833  0.56333 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.31417    0.06282   5.001 9.62e-06 ***
## treatB       0.36250    0.08885   4.080 0.000186 ***
## treatC       0.07833    0.08885   0.882 0.382739    
## treatD       0.22000    0.08885   2.476 0.017196 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2176 on 44 degrees of freedom
## Multiple R-squared:  0.3065, Adjusted R-squared:  0.2593 
## F-statistic: 6.484 on 3 and 44 DF,  p-value: 0.0009921
anova(mod)
## Analysis of Variance Table
## 
## Response: time
##           Df  Sum Sq  Mean Sq F value    Pr(>F)    
## treat      3 0.92121 0.307069  6.4836 0.0009921 ***
## Residuals 44 2.08388 0.047361                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Según lo anterior tenemos que:

El veneno III en comparación con el veneno I disminuye un 0.34125% en el tiempo que tarda en surtir efecto matando al animal y es significativo. Con el análisis de varianza realizado podemos ver que hay valores de veneno que presentan diferencias.

El tratamiento B en comparación con el tratamiento A incrementa un un 0.34125% el tiempo de superviviencia de las ratas y es significativo. Con el análisis de varianza realizado podemos ver que hay valores de tratamiento que presentan diferencias.

(time=0.4794 NA*treat)

##Prueba de Comparación Múltiple - Postanova

require (agricolae)
## Loading required package: agricolae
## Warning: package 'agricolae' was built under R version 4.1.1
mod=lm(time~poison,data=rats)
summary(mod)
## 
## Call:
## lm(formula = time ~ poison, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.31438 -0.15922 -0.03125  0.08594  0.69563 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.61750    0.05234  11.799 2.28e-15 ***
## poisonII    -0.07313    0.07401  -0.988    0.328    
## poisonIII   -0.34125    0.07401  -4.611 3.32e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2093 on 45 degrees of freedom
## Multiple R-squared:  0.3438, Adjusted R-squared:  0.3146 
## F-statistic: 11.79 on 2 and 45 DF,  p-value: 7.656e-05
compara=LSD.test(mod,"poison_")
compara
## NULL
mod=lm(time~treat,data=rats)
summary(mod)
## 
## Call:
## lm(formula = time ~ treat, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.38667 -0.15292 -0.01417  0.12833  0.56333 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.31417    0.06282   5.001 9.62e-06 ***
## treatB       0.36250    0.08885   4.080 0.000186 ***
## treatC       0.07833    0.08885   0.882 0.382739    
## treatD       0.22000    0.08885   2.476 0.017196 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2176 on 44 degrees of freedom
## Multiple R-squared:  0.3065, Adjusted R-squared:  0.2593 
## F-statistic: 6.484 on 3 and 44 DF,  p-value: 0.0009921
compara=LSD.test(mod,"treat_")
compara
## NULL
mod_1=lm(time~treat,data=rats)
summary(mod_1)
## 
## Call:
## lm(formula = time ~ treat, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.38667 -0.15292 -0.01417  0.12833  0.56333 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.31417    0.06282   5.001 9.62e-06 ***
## treatB       0.36250    0.08885   4.080 0.000186 ***
## treatC       0.07833    0.08885   0.882 0.382739    
## treatD       0.22000    0.08885   2.476 0.017196 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2176 on 44 degrees of freedom
## Multiple R-squared:  0.3065, Adjusted R-squared:  0.2593 
## F-statistic: 6.484 on 3 and 44 DF,  p-value: 0.0009921

Con los datos anteriores en emente tenemos que los venenos I y II son muy parecidos pues comparten intervalos indicados por la letra a. Una situación parecida ocurre con los tratamientos D y C al compartir ab y bc; esto por su parte indica la relación con los otros tratamientos A y B.

##veamos la variable poison

table(rats$poison)
## 
##   I  II III 
##  16  16  16
rats$poison_=as.factor(rats$poison)

rats$poison_
##  [1] I   I   I   I   I   I   I   I   I   I   I   I   I   I   I   I   II  II  II 
## [20] II  II  II  II  II  II  II  II  II  II  II  II  II  III III III III III III
## [39] III III III III III III III III III III
## Levels: I II III
mod_2=lm(time~poison_, data=rats)
summary(mod_2)
## 
## Call:
## lm(formula = time ~ poison_, data = rats)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.31438 -0.15922 -0.03125  0.08594  0.69563 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.61750    0.05234  11.799 2.28e-15 ***
## poison_II   -0.07313    0.07401  -0.988    0.328    
## poison_III  -0.34125    0.07401  -4.611 3.32e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2093 on 45 degrees of freedom
## Multiple R-squared:  0.3438, Adjusted R-squared:  0.3146 
## F-statistic: 11.79 on 2 and 45 DF,  p-value: 7.656e-05

el beta de poison_III indica que el efecto de 3% de poison es de 0.34125 menos respecto a Poison 1% y es significativo.

(time= 0.61750 + poison)