Skip to content

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

  1. 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

  1. Hacer un gráfico de dispersión de iris y agregarle una capa de geom_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

  1. En ciertos casos, nos puede interesar hacer un gráfico por grupo de observaciones (de acuerdo a alguna variable categórica).

    1. 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()