PARTE II: ANÁLISIS DE LA VARIANZA DE UN FACTOR

Problema 10 [Ver Script resuelto]

En Misiones coexisten tres especies de boyeros (Aves, Icteridae). Se midió la longitud del pico (en mm) para poner a prueba la hipótesis de que las especies de boyeros (Cacicus chrysopterus, C. haemorrhous, C. solitarius) difieren significativamente en el largo de sus picos. Los datos obtenidos se encuentran en el archivo BD_boyeros.txt:

Problema 11 [Ver script resuelto]

El cultivo de arveja (Pisum sativum L.), constituye actualmente un cultivo de alta importancia y gran demanda en el mercado nacional e internacional. Además de ser un cultivo rentable y de gran importancia alimenticia, tiene una gran capacidad de fijación simbiótica de nitrógeno atmosférico y como tal es una buena opción dentro de un plan de rotación de cultivos ya sea a campo abierto o bajo invernadero. Como los productores de esta especie están interesados en mejorar el rendimiento de su producción de arvejas, deciden promover una investigación sobre los efectos del agregado de diversos azúcares sobre el crecimiento. Para ello se hicieron crecer cortes de cotiledón de arveja en cultivos con auxina bajo la hipótesis previa de que se observará un crecimiento diferencial entre

  1. el control y los azúcares,
  2. azúcares puros y mezclas y
  3. azúcares monosacáridos y disacáridos.

Se obtuvieron los siguientes resultados expresados en unidades oculares (1 u.o. = 0,114 mm, datos en BD_arvejas.txt).

11.1.- Describa gráficamente los datos e identifique la unidad experimental, las variables respuesta y explicatoria, junto con sus tratamientos y cantidad de réplicas.

datos <- read.table("BD_arvejas.txt", header = T)
plot(datos)

Se realiza una prueba estadística ANOVA 1 factor (DCA) \(Y_{ij} = μ + α_i + ε_{ij}\).

11.2.- Enuncie las hipótesis previas planteadas en el enunciado, verifique la ortogonalidad y luego estudie la significación entre los tratamientos. Verificar los supuestos del modelo. ¿Qué conclusiones puede sacar?

AYUDA PARA PROBAR ORTOGONALIDAD:

Definan los contrastes como vectores “f_i”. Los elementos que compongan a dichos vectores serán los coeficientes de las hipótesis para cada contraste. Para definir a los vectores respeten el orden alfabético de los niveles del factor (1ero control, 2do fructosa, 3ro glucosa, 4to mezcla, 5to sacarosa).

En este ejemplo genérico el contraste f1 tiene como -> Ho: 2 μA= μB + μC => 2 μA- μB - μC = 0,
mientras que el f2 tiene como -> Ho: μB = μC => μB - μC = 0
f1 <- c(2,-1,-1)  
f2 <- c(0,1,-1)
sum(f1*f2)  # Este comando sirve para calcular el producto escalar, ya que multiplica elemento a 
elemento de cada vector y luego suma sus resultados.

Contrastes: C F G M S

  1. c vs f,g,m,s -> 4-1-1-1-1=0
  2. m vs f,g,s -> 0-1-1+3-1=0
  3. s vs f,g -> 0-1-1+0+2=0
f1 <- c(4,rep(-1,4))
f2 <- c(0,-1,-1,3,-1)
f3 <- c(0,-1,-1,0,2)
sum(f1*f2*f3)
## [1] 0

Supuestos del ANOVA:

modelo <- aov(datos$longitud~datos$tratamiento,datos)
res <- modelo$residuals
pre <- modelo$fitted.values
plot(modelo)

## normalidad analíticamente 
shapiro.test(res)
## 
##  Shapiro-Wilk normality test
## 
## data:  res
## W = 0.95835, p-value = 0.07579
0.07579<0.05 #no rechazo H_o
## [1] FALSE
## homcedasticidad analíticamente
library(car)
## Loading required package: carData
leveneTest(longitud~tratamiento,data=datos, center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  4  1.4541  0.232
##       45
0.232<0.05 #no rechazo H_o
## [1] FALSE

Prueba de ANOVA

anova(modelo)
## Analysis of Variance Table
## 
## Response: datos$longitud
##                   Df Sum Sq Mean Sq F value    Pr(>F)    
## datos$tratamiento  4 1045.0 261.250  71.035 < 2.2e-16 ***
## Residuals         45  165.5   3.678                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
2.2e-16 <0.05 #rechazo H_o
## [1] TRUE

Contrastes

## Paso 1 

# 1) control vs f+g+s+m  =>   Ho: 4control - 1f -1g - 1s - 1m = 0
# 2) mezclas vs puros    =>   Ho: 3m - 1f -1g - 1s = 0
# 3) di    vs. mono     =>    Ho: 2s - 1g - 1f = 0

# c / f / g / m / s → orden alfabetico

f1 <- c(4, -1, -1, -1, -1)
f2 <- c(0, -1, -1, 3, -1)
f3 <- c(0, -1, -1, 0, 2)

## Paso 2 (Comprobación de ortogonalidad) 
sum(f1*f2*f3) # da 0
## [1] 0
## Paso 3 (matriz con coeficientes) 
contrastes<-cbind(f1,f2, f3)
# La función "cbind" agrega los vectores como columnas de una matriz contrastes
contrastes
##      f1 f2 f3
## [1,]  4  0  0
## [2,] -1 -1 -1
## [3,] -1 -1 -1
## [4,] -1  3  0
## [5,] -1 -1  2
## Paso 4 Se agregan los contrastes a la base de datos
contrasts(datos$tratamiento)<-contrastes

## Paso 5 ANOVA
modelo <- aov(longitud ~ tratamiento, data = datos)
## Paso 6 contrastes ortogonales
summary(modelo, split = list(tratamiento = list("control vs. resto" = 1, "mezclas vs. puros" = 2, "disac. vs. mono" = 3)))
##                                  Df Sum Sq Mean Sq F value   Pr(>F)    
## tratamiento                       4 1045.0   261.3   71.03  < 2e-16 ***
##   tratamiento: control vs. resto  1  800.0   800.0  217.52  < 2e-16 ***
##   tratamiento: mezclas vs. puros  1   48.1    48.1   13.09 0.000749 ***
##   tratamiento: disac. vs. mono    1  190.8   190.8   51.88  5.1e-09 ***
## Residuals                        45  165.5     3.7                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  • Los azucares disminuyen significativamente (p<2x10-6) el crecimiento del coleoptile de arveja con respecto al control.
  • Entre los azúcares, la mezcla de G-F produce la mayor inhibición del crecimiento con respecto a los azúcares puros (p=0.000749) y el disacárido (sacarosa ) produce la menor inhibición (p:5.1x10-9) con respecto a monosacáridos.

11.3.- Calcule la magnitud de algún efecto que considere relevante.

medias <- tapply(datos$longitud, datos$tratamiento, mean)
medias
##  control fructosa  glucosa   mezcla sacarosa 
##     69.9     58.2     59.3     58.0     64.1
medias.azucares <- mean(medias[2:5])
medias.azucares
## [1] 59.9
  • Control-azúcares→ 69.9-59.9 = 10 → Los coleoptiles en presencia de azúcar crecen, en promedio, 10 uo menos que en cultivos sin azúcar → 14.3% menos en promedio.
  • Disacárido vs monosacárido → 64.1 vs 59.3+58.2/2 = 58.75 → 5.35 uo en promedio → 8.34% menos que en disacáridos

Problema 12

Se sabe que el dióxido de carbono tiene un efecto crítico en el crecimiento microbiológico. Cantidades pequeñas de \(CO_2\) estimulan el crecimiento de la mayoría de los organismos, mientras que altas concentraciones inhiben el crecimiento de la mayor parte de ellos. Este último efecto se utiliza comercialmente cuando se almacenan productos alimenticios perecederos. Se realizó un estudio para investigar el efecto de \(CO_2\) sobre la tasa de crecimiento del Pseudomonas fragii, un corruptor de alimentos. Se administró dióxido de carbono a cinco presiones atmosféricas diferentes. Se tomó como respuesta el cambio porcentual (crecimiento) en la masa celular, después de una hora. Se utilizaron diez cultivos en cada nivel y se obtuvieron los datos que se muestran en BD_tasa_crecimiento.txt.

datos <- read.table("BD_tasa_crecimiento.txt",header=T)

12.1.- Identificar la unidad experimental. Definir y clasificar la variable respuesta y la variable explicatoria.

12.2.- Representar gráficamente los datos.

plot(datos$crecimiento~as.factor(datos$presion_CO2))

12.3.- ¿Qué modelo estadístico podría usarse para analizar los datos? Plantear la hipótesis correspondiente y resolver.

Como se desea comparar si existen diferencias de meida entre los tratamientos se usa ANOVA de 1 factor, cuyo modelo se describe mediante la fórmula: \[y_{ij}=\mu+\alpha_i+\varepsilon_{ij}\] - \(y_{ij}\): tasa de crecimiento de la observación j del tratamiento i. - μ: tasa de crecimiento media poblacional de todos las observaciones. - \(ɑ_i\) : efecto del tratamiento (presión de CO_2) sobre la tasa de crecimiento (común a todos los individuos que recibieron ese tratamiento). - \(ε_{ij }\) : residuo, error aleatorio no explicado por el modelo.

Planteamos las hiótesis:

12.4.-¿Qué suposiciones se están haciendo? (explorar gráfica y analíticamente, cuando sea posible, la validez de los supuestos). Dar la conclusión estadística y biológica en base al valor p de la prueba.

  1. Muestras obtenidas de forma aleatoria.
  2. Observaciones o datos medidos independientes entre sí. (simplemente se asume)
  3. Normalidad de residuos [ \(ε_{ij } \sim N(0,2)\) ] : Los datos \(Y_i\) (o alternativamente los residuos) a registrar en las distintas unidades muestrales deben ser variables aleatorias con distribución normal (se verifica grafica y analiticamente).
  4. Homogeneidad de varianzas: todos los residuos deben tener una variabilidad constante de . (se verifica grafica y analiticamente)
datos$presion_CO2 <- as.factor(datos$presion_CO2)
modelo <- aov(datos$crecimiento~datos$presion_CO2)
plot(modelo)

r <- resid(modelo)
shapiro.test(r)
## 
##  Shapiro-Wilk normality test
## 
## data:  r
## W = 0.97205, p-value = 0.28
0.28<0.05 # no rechazo H_o 
## [1] FALSE
library(car)
leveneTest(datos$crecimiento~datos$presion_CO2, center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  4  0.4876 0.7447
##       45
0.7447<0.05 # no rechazo H_o
## [1] FALSE
summary(modelo)
##                   Df Sum Sq Mean Sq F value Pr(>F)    
## datos$presion_CO2  4  10637  2659.3   123.5 <2e-16 ***
## Residuals         45    969    21.5                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
2e-16<0.05 # rechazo H_o
## [1] TRUE

Conclusión: dado que el p-valor es menor al nivel de significación especificado (α:0.05) se rechaza Ho. Existe al menos un tratamiento que difiere significativamente del resto.

12.4.- Suponiendo que \(H_0\) es rechazada, ¿qué otras hipótesis podrían plantearse? Ponerlas a prueba e interpretar biológicamente los resultados.
Dado que se rechaza \(H_o\) se puede realizar la prueba de Tukey para analizar cuál de los tratamientos influye más.

Prueba de Tukey Contrastes de Tukey (10)

library(emmeans)
library(multcompView)
library(multcomp)
## Loading required package: mvtnorm
## Loading required package: survival
## Loading required package: TH.data
## Loading required package: MASS
## 
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
## 
##     geyser
# tukey 
Tukey <- emmeans(modelo, pairwise~presion_CO2, adjust ="tukey") 
Tukey
## $emmeans
##  presion_CO2 emmean   SE df lower.CL upper.CL
##  0             58.0 1.47 45     55.1     61.0
##  0.083         44.5 1.47 45     41.5     47.4
##  0.29          35.9 1.47 45     32.9     38.8
##  0.5           26.6 1.47 45     23.7     29.6
##  0.86          15.5 1.47 45     12.6     18.5
## 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast     estimate   SE df t.ratio p.value
##  0 - 0.083       13.58 2.08 45  6.544  <.0001 
##  0 - 0.29        22.18 2.08 45 10.689  <.0001 
##  0 - 0.5         31.41 2.08 45 15.137  <.0001 
##  0 - 0.86        42.50 2.08 45 20.481  <.0001 
##  0.083 - 0.29     8.60 2.08 45  4.144  0.0013 
##  0.083 - 0.5     17.83 2.08 45  8.593  <.0001 
##  0.083 - 0.86    28.92 2.08 45 13.937  <.0001 
##  0.29 - 0.5       9.23 2.08 45  4.448  0.0005 
##  0.29 - 0.86     20.32 2.08 45  9.793  <.0001 
##  0.5 - 0.86      11.09 2.08 45  5.344  <.0001 
## 
## P value adjustment: tukey method for comparing a family of 5 estimates
# grafico 
plot(Tukey, comparisons = T, horizontal = F)

# para tabla 
cld(Tukey$emmeans)
##  presion_CO2 emmean   SE df lower.CL upper.CL .group
##  0.86          15.5 1.47 45     12.6     18.5  1    
##  0.5           26.6 1.47 45     23.7     29.6   2   
##  0.29          35.9 1.47 45     32.9     38.8    3  
##  0.083         44.5 1.47 45     41.5     47.4     4 
##  0             58.0 1.47 45     55.1     61.0      5
## 
## Confidence level used: 0.95 
## P value adjustment: tukey method for comparing a family of 5 estimates 
## significance level used: alpha = 0.05

Se rechazan todas las hipótesis nulas entre contrastes dado que para todos el p valor es menor que el nivel de signficación. Es decir, existen evidencias singificativas para concluir que a diferntes preciones de CO_2, las tasas de crecimiento varían.

Problema 13 DUDA

Se estudió la capacidad de los nitropirenos,un grupo de contaminantes ambientales, para efectuar cambios en la capacidad reproductiva de la bacteria Salmonella sp,. que produce una enfermedad de transmisión alimentaria. Para ello se calculó la superficie cubierta por colonias de Salmonella (en cm) en placas tratadas con cuatro dosis del nitropireno 4NP. Los resultados obtenidos utilizando 28 placas fueron (Datos en BD_salmonella.txt):

datos <- read.table("BD_salmonella.txt",header = T)
datos$Dosis <- as.factor(datos$Dosis)

13.1.- Identificar y clasificar las variables del estudio.¿Qué análisis estadístico se podría utilizar para resolver la pregunta?

Análisis estadístico: ANOVA 1 factor

13.2.- Examinar los datos gráficamente y mediante estadística descriptiva.

library(psych)
## 
## Attaching package: 'psych'
## The following object is masked from 'package:car':
## 
##     logit
describeBy(datos$Superficie, datos$Dosis)
## 
##  Descriptive statistics by group 
## group: 0
##    vars n mean   sd median trimmed  mad min max range skew kurtosis   se
## X1    1 7 3.46 0.93    3.4    3.46 0.89 2.2   5   2.8 0.31    -1.33 0.35
## ------------------------------------------------------------ 
## group: 0.3
##    vars n  mean   sd median trimmed  mad min  max range skew kurtosis   se
## X1    1 7 10.23 1.98     10   10.23 2.08 7.8 13.4   5.6 0.36    -1.53 0.75
## ------------------------------------------------------------ 
## group: 1
##    vars n  mean   sd median trimmed  mad  min max range  skew kurtosis   se
## X1    1 7 21.89 3.09   22.6   21.89 2.67 17.6  26   8.4 -0.17    -1.74 1.17
## ------------------------------------------------------------ 
## group: 3
##    vars n  mean   sd median trimmed mad  min  max range  skew kurtosis   se
## X1    1 7 56.09 7.87   56.6   56.09 8.6 44.4 67.4    23 -0.04    -1.57 2.97
plot(datos)

13.3.- Verificar si se cumplen los supuestos de la prueba. Si se cumplen, resolver y dar la conclusión.

modelo <- lm(datos$Superficie~datos$Dosis)
plot(modelo) # no se cumple homocedasticidad

leveneTest(datos$Superficie~datos$Dosis,center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value    Pr(>F)    
## group  3   8.385 0.0005479 ***
##       24                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
0.0005479<0.05
## [1] TRUE

Probamos con una transoformación (¿¿??)

datos2 <- datos
datos2$Dosis <- as.factor(log(as.numeric(datos$Dosis), 10))
datos2$Superficie <- log(datos$Superficie, 10)
modelo2 <- lm(datos2$Superficie~datos2$Dosis)
plot(modelo2)

leveneTest(datos2$Superficie~datos2$Dosis, center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  3  1.0255 0.3989
##       24
0.3989<0.05 # no rechazo H_o
## [1] FALSE
shapiro.test(modelo2$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo2$residuals
## W = 0.99343, p-value = 0.9996
0.9996<0.05 # no rechazo H_o   
## [1] FALSE

Pasa los supuestos. Ahora planteamos la prueba.

anova(modelo)
## Analysis of Variance Table
## 
## Response: datos$Superficie
##             Df  Sum Sq Mean Sq F value    Pr(>F)    
## datos$Dosis  3 11486.4  3828.8  200.82 < 2.2e-16 ***
## Residuals   24   457.6    19.1                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
2.2e-16<0.05 # rechazo H_o. 
## [1] TRUE

Conclusión. Se concluye que al menos un log del promedio difiere del resto (existen diferencias significativas) pues el p.valor es menor al nivel de signficiación.

Problema 14

En un laboratorio de la facultad se quiso poner a prueba la efectividad de tres venenos contra cucarachas. Para ello se capturaron 24 cucarachas al azar y se asignaron aleatoriamente 8 a cada uno de tres tratamientos: veneno A, veneno B y veneno C. Se midió el tiempo de sobrevida de cada cucaracha como horas transcurridas desde que consumió el veneno hasta su muerte. Datos en BD_cucarachas.txt.

datos <- read.table("BD_cucarachas.txt",header = T)

14.1.- Definir la unidad experimental, la variable respuesta y la variable explicatoria con sus niveles.

14.2.- Indicar todos los supuestos del modelo y poner a prueba los que corresponda, planteando sus hipótesis.

modelo <- aov(datos$Sobrevida~datos$Veneno)
plot(modelo)

shapiro.test(resid(modelo))
## 
##  Shapiro-Wilk normality test
## 
## data:  resid(modelo)
## W = 0.94878, p-value = 0.255
0.255<0.05 # no rechazo H_o
## [1] FALSE
leveneTest(datos$Sobrevida~datos$Veneno, center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  2  0.1701 0.8447
##       21
0.8447<0.05 # no rechazo H_o
## [1] FALSE
#cumple los supuestos

14.3.- Plantear las hipótesis de la prueba estadística, resolver y concluir en términos estadísticos y biológicos.

summary(modelo)
##              Df Sum Sq Mean Sq F value   Pr(>F)    
## datos$Veneno  2  532.3   266.2   91.63 4.23e-11 ***
## Residuals    21   61.0     2.9                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
4.23e-11<0.05 #rechazo H_0
## [1] TRUE

Como se rechaza H_o, pues p.valor menor que el nivel de signficiación, hay evidencias para sostener que algún tratamiento hace difererir la media pobacional, por lo que no todas las medias pobacionales son iguales.

14.4- ¿Recomendarían un veneno? Realizar la prueba correspondiente y construir un gráfico que ilustre la respuesta.

Dado que ya rechazamos la H_o (no existe efecto debido al tratamiento con venenos sobre el tiempo de sobrevida (en horas)). Reealicemos ahora una prueba dee Tukey para comparar los tratamientos.

Tukey <- emmeans(modelo, pairwise~Veneno, adjust ="tukey") 
# grafico 
plot(Tukey, comparisons = T, horizontal = F)

# para tabla 
cld(Tukey$emmeans)
##  Veneno emmean    SE df lower.CL upper.CL .group
##  B        15.5 0.603 21     14.2     16.8  1    
##  A        17.2 0.603 21     16.0     18.5  1    
##  C        26.2 0.603 21     25.0     27.5   2   
## 
## Confidence level used: 0.95 
## P value adjustment: tukey method for comparing a family of 3 estimates 
## significance level used: alpha = 0.05

Recomendaría el veneno A o el veneno B, dado que no difiereen significativamente entre ellos.

Problema 15

Se sabe que los metales pesados pueden dañar el hepatopáncreas de la centolla (Lithodes santolla) , produciendo que excrete menos amonio. Se llevó a cabo un experimento a fin de evaluar la respuesta de la centolla a dosis sub-letales de distintos metales pesados. Para ello, ejemplares adultos fueron ubicados en recipientes individuales de vidrio con agua de mar artificial y se les aplicó uno de los siguientes tratamientos: a) Control, b) Cd; c) Cu; d) Zn, a 0,5 mg/l cada uno. A las 96 hs de exposición se midió el amonio excretado (en mg/litro). Los resultados se muestran en el archivo BD_Amonio.txt.

datos <- read.table("BD_Amonio.txt",header = T)

15.1. Identificar la unidad experimental, la variable respuesta, la variable explicatoria o factor y sus niveles, la cantidad de tratamientos y la cantidad de réplicas. ¿Cuál es el modelo? ¿Se trata de un estudio experimental u observacional?

15.2. Describir gráfica y estadísticamente las tres muestras.

plot(datos)

15.3.- Resolver comprobando previamente los supuestos.

  1. Muestras obtenidas de forma aleatoria.
  2. Observaciones o datos medidos independientes entre sí.
  3. Normalidad de residuos [ \(ε_{ij } \sim N(0,2)\) ] : Los datos \(Y_i\) (o alternativamente los residuos) a registrar en las distintas unidades muestrales deben ser variables aleatorias con distribución normal (se verifica grafica y analiticamente).
  4. Homogeneidad de varianzas: todos los residuos deben tener una variabilidad constante de . (se verifica grafica y analiticamente)
modelo <- aov(datos$amonio~datos$tratamiento, data=datos)
plot(modelo)

shapiro.test(modelo$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$residuals
## W = 0.97293, p-value = 0.328
0.328<0.05 # no rechazo H0
## [1] FALSE
leveneTest(datos$amonio~datos$tratamiento, center=mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  3  1.0279 0.3895
##       44
0.3895<0.05 # no rechazo 
## [1] FALSE
# Cumple los supuestos

Pureba de ANOVA

summary(modelo)
##                   Df Sum Sq Mean Sq F value   Pr(>F)    
## datos$tratamiento  3  262.8   87.60   32.07 3.81e-11 ***
## Residuals         44  120.2    2.73                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
3.81e-11<0.05
## [1] TRUE

15.4.- ¿Qué conclusión se puede sacar? De ser necesario, efectuar comparaciones. Existen evidencias signficativas para concluir que existen efectos debido a los tratamientos.

Tukey

Tukey <- emmeans(modelo, pairwise~tratamiento, adjust="tukey")
summary(Tukey)
## $emmeans
##  tratamiento emmean    SE df lower.CL upper.CL
##  Cd            8.43 0.477 44     7.47     9.39
##  control      14.47 0.477 44    13.51    15.43
##  Cu           12.47 0.477 44    11.51    13.44
##  Zn            9.79 0.477 44     8.83    10.75
## 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast     estimate    SE df t.ratio p.value
##  Cd - control    -6.03 0.675 44 -8.941  <.0001 
##  Cd - Cu         -4.04 0.675 44 -5.990  <.0001 
##  Cd - Zn         -1.36 0.675 44 -2.013  0.1988 
##  control - Cu     1.99 0.675 44  2.952  0.0251 
##  control - Zn     4.67 0.675 44  6.928  <.0001 
##  Cu - Zn          2.68 0.675 44  3.977  0.0014 
## 
## P value adjustment: tukey method for comparing a family of 4 estimates
plot(Tukey, comparison=T, horizontal=F)

cld(Tukey$emmeans)
##  tratamiento emmean    SE df lower.CL upper.CL .group
##  Cd            8.43 0.477 44     7.47     9.39  1    
##  Zn            9.79 0.477 44     8.83    10.75  1    
##  Cu           12.47 0.477 44    11.51    13.44   2   
##  control      14.47 0.477 44    13.51    15.43    3  
## 
## Confidence level used: 0.95 
## P value adjustment: tukey method for comparing a family of 4 estimates 
## significance level used: alpha = 0.05

15.5.- ¿Son igualmente tóxicos todos los metales pesados? Ordenar los metales según su toxicidad.

No. El Cd y el Zn no difieren signficativamente entre ellos en cuanto a toxicidad, mientras que estos sí difieren signficcativamente con el Cu.

Más tóxico a menos tóxico:(Zn y Cd), Cu.

Problema 16 [resuelto en campus]

En modelos animales se ha observado que la inhibición de la enzima convertidora de angiotensina (ECA) reduce la activación y acumulación de células inflamatorias y el estrés oxidativo, previniendo la formación de ateromas. Sin embargo, no existen ensayos clínicos en humanos.
Recientemente, unos investigadores efectuaron un ensayo con 90 pacientes mayores de 55 años de edad y alto riesgo de eventos cardiovasculares, con el propósito de averiguar si una droga inhibidora de ECA, Ramipril, era efectiva y qué dosis sería conveniente administrar. Los 90 pacientes fueron asignados aleatoriamente y en forma balanceada a alguna de tres dosis de Ramipril: 2.5 mg/día, 10 mg/día o 0 mg/día (placebo).

A los dos años de iniciado el tratamiento se midió, por ultrasonido, el espesor de la carótida (en micras/año) como indicador de ateroesclerosis (a mayor valor, mayor riesgo de ateroesclerosis). Los resultados se muestran en el archivo BD_ECA.txt).

datos <- read.table("BD_ECA.txt", header = T)
plot(datos)

16.1.- Identificar la unidad experimental, la variable respuesta, el factor, sus niveles y la cantidad de réplicas.

16.2.- ¿Se cumplen los supuestos del modelo? Fundamentar en forma gráfica o estadísticamente.

modelo <- lm(datos$espesor~datos$ramipril, data=datos)
plot(modelo)

shapiro.test(modelo$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo$residuals
## W = 0.97833, p-value = 0.1379
0.1379<0.05 #no rechazo
## [1] FALSE
library(car)
leveneTest(datos$espesor~datos$ramipril, data=datos, center = mean)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  2   1.234 0.2962
##       87
0.2962<0.05 #no rechazo 
## [1] FALSE

16.3.- ¿Es efectivo el uso del Ramipril para reducir el crecimiento del espesor de la carótida? ¿Qué dosis recomendaría usar?

anova(modelo)
## Analysis of Variance Table
## 
## Response: datos$espesor
##                Df Sum Sq Mean Sq F value    Pr(>F)    
## datos$ramipril  2 2105.3 1052.65   19.02 1.404e-07 ***
## Residuals      87 4815.0   55.34                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
1.404e-07<0.05 # rechazo 
## [1] TRUE
library(emmeans)
library(multcompView)
library(multcomp)
Tukey <- emmeans(modelo, pairwise~ramipril, adjust ="tukey")
Tukey
## $emmeans
##  ramipril emmean   SE df lower.CL upper.CL
##  10         18.6 1.36 87     15.9     21.3
##  2.5        12.8 1.36 87     10.1     15.5
##  placebo    24.6 1.36 87     21.9     27.3
## 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast      estimate   SE df t.ratio p.value
##  10 - 2.5          5.79 1.92 87  3.012  0.0094 
##  10 - placebo     -6.06 1.92 87 -3.155  0.0062 
##  2.5 - placebo   -11.85 1.92 87 -6.167  <.0001 
## 
## P value adjustment: tukey method for comparing a family of 3 estimates
cld(Tukey$emmeans)
##  ramipril emmean   SE df lower.CL upper.CL .group
##  2.5        12.8 1.36 87     10.1     15.5  1    
##  10         18.6 1.36 87     15.9     21.3   2   
##  placebo    24.6 1.36 87     21.9     27.3    3  
## 
## Confidence level used: 0.95 
## P value adjustment: tukey method for comparing a family of 3 estimates 
## significance level used: alpha = 0.05
plot(Tukey, comparisons = T, horizontal = F)

plot(Tukey)