10/7/2019
Retomemos el ejemplo
Tema
Efecto de un herbicida sobre el desarrollo radicular de plántulas
Variable medida: Longitud de la raíz de la plántula
Tratamientos: Distintas dosis de herbicida
Mediciones: al día 3 y 5.
Muestreo aleatorio en cada vez.
Cálculo sugerido por la bibliografía
\(ln(\frac {l_{ij}} {lmean_{i0} - l_{ij}} + 150)\)
Donde \(l_{ij}\) es la longitud de la raíz de la plántula para el día \(i\) del tratamiento \(j\).
Sigamos con el ejemplo
# Cargo librerías
library(knitr)
library(kableExtra)
library(readxl)
# Cargo datos
datos <-read_xlsx("./datos/datos_limpios.xlsx")
datos <- as.data.frame(datos) # redefino el objeto como un dataframe
# Separo la planilla en función del día de medición
datos3 <- datos[datos$dia==3,] # datos para el día 3
datos5 <- datos[datos$dia==5,] # datos para el día 5
# Calculo promedio de controles
lmean0_3 <- mean(datos3[datos3$tratamiento=="0mg_ml",2])
lmean0_5 <- mean(datos5[datos5$tratamiento=="0mg_ml",2])
# Calculo el índice para cada día de medición datos3$indice<-log((datos3$longitud/(lmean0_3-datos3$longitud))+150) datos5$indice<-log((datos5$longitud/(lmean0_5-datos5$longitud))+150) # elimino los valores para el tratamiento control # asignandoles el valor NA (dato faltante) datos3[datos3$tratamiento=="0mg_ml","indice"] <- NA datos5[datos5$tratamiento=="0mg_ml","indice"] <- NA # concateno los resultados del índice # lo asigno a una variable nueva "indice" en el dataframe original datos$indice <- c(datos3$indice,datos5$indice)
| tratamiento | longitud | dia | indice |
|---|---|---|---|
| 0mg_ml | 15.70 | 3 | NA |
| 0mg_ml | 21.95 | 3 | NA |
| 0mg_ml | 14.60 | 3 | NA |
| 0mg_ml | 19.55 | 3 | NA |
| 0mg_ml | 13.87 | 3 | NA |
| 0mg_ml | 14.83 | 3 | NA |
Exploración
datos$indice datos$tratamiento plot(datos$tratamiento,datos$indice)
¡ERROR!
Factor
Variables
- Numérica
## [1] 15.70 21.95 14.60 19.55 13.87 14.83
- Categórica
## [1] 0mg_ml 0mg_ml 0mg_ml 0mg_ml 0mg_ml 0mg_ml ## Levels: 0,0005mg_ml 0,005mg_ml 0,1mg_ml 0mg_ml
- Ordinal
¿Qué es un factor?
- Una variable categórica.
- Tiene niveles
- El texto pasa a ser una etiqueta de cada nivel
Tratamiento como factor
datos$tratamiento<- as.factor(datos$tratamiento) datos$tratamiento
## [1] 0mg_ml 0mg_ml 0mg_ml 0mg_ml 0mg_ml ## [6] 0mg_ml 0mg_ml 0mg_ml 0,0005mg_ml 0,0005mg_ml ## [11] 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml ## [16] 0,0005mg_ml 0,005mg_ml 0,005mg_ml 0,005mg_ml 0,005mg_ml ## [21] 0,005mg_ml 0,005mg_ml 0,005mg_ml 0,005mg_ml 0,1mg_ml ## [26] 0,1mg_ml 0,1mg_ml 0,1mg_ml 0,1mg_ml 0,1mg_ml ## [31] 0,1mg_ml 0mg_ml 0mg_ml 0mg_ml 0mg_ml ## [36] 0mg_ml 0mg_ml 0mg_ml 0mg_ml 0,0005mg_ml ## [41] 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml 0,0005mg_ml ## [46] 0,0005mg_ml 0,0005mg_ml 0,005mg_ml 0,005mg_ml 0,005mg_ml ## [51] 0,005mg_ml 0,005mg_ml 0,005mg_ml 0,005mg_ml 0,1mg_ml ## [56] 0,1mg_ml 0,1mg_ml 0,1mg_ml 0,1mg_ml 0,1mg_ml ## Levels: 0,0005mg_ml 0,005mg_ml 0,1mg_ml 0mg_ml
plot(datos$tratamiento,datos$indice)
library(ggplot2) ggplot(data=datos,aes(x=tratamiento,y=indice))+ geom_point()+ stat_summary(fun.data = mean_se,geom="pointrange",color="red")+ facet_grid(dia~.)
¿qué análisis iría?
ACLARACIÓN
A continuación se muestran ejemplos de cómo se implementan algunos análisis básicos en R con el fin de familiarizarse con la dinámica de interacción con R. En ninguno de los procedimientos que vemos a continuación se evaluaron supuestos ni se hizo un análisis exploratorio para asegurarse que la herramienta utilizada sea la adecuada. Se muestran los procedimientos aplicados a los datos de ejemplo solo a título ilustrativo.
ANOVA
Día 3
modelo <-lm(data=datos3, indice ~ tratamiento) modelo
## ## Call: ## lm(formula = indice ~ tratamiento, data = datos3) ## ## Coefficients: ## (Intercept) tratamiento0,005mg_ml tratamiento0,1mg_ml ## 5.02817 0.00557 0.00252
Uno debe dialogar con R
modelo$coefficients
## (Intercept) tratamiento0,005mg_ml tratamiento0,1mg_ml ## 5.028171636 0.005569826 0.002519710
anova(modelo)
## Analysis of Variance Table ## ## Response: indice ## Df Sum Sq Mean Sq F value Pr(>F) ## tratamiento 2 0.0001244 0.00006222 0.1369 0.8729 ## Residuals 20 0.0090924 0.00045462
coef(modelo)
## (Intercept) tratamiento0,005mg_ml tratamiento0,1mg_ml ## 5.028171636 0.005569826 0.002519710
summary(modelo)
## ## Call: ## lm(formula = indice ~ tratamiento, data = datos3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -0.046379 -0.011107 -0.004740 0.007265 0.048042 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 5.028172 0.007538 667.006 <0.0000000000000002 *** ## tratamiento0,005mg_ml 0.005570 0.010661 0.522 0.607 ## tratamiento0,1mg_ml 0.002520 0.011035 0.228 0.822 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 0.02132 on 20 degrees of freedom ## (8 observations deleted due to missingness) ## Multiple R-squared: 0.0135, Adjusted R-squared: -0.08515 ## F-statistic: 0.1369 on 2 and 20 DF, p-value: 0.8729
par(mfrow=c(2,2)) plot(modelo)
Día 5
modelo <-lm(data=datos5, indice ~ tratamiento) anova(modelo)
## Analysis of Variance Table ## ## Response: indice ## Df Sum Sq Mean Sq F value Pr(>F) ## tratamiento 2 0.0036165 0.00180827 2.3598 0.123 ## Residuals 18 0.0137931 0.00076628
Tratamiento como variable numérica
## [1] 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0005 0.0005 ## [11] 0.0005 0.0005 0.0005 0.0005 0.0005 0.0005 0.0050 0.0050 0.0050 0.0050 ## [21] 0.0050 0.0050 0.0050 0.0050 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 ## [31] 0.1000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0005 ## [41] 0.0005 0.0005 0.0005 0.0005 0.0005 0.0005 0.0005 0.0050 0.0050 0.0050 ## [51] 0.0050 0.0050 0.0050 0.0050 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000
plot(datos$tratamiento,datos$indice)
library(ggplot2) ggplot(data=datos,aes(x=tratamiento,y=indice))+ geom_point()+ stat_smooth(method = "lm", color="red")+ facet_grid(dia~.)
¿Qué análisis iría?
Regresión
modelo <-lm(data=datos3, indice ~ tratamiento) modelo
## ## Call: ## lm(formula = indice ~ tratamiento, data = datos3) ## ## Coefficients: ## (Intercept) tratamiento ## 5.0308936 -0.0005491
anova(modelo)
## Analysis of Variance Table ## ## Response: indice ## Df Sum Sq Mean Sq F value Pr(>F) ## tratamiento 1 0.0000000 0.00000001 0 0.9956 ## Residuals 21 0.0092169 0.00043890
coef(modelo)
## (Intercept) tratamiento ## 5.0308935984 -0.0005491239
¿Correlación?
#cor(datos$longitud,datos$indice)
cor.test(datos$longitud,datos$indice,
na.rm=T)
## ## Pearson's product-moment correlation ## ## data: datos$longitud and datos$indice ## t = 4.2279, df = 42, p-value = 0.0001246 ## alternative hypothesis: true correlation is not equal to 0 ## 95 percent confidence interval: ## 0.2978161 0.7255715 ## sample estimates: ## cor ## 0.5463872
¡Panmixis de correlaciones!
Si quiero correr correlaciones entre mucha variables de un dataframe.
# datos ejemplo de R iris datos<-iris[-5] datos
## Sepal.Length Sepal.Width Petal.Length Petal.Width ## 1 5.1 3.5 1.4 0.2 ## 2 4.9 3.0 1.4 0.2 ## 3 4.7 3.2 1.3 0.2 ## 4 4.6 3.1 1.5 0.2 ## 5 5.0 3.6 1.4 0.2 ## 6 5.4 3.9 1.7 0.4 ## 7 4.6 3.4 1.4 0.3 ## 8 5.0 3.4 1.5 0.2 ## 9 4.4 2.9 1.4 0.2 ## 10 4.9 3.1 1.5 0.1 ## 11 5.4 3.7 1.5 0.2 ## 12 4.8 3.4 1.6 0.2 ## 13 4.8 3.0 1.4 0.1 ## 14 4.3 3.0 1.1 0.1 ## 15 5.8 4.0 1.2 0.2 ## 16 5.7 4.4 1.5 0.4 ## 17 5.4 3.9 1.3 0.4 ## 18 5.1 3.5 1.4 0.3 ## 19 5.7 3.8 1.7 0.3 ## 20 5.1 3.8 1.5 0.3 ## 21 5.4 3.4 1.7 0.2 ## 22 5.1 3.7 1.5 0.4 ## 23 4.6 3.6 1.0 0.2 ## 24 5.1 3.3 1.7 0.5 ## 25 4.8 3.4 1.9 0.2 ## 26 5.0 3.0 1.6 0.2 ## 27 5.0 3.4 1.6 0.4 ## 28 5.2 3.5 1.5 0.2 ## 29 5.2 3.4 1.4 0.2 ## 30 4.7 3.2 1.6 0.2 ## 31 4.8 3.1 1.6 0.2 ## 32 5.4 3.4 1.5 0.4 ## 33 5.2 4.1 1.5 0.1 ## 34 5.5 4.2 1.4 0.2 ## 35 4.9 3.1 1.5 0.2 ## 36 5.0 3.2 1.2 0.2 ## 37 5.5 3.5 1.3 0.2 ## 38 4.9 3.6 1.4 0.1 ## 39 4.4 3.0 1.3 0.2 ## 40 5.1 3.4 1.5 0.2 ## 41 5.0 3.5 1.3 0.3 ## 42 4.5 2.3 1.3 0.3 ## 43 4.4 3.2 1.3 0.2 ## 44 5.0 3.5 1.6 0.6 ## 45 5.1 3.8 1.9 0.4 ## 46 4.8 3.0 1.4 0.3 ## 47 5.1 3.8 1.6 0.2 ## 48 4.6 3.2 1.4 0.2 ## 49 5.3 3.7 1.5 0.2 ## 50 5.0 3.3 1.4 0.2 ## 51 7.0 3.2 4.7 1.4 ## 52 6.4 3.2 4.5 1.5 ## 53 6.9 3.1 4.9 1.5 ## 54 5.5 2.3 4.0 1.3 ## 55 6.5 2.8 4.6 1.5 ## 56 5.7 2.8 4.5 1.3 ## 57 6.3 3.3 4.7 1.6 ## 58 4.9 2.4 3.3 1.0 ## 59 6.6 2.9 4.6 1.3 ## 60 5.2 2.7 3.9 1.4 ## 61 5.0 2.0 3.5 1.0 ## 62 5.9 3.0 4.2 1.5 ## 63 6.0 2.2 4.0 1.0 ## 64 6.1 2.9 4.7 1.4 ## 65 5.6 2.9 3.6 1.3 ## 66 6.7 3.1 4.4 1.4 ## 67 5.6 3.0 4.5 1.5 ## 68 5.8 2.7 4.1 1.0 ## 69 6.2 2.2 4.5 1.5 ## 70 5.6 2.5 3.9 1.1 ## 71 5.9 3.2 4.8 1.8 ## 72 6.1 2.8 4.0 1.3 ## 73 6.3 2.5 4.9 1.5 ## 74 6.1 2.8 4.7 1.2 ## 75 6.4 2.9 4.3 1.3 ## 76 6.6 3.0 4.4 1.4 ## 77 6.8 2.8 4.8 1.4 ## 78 6.7 3.0 5.0 1.7 ## 79 6.0 2.9 4.5 1.5 ## 80 5.7 2.6 3.5 1.0 ## 81 5.5 2.4 3.8 1.1 ## 82 5.5 2.4 3.7 1.0 ## 83 5.8 2.7 3.9 1.2 ## 84 6.0 2.7 5.1 1.6 ## 85 5.4 3.0 4.5 1.5 ## 86 6.0 3.4 4.5 1.6 ## 87 6.7 3.1 4.7 1.5 ## 88 6.3 2.3 4.4 1.3 ## 89 5.6 3.0 4.1 1.3 ## 90 5.5 2.5 4.0 1.3 ## 91 5.5 2.6 4.4 1.2 ## 92 6.1 3.0 4.6 1.4 ## 93 5.8 2.6 4.0 1.2 ## 94 5.0 2.3 3.3 1.0 ## 95 5.6 2.7 4.2 1.3 ## 96 5.7 3.0 4.2 1.2 ## 97 5.7 2.9 4.2 1.3 ## 98 6.2 2.9 4.3 1.3 ## 99 5.1 2.5 3.0 1.1 ## 100 5.7 2.8 4.1 1.3 ## 101 6.3 3.3 6.0 2.5 ## 102 5.8 2.7 5.1 1.9 ## 103 7.1 3.0 5.9 2.1 ## 104 6.3 2.9 5.6 1.8 ## 105 6.5 3.0 5.8 2.2 ## 106 7.6 3.0 6.6 2.1 ## 107 4.9 2.5 4.5 1.7 ## 108 7.3 2.9 6.3 1.8 ## 109 6.7 2.5 5.8 1.8 ## 110 7.2 3.6 6.1 2.5 ## 111 6.5 3.2 5.1 2.0 ## 112 6.4 2.7 5.3 1.9 ## 113 6.8 3.0 5.5 2.1 ## 114 5.7 2.5 5.0 2.0 ## 115 5.8 2.8 5.1 2.4 ## 116 6.4 3.2 5.3 2.3 ## 117 6.5 3.0 5.5 1.8 ## 118 7.7 3.8 6.7 2.2 ## 119 7.7 2.6 6.9 2.3 ## 120 6.0 2.2 5.0 1.5 ## 121 6.9 3.2 5.7 2.3 ## 122 5.6 2.8 4.9 2.0 ## 123 7.7 2.8 6.7 2.0 ## 124 6.3 2.7 4.9 1.8 ## 125 6.7 3.3 5.7 2.1 ## 126 7.2 3.2 6.0 1.8 ## 127 6.2 2.8 4.8 1.8 ## 128 6.1 3.0 4.9 1.8 ## 129 6.4 2.8 5.6 2.1 ## 130 7.2 3.0 5.8 1.6 ## 131 7.4 2.8 6.1 1.9 ## 132 7.9 3.8 6.4 2.0 ## 133 6.4 2.8 5.6 2.2 ## 134 6.3 2.8 5.1 1.5 ## 135 6.1 2.6 5.6 1.4 ## 136 7.7 3.0 6.1 2.3 ## 137 6.3 3.4 5.6 2.4 ## 138 6.4 3.1 5.5 1.8 ## 139 6.0 3.0 4.8 1.8 ## 140 6.9 3.1 5.4 2.1 ## 141 6.7 3.1 5.6 2.4 ## 142 6.9 3.1 5.1 2.3 ## 143 5.8 2.7 5.1 1.9 ## 144 6.8 3.2 5.9 2.3 ## 145 6.7 3.3 5.7 2.5 ## 146 6.7 3.0 5.2 2.3 ## 147 6.3 2.5 5.0 1.9 ## 148 6.5 3.0 5.2 2.0 ## 149 6.2 3.4 5.4 2.3 ## 150 5.9 3.0 5.1 1.8
library(psych) corr.test(datos)
## Call:corr.test(x = datos) ## Correlation matrix ## Sepal.Length Sepal.Width Petal.Length Petal.Width ## Sepal.Length 1.00 -0.12 0.87 0.82 ## Sepal.Width -0.12 1.00 -0.43 -0.37 ## Petal.Length 0.87 -0.43 1.00 0.96 ## Petal.Width 0.82 -0.37 0.96 1.00 ## Sample Size ## [1] 150 ## Probability values (Entries above the diagonal are adjusted for multiple tests.) ## Sepal.Length Sepal.Width Petal.Length Petal.Width ## Sepal.Length 0.00 0.15 0 0 ## Sepal.Width 0.15 0.00 0 0 ## Petal.Length 0.00 0.00 0 0 ## Petal.Width 0.00 0.00 0 0 ## ## To see confidence intervals of the correlations, print with the short=FALSE option