Taller IFIBYNE (M3.C)¶
Lo que sigue es el .rmd de la guía.
--- title: "\[ATR\] Taller IFIBYNE (M3.C)" author: "Nicolás Méndez, Luis de Haro" date: "2019/8/13" output:
pdf_document:
toc: true
toc_depth: 3
number_sections: true
urlcolor: blue ---
Tip del módulo: si tienen alguna duda, pidan ayuda :)
{r setup, include=FALSE} # Opciones globales knitr::opts_chunk\$set(echo = TRUE) knitr::opts_chunk\$set(out.width='90%', dpi=300) knitr::opts_chunk\$set(fig.align='center', dev = "png") options(max.print=50) library("tidyverse") library("gridExtra") library("reshape2") library("gplots")
Aprender a caminar corriendo (2/2)¶
Como la sintaxis se repite en todos los gráficos, la idea es tomar la "anatomía" de ggplot como referencia y hacer cualquier gráfico.
Los organizamos en varias categorías: 1D, 2D, 3D y otras.
Good luck!
Preparar entorno¶
-
Cargar la librería de ggplot. 2. Usar el dataset "iris".
iris
es otro dataset de ejemplo: una tabla de mediciones del ancho y largo de sépalos y pétalos en diferentes especies.
```{r cars, message=FALSE} library(tidyverse) # Incluye ggplot2
Usaremos datos de ejemplo en iris ```¶
Flashback¶
-
¿Qué columnas tiene el data frame
iris
? -
¿Qué tipo de datos contiene cada columna?
¿Cuales son contínuas? ¿Alguna debería considerarse discreta o categórica?
-
Resumir las variables agrupando por
Species
usando la media:En base R es:
aggregate(. ~ Species, iris, FUN = mean) %>% melt()
O en tidycosas:
```{r} # A veces no es más elegante iris_tidy_summary <- iris %>%
gather(key = "attribute", value = "measurement",-Species) %>%
group_by(Species, attribute) %>%
summarise(media = mean(measurement),
varianza = var(measurement))
```{r, include = F} iris_means <- aggregate(. \~ Species, iris, FUN = mean) %>% melt(id.vars = "Species") iris_melt <- melt(iris, id.vars = "Species") ```
## Gráficos 2D
Acá incluímos gráficos con dos ejes para variables contínuas.
Finally!
### geom_point()
Ejemplos:
1. Hacer un gráfico de dispersión con `geom_point()` usando los datos de `iris`.
Busquen ayuda de nuestra parte, [acá](http://www.sthda.com/english/wiki/print.php?id=188) y [más allá](https://www.google.com).
2. Intenten colorear los puntos por la variable categórica `Species` poniendo un `aes()`.
3. Cambiar las variables graficadas. ¿Hay alguna otra combinación separe mejor las especies?
```{r, out.width="65%", eval = F} ggplot(iris, aes(Sepal.Length, Petal.Length)) + # Elegir dos variables
geom_point()
ggplot(iris, aes(Sepal.Length, Sepal.Width)) +
geom_point(aes(color = Species)) # Colorear por especie
Ejercicio 1¶
-
Hacer un gráfico de dispersión de
iris
y agregarle una capa degeom_smooth
con regresión lineal (una por especie).Hay otros métodos de regresión; consultar con
?geom_smooth
.
```{r, echo = F, out.width="65%"} ggplot(iris, aes(Sepal.Length, Sepal.Width, color = Species)) +
geom_jitter(width = .3, height = .3) +
geom_smooth(method = lm, aes(fill = Species)) +
theme_minimal()
- Como usamos `lm`, es relevante conocer la ecuación de la recta, el R2, etc.
<!-- -->
Para incluirlos en gráfico, hay que calcularlos con herramientas estadísticas y [agregarlos después](https://stackoverflow.com/questions/7549694/adding-regression-line-equation-and-r2-on-graph).
2. Para generar el siguiente gráfico puede ser útil:
- Desheredar los `aes()` en `geom_smooth` con `inherit_aes = FALSE`, y luego especificar `aes(Sepal.Length, Sepal.Width)`.
<!-- -->
- O simplemente poner `color = "black` en `geom_smooth()`. Esto sobreescribe el `aes()` original, desagrupando los datos.
```{r, echo = F, , out.width="65%"} ggplot(iris, aes(Sepal.Length, Sepal.Width, color = Species)) +
geom_jitter(width = .3, height = .3) +
geom_smooth(method = lm, mapping = aes(Sepal.Length, Sepal.Width), inherit.aes = F, color = "Black") +
theme_minimal()
Bonus: paradoja estadística¶
-
Dado el gráfico anterior y el siguiente, responda:
¿El ancho del sépalo aumenta o disminuye con su longitud?
¿Se puede resolver esa pregunta solo con estadística?
rta.
Bonus: todos contra todos¶
Sería interesante visualizar gráficos de todas las variables contra todas las variables.
Y sería más genial hacerlo en dos líneas.
Hay extensiones a ggplot, como [GGally](https://ggobi.github.io/ggally/#ggallyggpairs) que facilitan mucho la construcción de gráficos como esos.
Una lista de extensiones y ejemplos está disponible en: http://www.ggplot2-exts.org/gallery/
A continuación, un gráfico "todos contra todos" del dataset iris.
library("GGally")
ggpairs(iris, columns = 1:4,
mapping = aes(colour=Species),
diag = list(mapping = aes(alpha = .5)),
axisLabels = "show",
title = "ggpairs del iris dataset",
columnLabels = colnames(iris[, 1:4])) +
theme_minimal()
Desde acá toma la posta el compañero Luis :) con sus metabolitos. Gracias!
\newpage
-
En ciertos casos, nos puede interesar hacer un gráfico por grupo de observaciones (de acuerdo a alguna variable categórica).
- Cargamos los datos
read.table
indicando a la función que las columnas están separadas por tabs y que la tabla posee nombres de columna
```{r, message=F, warning=F, results=F} ### Leemos la tabla GCMS_data <- read.table("data/GCMS_Fruits_Landraces.txt", header = T, sep = "\t")
¿Se cargaron todas las filas y columnas? dim(GCMS_data) str(GCMS_data) ```¶
Tenemos datos de 78 metabolitos, para 16 genotipos distintos de tomate, en dos condiciones distintas de crecimiento (n = 4).
2. Graficamos algunos de los metabolitos
```{r, message=F, warning=F, results="hold"} ### Elegimos el metabolito "Adenine" y lo graficamos según el genotipo y en ambos tratamientos ggplot(GCMS_data, aes(x=Genotype, y=Adenine, fill=Condition)) +
geom_boxplot()
podemos cambiar el metabolito (y=XXX) para ver distintos boxplots ggplot(GCMS_data, aes(x=Genotype, y=Pyroglutamic.acid, fill=Condition)) +¶
geom_boxplot()
``` Dado que nos interesa conocer el efecto del tratamiento, si tuvieramos una escala independiente para cada genotipo podríamos visualizar mejor el cambio.
3. Una forma de lograr esto sin transformar los datos, es agregar las facets: `facet_wrap` o `facet_grid`
En el siguiente ejemplo `scales = "free"` implica que la escala de cada facet es independiente de las otras.
`theme(axis.text.y=element_text(size=6))` es para cambiar el tamaño del texto.
```{r, message=F, warning=F} ggplot(data = GCMS_data, aes(x=Condition, y=Pyroglutamic.acid )) +
geom_boxplot(aes(fill=Condition)) +
facet_wrap( ~ Genotype, scales="free") +
theme(axis.text.y=element_text(size=6))
### Heatmaps
¿Qué pasa si queremos graficar el cambio de todos los metabolitos juntos? ¿Hacemos 78 gráficos distintos?
Podemos simplificar la representación a través de un heatmap, coloreando según el fold change (veces de cambio) de las medias entre condiciones.
Lamentablemente, el soporte de ggplot para heatmaps es limitado. Así que vamos a usar otra librería, llamada "gplot".
0. Instalar gplot
```{r, eval = F} install.packages("gplot")
library(gplot) ```
1. Preparar los datos
```{r, message=F, warning=F} ### Calculamos la media para todos los metabolitos para cada genotipo y cada tratamiento x_means <- GCMS_data %>% group_by(Genotype, Condition) %>%
summarise_all(funs(mean))
### Removemos la columna que indica el número de réplica x_means <- subset(x_means, select=-c(Replicate))
### Separamos los tratamientos x_means_HT <- subset(x_means, Condition == "HT") x_means_LT <- subset(x_means, Condition == "LT")
### Calculamos el logaritmo del ratio x_log_ratio <- log(x_means_HT\[, 3:ncol(x_means_HT)\] / x_means_LT\[, 3:ncol(x_means_LT)\], 2) # FC
### Renombramos las filas (para la función heatmap) row.names(x_log_ratio) <- x_means_HT\$Genotype ```
2. Construimos el heatmap.
```{r} ### Elegimos la paleta de colores que m?s nos guste my_palette <- colorRampPalette(c("blue", "white", "red"))(n = 50) ```
```{r, include = F} # png(filename=paste("heatmap", ".png", sep = ""), # type="cairo", # units="cm", # width=50, # height=50, # pointsize=8, # res=300)
# install.packages("gplots") ```
```{r, fig.height=20, fig.width=20} # Necesitamos instalar y cargar el paquete "gplot" (si, con UNA g). heatmap.2(
t(x_log_ratio),
scale = "none",
trace = "none",
na.color = "dark grey",
#cellnote=t(p),
#notecol="black",
#notecex=1.5,
col = my_palette,
margins = c(8, 25),
cexRow = 1,
cexCol = 1.7,
key.xlab = "log2(Fold Change)",
tracecol="black"
#RowSideColors = sample(brewer.pal(8, name = "Accent"), ncol(x), replace = T),
#ColSideColors = sample(brewer.pal(8, name = "Accent"), nrow(x), replace = T)
) ```
```{r, include = F} # dev.off() ```
#### Bonus: resultados de la encuesta
Con gráficos tipo "alluvial" o "sankey"\[\^alluvial\]. \[\^alluvial\]: Hay bastante material para customizar esto \[acá\](<http://corybrunson.github.io/ggalluvial/>), \[acá\](<http://corybrunson.github.io/ggalluvial/articles/ggalluvial.html>), \[acá\](<https://cran.r-project.org/web/packages/ggalluvial/ggalluvial.pdf>), \[acá\](<https://www.datisticsblog.com/2018/10/intro_easyalluvial/>)
```{r, fig.height=10, fig.width=10, out.width="95%", message=F, warning=F} library("tidyverse") library("ggalluvial")
# Cargar y amasar datos respuestas <- read_csv("data/encuesta/respuestas.csv") %>%
mutate(timestamp = as.POSIXct(timestamp,
tz = "GMT-3",
format = " %Y/%m/%d %I:%M:%S %p")) # 2019/07/29 2:40:37 pm GMT-3
data_long_a <- respuestas %>%
group_by(exposicion, pref_modulo_5, pref_modulo_6) %>%
count() %>%
mutate(fill_stat = factor(
exposicion,
levels = c(
"Ya no tengo (tanto) miedo",
"Lo vi de lejos",
"Debo mi código a Stackoverflow"
)
)) %>%
ungroup %>%
arrange(fill_stat) %>%
mutate(subject = seq(1, n())) %>%
gather(key, value,-n ,-subject,-fill_stat) %>%
mutate(key = factor(key, levels = c(
"exposicion", "pref_modulo_5", "pref_modulo_6"
))) %>%
arrange(key, fill_stat)
#Graficar data_long_a %>%
ggplot(aes(
x = key,
y = n,
label = str_wrap(value, width = 10),
stratum = value,
alluvium = subject
)) +
geom_alluvium(aes(fill = fill_stat), aes.bind = TRUE) +
geom_stratum() +
geom_text(stat = "stratum", size = 3) +
scale_fill_manual(values = c("blue", "red", "yellow", "blue")) +
theme_void()