Скачать презентацию Cómo analizar los datos crudos de microarrays Guía Скачать презентацию Cómo analizar los datos crudos de microarrays Guía

0569b0ceb879437876eae1507c42686b.ppt

  • Количество слайдов: 99

¿Cómo analizar los datos crudos de microarrays? Guía práctica de análisis de datos de ¿Cómo analizar los datos crudos de microarrays? Guía práctica de análisis de datos de microarrays Juan Pablo Fededa y Carlos Rocco Análisis Exploratorio y Confirmatorio de Datos de Experimentos de Microarrays Primer Cuatrimestre 2006 Instituto de Cálculo y Departamento de Matemática. Facultad de Ciencias Exactas y Naturales Universidad de Buenos Aires 1

IMPORTANTE: Esta es una guía acerca de cómo analizar los datos que se obtienen IMPORTANTE: Esta es una guía acerca de cómo analizar los datos que se obtienen del scanneado de un microarray de dos colores. No contiene información acerca de cómo construir, ensayar, diseñar un experimento de y/ó scannear microarrays. El objetivo final es ayudar a un biólogo a que, a partir del dato original de scanneado de un microarray, obtenga una lista de genes diferencialmente expresados en las dos condiciones confrontadas. Tampoco se hablará de data mining y generación de patrones de comportamiento génico. El análisis de datos que promueve esta guía involucra el uso de R, un software libre de estadística sobre el cuál se fueron construyendo múltiples herramientas para analizar microarrays. Se asume que el lector tiene un background matemático básico y que entiende los conceptos fundamentales de la tecnología de microarrays. 2

Microarrays…………… Fabricación de un microarray FLUJO DE TRABAJO u$s Diseño del experimento Preparación de Microarrays…………… Fabricación de un microarray FLUJO DE TRABAJO u$s Diseño del experimento Preparación de la muestra - hibridación Scanning Análisis de datos - lista genes diferencialmente expresados Data mining – clustering 3

Microarrays…………… FLUJO DE TRABAJO Fabricación de un microarray Diseño del experimento Preparación de la Microarrays…………… FLUJO DE TRABAJO Fabricación de un microarray Diseño del experimento Preparación de la muestra - hibridación Scanning Análisis de datos - lista genes diferencialmente expresados Data mining – clustering (……. . 2007!!!!) 4

FLUJO DE TRABAJO Análisis de Datos de Microarrays Obtención de los datos crudos – FLUJO DE TRABAJO Análisis de Datos de Microarrays Obtención de los datos crudos – Ingreso de los datos al R – carga de paquetes de análisis Gráficos de estado del microarray – toma de decisión ( sigo con el análisis o dejo? ) Normalizaciónes Gráficos de normalización – toma de decisión ( sigo con el análisis, realizo otras normalizaciones, dejo? ) Identificación de genes diferencialmente expresados Gráficos y tablas de genes diferencialmente expresados 5

FLUJO DE TRABAJO Análisis de Datos de Microarrays Obtención de los datos crudos – FLUJO DE TRABAJO Análisis de Datos de Microarrays Obtención de los datos crudos – Ingreso de los datos al R – carga de paquetes de análisis Gráficos de estado del microarray – toma de decisión ( sigo con el análisis o dejo? ) Normalizaciónes R Gráficos de normalización – toma de decisión ( sigo con el análisis, realizo otras normalizaciones, dejo? ) Identificación de genes diferencialmente expresados Gráficos y tablas de genes diferencialmente expresados 6

Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de los datos crudos Ingreso de los datos al R 7

-1: Obtener el software R !!!! de http: //cran. r-project. org/ 8 -1: Obtener el software R !!!! de http: //cran. r-project. org/ 8

-1: Obtener el software R !!!! de http: //cran. r-project. org/ 9 -1: Obtener el software R !!!! de http: //cran. r-project. org/ 9

-1: Obtener el software R !!!! de http: //cran. r-project. org/ 10 -1: Obtener el software R !!!! de http: //cran. r-project. org/ 10

0: Obtener ( bajar) los siguientes paquetes de trabajo de R en http: //www. 0: Obtener ( bajar) los siguientes paquetes de trabajo de R en http: //www. bioconductor. org/packages/bioc/1. 8/index. html : Affyaffyio, array. Quality, Biobase, colorspace, convert, grid. Base, hexbin, limma, marray, RColor. Brewer, vsn 11

Una vez hecho esto instalar el R y dentro de la consola del R Una vez hecho esto instalar el R y dentro de la consola del R instalar los paquetes bajados en archivos zip : NOTA: Se pueden actualizar e instalar desde R (en otra de las opciones de 12 package).

Una vez hecho esto instalar el R y dentro de la consola del R Una vez hecho esto instalar el R y dentro de la consola del R instalar los paquetes bajados en archivos zip : 13

Una vez hecho esto instalar el R y dentro de la consola del R Una vez hecho esto instalar el R y dentro de la consola del R instalar los paquetes bajados en archivos zip : 14

Antes que nada: vamos a trabajar escribiendo nuestras ordenes en un script, para lo Antes que nada: vamos a trabajar escribiendo nuestras ordenes en un script, para lo cuál primero hay que abrir un script: 15

Es más fácil trabajar con las ventanas de la consola y del script (Editor) Es más fácil trabajar con las ventanas de la consola y del script (Editor) en paralelo ( tile): 16

Entonces…. . esto funciona escribiendo la orden en el script y luego se utiliza Entonces…. . esto funciona escribiendo la orden en el script y luego se utiliza el botón run line para correrla en la consola: SCRIPT CONSOLA NOTA: run line se puede hacer con F 5. 17

La orden se corre en la consola: 18 La orden se corre en la consola: 18

Cargamos el paquete marray para empezar a trabajar (al cargar marray se carga automaticamente Cargamos el paquete marray para empezar a trabajar (al cargar marray se carga automaticamente limma ) 19

Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de los datos crudos Ingreso de los datos al R 20

OK; los microarrays sobre los cuáles vamos a trabajar como ejemplo corresponden al siguiente OK; los microarrays sobre los cuáles vamos a trabajar como ejemplo corresponden al siguiente paper: 21

A los datos originales los bajamos de la base de datos de microarrays de A los datos originales los bajamos de la base de datos de microarrays de stanford: http: //genome-www 5. stanford. edu/ 22

Vamos a la carpeta SMD, en donde están los datos crudos: 23 Vamos a la carpeta SMD, en donde están los datos crudos: 23

Particularmente vamos a analizar los datos de los microarrays correspondientes a los experimentos de Particularmente vamos a analizar los datos de los microarrays correspondientes a los experimentos de inmunoprecipitación de PUMILIO en moscas adultas: 24

Vamos a utilizar los datos originales del scanner (ORI DATA), en este caso son Vamos a utilizar los datos originales del scanner (ORI DATA), en este caso son 4 experimentos de hibridación: expt. ID 52987, 53530, 54253 y 54596 25

Utilizamos (bajamos) los archivos. gpr (genepix): 26 Utilizamos (bajamos) los archivos. gpr (genepix): 26

NOTA: Diferentes scanners generan diferentes tipos de archivos con diferentes extensiones; los paquetes de NOTA: Diferentes scanners generan diferentes tipos de archivos con diferentes extensiones; los paquetes de análisis de datos que utilizaremos leen la mayoría de estas extensiones. Agilent Arrayvision Bluefuse Genepix spot. close. open genepix. custom genepix. median Imagene Smd Quantarray Scanarrayexpress smd. old Spot El procedimiento para ingresar estas files en R es el mismo que el que veremos a continuación. 27

Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de Obtención de R y paquetes de análisis Carga de paquetes de análisis Obtención de los datos crudos Ingreso de los datos al R 28

Volviendo a R, cambiamos de directorio en el cuál están nuestros archivos originales salidos Volviendo a R, cambiamos de directorio en el cuál están nuestros archivos originales salidos del scanner (ej. : . gpr para genepix): 29

Cambiamos de directorio en el cuál están nuestros archivos originales salidos del scanner (ej. Cambiamos de directorio en el cuál están nuestros archivos originales salidos del scanner (ej. : . gpr para genepix): 30

Cambiamos de directorio en el cuál están nuestros archivos originales salidos del scanner (ej. Cambiamos de directorio en el cuál están nuestros archivos originales salidos del scanner (ej. : . gpr para genepix): 31

Empezamos a trabajar en R: dirección. experimentospumilio <- Empezamos a trabajar en R: dirección. experimentospumilio <-"D: //juan//análisis exploratorio y confirmatorio de experimentos de microarrays//curso microarrays//pumilio//Experiment Sets Name Organism Type Created By Options“ Esta es la otra opción para indicar la carpeta en donde están nuestros datos (Por defecto lee en la carpeta de trabajo). setwd(dirección. experimentospumilio) R va a leer los datos desde esta dirección seteada. archivosexperimentospumilio <- c("52987. gpr", "53530. gpr", "54253. gpr", "54596. gpr") Generamos un objeto marray. Raw con los archivos de los arrays que vamos a analizar. crudosexperimentospumilio<-read. Gene. Pix(archivosexperimentospumilio) Generamos un objeto que NO ignore las estimaciones de las intensidades del background para los canales verde ( G ) y rojo ( R ). crudosexperimentospumilio 2<-read. Gene. Pix(archivosexperimentospumilio, name. Gb= NULL, name. Rb= NULL) Generamos un objeto que ignore las estimaciones de las intensidades del background para los canales verde ( G ) y rojo ( R ) fijando el valor NULL para los argumentos name. Gb y name. Rb. 32

33 33

Estructura de los datos que estamos analizando en el paquete marray: slo ts etceteras Estructura de los datos que estamos analizando en el paquete marray: slo ts etceteras Gb (background Cy 3) Rb (background Cy 5) Gf (foreground Cy 3) microarreglos genes Rf (foreground Cy 5) 34

summary(crudosexperimentospumilio) Nos muestra un resumen de las características del objeto crudosexperimentospumilio que generamos: Son summary(crudosexperimentospumilio) Nos muestra un resumen de las características del objeto crudosexperimentospumilio que generamos: Son datos sin normalizar, es clase marray, contiene los datos de los 4 microarrays 17664 spots organizados en 8 filas x 4 columnas de subconjuntos ( grids ), cada grid contiene 23 filas x 24 columnas de spots 35

slot. Names (crudosexperimentospumilio) Nos muestra la estructura de slots del objeto crudosexperimentopumilio. slot. Names slot. Names (crudosexperimentospumilio) Nos muestra la estructura de slots del objeto crudosexperimentopumilio. slot. Names (crudosexperimentospumilio 2) Nos muestra la estructura de slots del objeto crudosexperimentopumilio 2. Ambos objetos (crudosexperimentospumilio, crudosexperimentospumilio 2 ) tienen la misma estructura de slots. Los primeros 5 ("ma. Rf" "ma. Gf" "ma. Rb" "ma. Gb" "ma. W“) son las matrices que contienen la información cuantitativa básica extraída de los archivos. gpr. Los demás están asociados con los archivos de la estructura de los arreglos (layout) y las 36 anotaciones.

crudosexperimentospumilio@ma. Rf Objeto con las intensidades del foreground para el canal rojo. crudosexperimentospumilio@ma. Rb [email protected] Rf Objeto con las intensidades del foreground para el canal rojo. [email protected] Rb Objeto con las intensidades del background para el canal rojo. [email protected] Gf Objeto con las intensidades del foreground para el canal verde. [email protected] Gb Objeto con las intensidades del background para el canal verde. [email protected] Layout Objeto con la estructura o geometría del arreglo. [email protected] Gnames Objeto con los nombres de los genes. [email protected] Targets Objeto con información acerca de que muestras fueron hibridadas en cada 37 arreglo.

crudosexperimentospumilio@ma. Gf[1: 5, ] Generamos un objeto con las intensidades del foreground para el [email protected] Gf[1: 5, ] Generamos un objeto con las intensidades del foreground para el canal verde para las 5 primeras filas (genes) para todos los microarrays En este caso el microarray de la columna 1 es 52987 [1: 5, 1 ] microarrays genes Recuerden la estructura de datos de marray!!! 38

Existen métodos específicos para objetos de la clase marray. Raw (como crudosexperimentospumilio), que luego Existen métodos específicos para objetos de la clase marray. Raw (como crudosexperimentospumilio), que luego se utilizaran para tomar decisiones acerca de los datos analizados: ma. A : matriz de log 2 de intensidades (con corrección por background) ma. M : matriz de log 2 de cocientes (con corrección por background) ma. LR : matriz de log 2 de ( intensidades - background) para el canal rojo ma. LG: matriz de log 2 de ( intensidades - background) para el canal verde LR = log 2 ( Rf -Rb ) LG =log 2 ( Gf - Gb ) A = 0. 5 x ( LR + LG ) M= LR - LG 39

Los datos se transforman aplicando log 2 para que las intensidades se distribuyan en Los datos se transforman aplicando log 2 para que las intensidades se distribuyan en forma aproximadamente simétrica. Esto mejora la visualización gráfica de los datos. También se intenta que se reduzca la relación entre la intensidad y la varianza que aparece cuando se realizan experimentos con replicaciones. La mayoría de los paquetes de análisis de microarrays utilizan log 2. 40

Gráficos espaciales del estado del microarray 41 Gráficos espaciales del estado del microarray 41

image(crudosexperimentospumilio[, 1], xvar= image(crudosexperimentospumilio[, 1], xvar="ma. M“) Generamos un gráfico de los valores de M del primer microarray ( [ , 1 ] ) de crudosexperimentospumilio: valores altos de M (cy 5 > cy 3) valores bajos de M (cy 5 < cy 3) 42

image(crudosexperimentospumilio[, 1], xvar= image(crudosexperimentospumilio[, 1], xvar="ma. Rb“) Generamos un gráfico de los valores del background para el canal rojo del primer microarray ( [ , 1 ] ) de crudosexperimentospumilio: valores altos de background en el canal rojo valores bajos de background en el canal rojo 43

par(mfrow=c(2, 2)) Generamos un espacio gráfico de 2 x 2 gráficos. image (crudosexperimentospumilio[, 1], par(mfrow=c(2, 2)) Generamos un espacio gráfico de 2 x 2 gráficos. image (crudosexperimentospumilio[, 1], xvar="ma. Spot. Col", bar=FALSE) Generamos un gráfico del vector columna para cada spot en el primer microarray ( [ , 1 ] ) de crudosexperimentospumilio ( bar=FALSE : para que no dibuje la barra de colores). image (crudosexperimentospumilio[, 1], xvar="ma. Print. Tip", bar=FALSE) Generamos un gráfico del vector print-tip para cada spot en el primer microarray ( [ , 1 ] ) de crudosexperimentospumilio image (crudosexperimentospumilio[, 1], xvar="ma. Controls", col=heat. colors(10), bar=FALSE) Generamos un gráfico espacial que indica en donde se encuentran los spots control en el primer array de ( [ , 1 ] ) de crudosexperimentospumilio, indicados en rojo image (crudosexperimentospumilio[, 1], xvar="ma. Plate", bar=FALSE) Generamos un gráfico de los valores de M del primer microarray ( [ , 1 ] ) de crudosexperimentospumilio 44

45 45

boxplot(crudosexperimentospumilio[, 1], xvar = boxplot(crudosexperimentospumilio[, 1], xvar = "ma. Print. Tip", yvar = "ma. M", main="arreglo 52987. gpr") Podemos generar un gráfico boxplot que nos muestre cómo varía M en función del printtip para el primer microarray ( [ , 1 ] ) de crudosexperimentospumilio: 46

boxplot(crudosexperimentospumilio, yvar= boxplot(crudosexperimentospumilio, yvar="ma. M") Podemos generar un gráfico boxplot que nos muestre cómo varía M en función del array para crudosexperimentospumilio: 47

library (array. Quality) Cargamos el paquete array. Quality. ma. Quality. Plots(crudosexperimentospumilio [, 1]) Generamos library (array. Quality) Cargamos el paquete array. Quality. ma. Quality. Plots(crudosexperimentospumilio [, 1]) Generamos un gráfico ma. Quality. Plots para evaluar a priori la calidad delmicroarray: 48

1. MA-plot de los datos crudos sin sustracción de background. Cada línea coloreada representa 1. MA-plot de los datos crudos sin sustracción de background. Cada línea coloreada representa la curva loess para cada grupo de print-tip. Estas curvas representan que cantidad de normalización se deberá realizar a cada grupo de datos. Mientras más separadas las curvas de la línea 0 de M, mayor normalización será necesaria. Mientras mayor la normalización necesaria, menor la calidad del experimento. 49

2. MA-plot usando el paquete hexbin que destaca la densidad de puntos, de los 2. MA-plot usando el paquete hexbin que destaca la densidad de puntos, de los datos normalizados; amarillo = alta densidad de puntos, azul = baja densidad. Este diagrama nos indica como quedó la normalización por grupo de print-tip (default). 50

3. Gráfico espacial de los valores de M crudos sin sustracción del fondo. Amarillo 3. Gráfico espacial de los valores de M crudos sin sustracción del fondo. Amarillo = valores bajos de M, azul = valores altos de M, blanco = valores faltantes. Esto nos permite visualizar una hibridación despareja. 51

4. Gráfico espacial de los rangos de los valores de M normalizados. Por defecto, 4. Gráfico espacial de los rangos de los valores de M normalizados. Por defecto, se utiliza la normalización loess por print-tip-group. Además, los puntos marcados (flagged) son señalados por medio de un recuadro negro. Este tipo de representación gráfica permite verificar que la normalización ha quitado efectos espaciales. 52

5. Gráfico espacial de los valores de A crudos sin sustracción del fondo. El 5. Gráfico espacial de los valores de A crudos sin sustracción del fondo. El color representa la intensidad de la señal. 53

6. Histogramas del cociente señal ruido (signal -to-noise ratio) (S 2 N) para cada 6. Histogramas del cociente señal ruido (signal -to-noise ratio) (S 2 N) para cada canal Cy 5 y Cy 3. S 2 N = log 2(intensidad del spot / intensidad del fondo) En la parte superior de cada histograma se imprimen la media y varianza de la señal. Además, se superponen curvas de densidad de SNR estratificado por diversos tipos de control (estado). Los esquemas de color de las curvas de densidad están dados la tabla 1. El SNR es un buen indicador para los problemas del tinte. Las líneas de densidad de controles negativos (por ej. secuencias artificiales que han sido diseñadas para no presentar homología con ningún genoma) y vacíos deben estar cercanas, casi superpuestas. 54

7. Diagrama de punto de los valores M normalizados para los controles. Los controles 7. Diagrama de punto de los valores M normalizados para los controles. Los controles con más de 3 réplicas se representan en el eje Y, el esquema de colores se representa en la tabla 1. Los valores de los controles M deben estar concentrados y cerca de 0. 55

8. Diagrama del punto de A, sin la sustracción del fondo. Los controles con 8. Diagrama del punto de A, sin la sustracción del fondo. Los controles con más de 3 réplicas se representan en el eje Y, el esquema de colores se representa en la tabla 1. La intensidad de los controles positivos debe estar en la región de alta intensidad, los controles negativos y vacíos deben estar en la región de intensidad más baja. Los rangos de los controles positivos y los negativos/vacíos deben estar separados. 56

Normalizaciónes Gráficos de normalización 57 Normalizaciónes Gráficos de normalización 57

Normalizaciones de datos : ¿ POR QUE NORMALIZAMOS? Para analizar los datos es necesario Normalizaciones de datos : ¿ POR QUE NORMALIZAMOS? Para analizar los datos es necesario reducir el error sistemático. Como suponemos que sólo una pequeña cantidad de genes se expresa diferencialmente, las distribuciones de intensidades deben ser uniformes espacialmente dentro del arreglo y entre arreglos. Entonces… normalizaremos primero dentro de un microarray utilizando la función ma. Norm y luego, si fuese necesario, entre arreglos con la función ma. Norm. Scale. ¿Cómo nos damos cuenta si es necesario? Luego de la normalización dentro de cada microarray, vamos a realizar diferentes gráficos que nos van a indicar si la escala de los datos entre microarrays es diferente o no. Si es muy diferente, tendremos que normalizar por escala. Recuerden que por cada paso de normalización estamos perdiendo información. 58

Normalización dentro de cada microarray: experimentospumilionorm <- ma. Norm(crudosexperimentospumilio, norm= Normalización dentro de cada microarray: experimentospumilionorm <- ma. Norm(crudosexperimentospumilio, norm= "p") Generamos un objeto que contiene los datos normalizados por grupo de printtip ( norm = “p” : normaliza por intensidad DENTRO de cada aguja ) de crudosexperimentospumilio. summary (experimentospumilionorm) Nos muestra un resumen de las características del objeto experimentospumilionorm que generamos. summary (crudosexperimentospumilio) Nos muestra nuevamente un resumen de las características del objeto crudosexperimentospumilio. 59

Normalización dentro de cada microarray: Podemos comparar la parte estadística de crudosexperimentospumilio (sin normalizar) Normalización dentro de cada microarray: Podemos comparar la parte estadística de crudosexperimentospumilio (sin normalizar) con la misma de experimentospumilionorm ( mismos datos normalizados por print-tip group): antes después Las medias y las medianas quedaron cerca de cero en todos los arrays del experimento!!!! 60

Normalización dentro de cada microarray: veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un espacio gráfico Normalización dentro de cada microarray: veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un espacio gráfico de 2 gráficos independientes. boxplot(crudosexperimentospumilio[, 1], xvar = "ma. Print. Tip", yvar = "ma. M", main="crudosexperimentospumilio") Generamos un gráfico boxplot MA ( yvar = "ma. M" ) de los datos sin normalizar de crudosexperimentospumilio, divididos por grupo de print-tip ( xvar = "ma. Print. Tip“ ) del primer microarray ( [ , 1 ] ). boxplot(experimentospumilionorm[, 1], xvar = "ma. Print. Tip", yvar = "ma. M", main="experimentospumilionorm") Generamos un gráfico boxplot MA ( yvar = "ma. M" ) de los datos normalizados por print-tip de experimentospumilionorm, divididos por grupo de print-tip ( xvar = "ma. Print. Tip“ ) del primer microarray ( [ , 1 ] ). 61

Normalización dentro de cada microarray: veámoslo en gráficos Antes de normalizar por Print-tip Hay Normalización dentro de cada microarray: veámoslo en gráficos Antes de normalizar por Print-tip Hay variaciones entre las medias o medianas, en los spoteados por diferentes print-tips. La mayoría de los valores está debajo de cero Luego de normalizar por Print-tip Transformamos los datos de manera que las medias o medianas de cada print tip son iguales ( 0) !!!! 62

Normalización dentro de cada microarray: veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un espacio gráfico Normalización dentro de cada microarray: veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un espacio gráfico de 2 gráficos independientes. boxplot(crudosexperimentospumilio, yvar = "ma. M", main = "pre normalización") Generamos un gráfico boxplot MA ( yvar = "ma. M" ) de los datos sin normalizar de crudosexperimentospumilio, dividos por microarray. boxplot(experimentospumilionorm, yvar = "ma. M", main = "post normalización") Generamos un gráfico boxplot MA ( yvar = "ma. M" ) de los datos normalizados de experimentospumilionorm, divididos por microarray. 63

Normalización dentro de cada microarray: veámoslo en gráficos Antes de normalizar Hay variaciones entre Normalización dentro de cada microarray: veámoslo en gráficos Antes de normalizar Hay variaciones entre las medias o medianas en los datos de los diferentes microarrays. La mayoría de los valores está debajo de cero Luego de normalizar Transformamos los datos de manera que las medias o medianas de todos los microarrays son similares ( 0) !!!! 64

Normalización dentro de cada microarray: veámoslo en gráficos ( MA plot ) par(mfrow=c(2, 1)) Normalización dentro de cada microarray: veámoslo en gráficos ( MA plot ) par(mfrow=c(2, 1)) Generamos un espacio de 2 gráficos. ma. Plot(crudosexperimentospumilio, legend. func = "NULL")Generamos un gráfico MA plot de los datos sin normalizar de crudosexperimentospumilio; las líneas de colores representan las curvas loess para cada print tip group. ma. Plot(experimentospumilionorm, legend. func = "NULL") Generamos un gráfico MA plot de los datos normalizados de experimentospumilionorm; las líneas de colores representan las curvas loess para cada print tip group. 65

Normalización dentro de cada microarray: veámoslo en gráficos ( MA plot ) Antes de Normalización dentro de cada microarray: veámoslo en gráficos ( MA plot ) Antes de normalizar Luego de normalizar 66

Normalización entre microarrays ( por escala ) : experimentospumilionormescala <- ma. Norm. Scale(experimentospumilionorm) Generamos Normalización entre microarrays ( por escala ) : experimentospumilionormescala <- ma. Norm. Scale(experimentospumilionorm) Generamos un objeto que contiene los datos normalizados por escala ( función ma. Norm. Scale ) a partir de objeto de datos prenormalizados por printtip group experimentospumilionorm. par(mfrow=c(2, 1)) Generamos un espacio gráfico de 2 gráficos. boxplot(experimentospumilionorm, yvar = "ma. M", main = "pre normalización") Generamos un gráfico boxplot de M ( yvar = "ma. M" ) de los datos sin normalizar por escala de experimentospumilionorm, dividos por microarray. boxplot(experimentospumilionormescala, yvar = "ma. M", main = "post normalización") Generamos un gráfico boxplot de M ( yvar = "ma. M" ) de los datos normalizados por escala de experimentospumilionormescala, divididos por microarray. 67

Normalización entre microarrays ( por escala ): veámoslo en gráficos Antes de normalizar por Normalización entre microarrays ( por escala ): veámoslo en gráficos Antes de normalizar por escala Hay variaciones entre las distribuciones de los datos de los diferentes microarrays Luego de normalizar por escala, las cajas de los box-plots tienen anchos casi idénticos!!!! 68

Normalización entre microarrays ( por cuantiles ) : library(convert) Necesitamos cargar el paquete convert Normalización entre microarrays ( por cuantiles ) : library(convert) Necesitamos cargar el paquete convert para transformar los objetos que obtuvimos en el paquete marray ( crudosexperimentospumilio, experimentospumilionorm ) a objetos legibles en el paquete limma, sobre el cuál realizaremos las normalizaciones por cuantiles. crudosexperimentospumilio. limma <- as(crudosexperimentospumilio, "RGList") Creamos un objeto crudosexperimentospumilio. limma que contiene los datos del objeto crudosexperimentospumilio de marray. experimentospumilionorm. limma <- as(experimentospumilionorm, "MAList") Creamos un objeto experimentospumilionorm. limma que contiene los datos del objeto experimentospumilionorm de marray. experimentospumilionormescala. limma

Normalización entre microarrays ( por cuantiles ): veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un Normalización entre microarrays ( por cuantiles ): veámoslo en gráficos par(mfrow=c(2, 1)) Generamos un espacio gráfico de 2 gráficos independientes. boxplot(experimentospumilionormescala. limma$M~col(experimentospumilionorm escala. limma$M), names=colnames(experimentospumilionormescala. limma$M), main="norm. por escala", col="red") Generamos un gráfico boxplot de M de los datos normalizados por escala experimentospumilionormescala. limma, dividos por microarray. boxplot(experimentospumilionormporquantiles. limma$M~col(experimentospumili onormporquantiles. limma$M), names=colnames(experimentospumilionormporqua ntiles. limma$M), main="norm. por quantiles", col="red") Generamos un gráfico boxplot de M de los datos normalizados por escala experimentospumilionormporquantiles. limma, dividos por microarray. 70

Normalización entre microarrays ( por cuantiles ): veámoslo en gráficos Antes de normalizar por Normalización entre microarrays ( por cuantiles ): veámoslo en gráficos Antes de normalizar por cuantiles Hay variaciones entre las distribuciones de los datos de intensidad de los diferentes microarrays Luego de normalizar por cuantiles Las diferencias entre las distribuciones se reducen!!!! 71

Podemos también normalizar enteramente en el paquete limma: Para ello generaremos un archivo. txt Podemos también normalizar enteramente en el paquete limma: Para ello generaremos un archivo. txt con información acerca de los experimentos que estamos analizando, que tenga la siguiente estructura de datos: Slide. Number 075 -dm 0603 B 2 079 -dm 0603 B 2 085 -dm 0603 B 2 097 -dm 0603 B 2 File. Name 52987. gpr 53530. gpr 54253. gpr 54596. gpr Cy 3 RNA total Cy 5 RNA IP Date 2004/06/11 2004/07/06 2004/07/27 2004/08/09 Llamaremos a este archivo Pumilio. Sample. txt , más tarde nos servirá para poder trabajar con los datos. 72

Normalización entre microarrays ( en limma ): library(limma) Cargamos el paquete limma para trabajar Normalización entre microarrays ( en limma ): library(limma) Cargamos el paquete limma para trabajar sobre el mismo. targetspumilio<-read. Targets("Pumilio. Sample. txt") Generamos un objeto targetspumilio que es un dataframe, generado por la función read. Targets, con los datos de Pumilio. Sample. txt. f <- function(x) as. numeric(x$Flags > -99) Generamos una función f que se utilizará para fijar un filtro de manera que cualquier spot marcado con ( flag de) -99 o menos tenga peso cero. RGpumilio <- read. maimages(targetspumilio$File. Name, source="genepix", wt. fun=f) Los nombres de los archivos guardados en targetspumilio permiten identificar la información de intensidades contenida en los archivos de salida del procesamiento de imágenes, en este caso son archivos con extensión. gpr 73

experimentospumiliocorregidospor. RMA <- background. Correct(RGpumilio, method= experimentospumiliocorregidospor. RMA <- background. Correct(RGpumilio, method="rma") Generamos un objeto experimentospumiliocorregidospor. RMA que contiene los datos de la corrección por “rma”, que es la recomendada por la guía de usuarios de limma para los datos de Gene. Pix. experimentospumilionormalizadosporprinttip

Normalización entre microarrays ( en limma ): veámoslo en gráficos par(mfrow=c(2, 2)) Generamos espacio Normalización entre microarrays ( en limma ): veámoslo en gráficos par(mfrow=c(2, 2)) Generamos espacio de 2 x 2 gráficos. plot. Densities(experimentospumiliocorregidospor. RMA) text(15, 0. 8, "RMA") Generamos un gráfico de las densidades en función de la intensidad de experimentospumiliocorregidospor. RMA. plot. Densities(experimentospumilionormalizadosporprinttip) text(15, 0. 45, "print tip") Generamos un gráfico de las densidades en función de la intensidad de experimentospumilionormalizadosporprinttip. plot. Densities(experimentospumilionormalizadosporprinttipycuantiles) text(15, 0. 25, "cuantiles") Generamos un gráfico de las densidades en función de la intensidad de experimentospumilionormalizadosporprinttipycuantiles. plot. Densities(experimentospumilionormalizadosporvsn) text(15, 0. 12, "vsn") Generamos un gráfico de las densidades en función de la intensidad de experimentospumilionormalizadosporprinttipyvsn. 75

Nos quedamos con esta normalización: Print-tip + cuantiles 76 Nos quedamos con esta normalización: Print-tip + cuantiles 76

Identificación de genes diferencialmente expresados Gráficos y tablas de genes diferencialmente expresados 77 Identificación de genes diferencialmente expresados Gráficos y tablas de genes diferencialmente expresados 77

OK, ahora que ya tenemos los datos normalizados, necesitamos saber en que canales ( OK, ahora que ya tenemos los datos normalizados, necesitamos saber en que canales ( Cy 3 – Cy 5 ) está hibridada cuál de las dos muestras ( RNA total – RNA inmunoprecipitado con pumilio ) en cada uno de los microarrays que analizamos, clickeamos el botón view de la página de Stanford de la cuál bajamos nuestros datos: 78

RNA total Cy 3 = canal 1 RNA IP Cy 5 = canal 2 RNA total Cy 3 = canal 1 RNA IP Cy 5 = canal 2 79

Exp ID Cy 3 Cy 5 52987 RNA total RNA IP 53530 RNA total Exp ID Cy 3 Cy 5 52987 RNA total RNA IP 53530 RNA total RNA IP 54253 RNA total RNA IP 54596 RNA total RNA IP Aquí observamos un gran error que hicieron nuestros amigos de Stanford en cuanto al diseño experimental: NO hicieron dye swapping en este experimento, RECUERDEN REALIZAR DYE SWAPPING EN SUS EXPERIMENTOS!!!!! 80

Podemos adoptar el siguiente modelo lineal para nuestro experimento: o, equivalentemente: para un dye-swap: Podemos adoptar el siguiente modelo lineal para nuestro experimento: o, equivalentemente: para un dye-swap: matriz de diseño Cy 3 La matriz de diseño para un experimento con dye-swap sería: Cy 5 RNA IP RNA total 1 RNA total RNA IP -1 81

Si asignamos β = E(log 2 RNA total. Cy 5/RNA IPCy 3) nos queda Si asignamos β = E(log 2 RNA total. Cy 5/RNA IPCy 3) nos queda una grilla de este tipo: -M M Exp ID Cy 3/Cy 5/Cy 3 52987 -1 1 53530 -1 1 54253 -1 1 54596 -1 1 Con estos datos creamos una matriz de diseño del experimento que será la siguiente : design = c(-1, -1, -1) para estimar la relación RNA_total/RNA_IP (el inverso de M, ya que aquí marcamos con Cy 5 a RNA IP). 82

pumiliofit 1 <-lm. Fit(experimentospumilionormalizadosporprinttipycuantiles, design=c(-1, -1, -1)) Generamos un objeto pumiliofit 1 que contiene pumiliofit 1 <-lm. Fit(experimentospumilionormalizadosporprinttipycuantiles, design=c(-1, -1, -1)) Generamos un objeto pumiliofit 1 que contiene los datos de M ajustados de cada gen (previamente normalizados) con la función lm. Fit. Ésta, por defecto realiza un ajuste por cuadrados mínimos ordinario. Tiene opciones para incluir una estructura de covarianzas de los errores y un ajuste robusto. La matriz design indica cuáles arreglos son dye-swaps. pumiliofit <- e. Bayes(pumiliofit 1) La función e. Bayes calcula el estadístico t moderado (Lönnstedt and Speed 2001) y un p-valor, para cada gen. summary(pumiliofit) Nos muestra las características de pumiliofit. 83

Tablas de genes diferencialmente expresados: tablatop 10 Bpumilio <top. Table(pumiliofit, coef=1, number=10, genelist=pumiliofit$genes, adjust. Tablas de genes diferencialmente expresados: tablatop 10 Bpumilio

85 85

tablatop 17000 Bpumilio <top. Table(pumiliofit, coef=1, number=17000, genelist=pumiliofit$genes, adjust. method= tablatop 17000 Bpumilio

87 87

Puedo ver en un gráfico MA cuáles son los genes top score de expresión Puedo ver en un gráfico MA cuáles son los genes top score de expresión diferencial, con respecto a diferentes criterios de selección ( M, t o B): M<-pumiliofit$coefficients Generamos un objeto M que contenga los datos del espacio coefficients ( valores M ) de pumiliofit. A<-pumiliofit$Amean Generamos un objeto A que contenga los datos del espacio Amean ( valores A ) de pumiliofit. plot(A, M, pch=". ", col="lightblue", cex=3, main="Valores de M ") Generamos un gráfico de M vs. A. points(tablatop 10 Bpumilio$A[1], tablatop 10 Bpumilio$M[1], pch=“? ", col="blue", cex=2) ; text(12. 27, 4. 738, "CG 2229") Generamos un punto en el gráfico anterior, que representa al gen con mejor score de tablatop 10 Bpumilio, ubicándolo espacialmente por su A ( 12. 27 ) y M ( 4. 738 ). points(tablatop 10 Mpumilio$A[1], tablatop 10 Mpumilio$M[1], pch=“? ", col="blue", cex=2 ); text(9. 206855, -6. 772648, "CG 7979") Generamos un punto en el gráfico anterior, que representa al gen con mejor score de tablatop 10 Mpumilio, ubicándolo espacialmente por su A (9. 206855 ) y M (- 88 6. 772648 ).

89 89

También puedo buscar que datos estadísticos de expresión diferencial presenta un gen x del También puedo buscar que datos estadísticos de expresión diferencial presenta un gen x del cuál conozco su nombre o ID: genquemeinteresa

También podemos intersectar los 100 genes con mayores valores de M, B y t, También podemos intersectar los 100 genes con mayores valores de M, B y t, los cuáles van a ser los mayores candidatos a estar diferencialmente expresados: ordinary. t <- pumiliofit$coef / pumiliofit$stdev. unscaled / pumiliofit$sigma Generamos un objeto que contiene los datos de coef, stdev. unscaled y sigma de pumiliofit. ord. M <- order(abs(pumiliofit$coefficients), decreasing=TRUE) Generamos un objeto que contiene los datos ordenados de M (coefficients) en orden decreciente (valores absolutos). ord. B <- order(pumiliofit$lods, decreasing=TRUE) Generamos un objeto que contiene los datos ordenados de B (lods) en orden decreciente. ordt <- order(abs(ordinary. t), decreasing=TRUE) Generamos un objeto que contiene los datos ordenados de t (ordinary. t) en orden decreciente (valores absolutos). vector <- rep(0, 17664) Generamos un vector de tantos valores como spots tengael array ( 0 a 17664 en éste caso). top 100 M <- ord. M[1: 100] Generamos un objeto que contenga los 100 primeros valores de ord. M. top 100 B <- ord. B[1: 100] Generamos un objeto que contenga los 100 primeros valores de ord. B. 91 top 100 t <- ordt[1: 100] Generamos un objeto que contenga los 100 primeros valores de ordt.

M <- vector Asigno a M el objeto vector. B <- vector Asigno a M <- vector Asigno a M el objeto vector. B <- vector Asigno a B el objeto vector. t <- vector Asigno a t el objeto vector. B[top 100 B]<-1 A los elementos de top 100 B que estén en el vector B, les asigno un valor de 1. t[top 100 t]<-1 A los elementos de top 100 t que estén en el vector t, les asigno un valor de 1. M[top 100 M]<-1 A los elementos de top 100 M que estén en el vector M, les asigno un valor de 1. interseccion <- (B*t*M)/(B*t*M) Genero un objeto interseccion que contiene los datos resultantes de dicha operación entre vectores. interseccion Visualizamos intersección. 92

gen 10020 gen 10328 93 gen 10020 gen 10328 93

¿Cómo se qué gen es el 10020? Names<-pumiliofit$genes$Name Generamos un objeto Names que contiene ¿Cómo se qué gen es el 10020? Names<-pumiliofit$genes$Name Generamos un objeto Names que contiene los datos de nombres de genes de pumiliofit. diferencial 10020<-Names[10020] Generamos un objeto que contenga el nombre del gen 10020. diferencial 10020 Visualizamos a diferencial 10020. 94

¿Cuáles son los genes de la intersección entre t, B y M? : pumiliointerseccion ¿Cuáles son los genes de la intersección entre t, B y M? : pumiliointerseccion 2 <- data. frame (interseccion, Names) Generamos un objeto pumiliointerseccion 2 que contiene dos columnas de datos, en la primera interseccion y en la segunda names. subset (pumiliointerseccion 2 , interseccion = =1) Visualizamos los elementos de pumiliointerseccion 2 que tienen un valor de la columna correspondiente a interseccion = 1 ( recordar que los elementos de intersección entre t, B y M tienen valor =1). 95

subset (pumiliointerseccion 2, Names==“Vha 16”) Exploramos los elementos de pumiliointerseccion 2 que contienen datos subset (pumiliointerseccion 2, Names==“Vha 16”) Exploramos los elementos de pumiliointerseccion 2 que contienen datos acerca del gen de nombre “Vha 16”. Existen varios spots del gen Vha 16, uno de los cuáles está en la intersección de los estadísticos t, B y M. 96

Graficamos los genes de la intersección sobre MA y t. A: M 1 <- Graficamos los genes de la intersección sobre MA y t. A: M 1 <- pumiliofit$coefficients Generamos un objeto M 1 que contenga los datos de M (coefficients) de pumiliofit. A <- pumiliofit$Amean Generamos un objeto A que contenga los datos de A (Amean) de pumiliofit. par(mfrow=c(1, 2)) Generamos un espacio de 2 gráficos. plot(A, M 1, pch='. ', col="green") Generamos un gráfico de MA con los objetos M 1 y A. points(A*interseccion, M 1*interseccion, pch="*") Generamos puntos * en los datos correspondientes a la intersección. points(tablatop 10 Mpumilio$A[1], tablatop 10 Mpumilio$M[1], pch="*", col="blue", ce x=2); text(9. 206855, -6. 772648, "CG 7979“) Generamos un punto * azul en el dato correspondiente al primer gen de tablatop 10 Mpumilio. plot(A, ordinary. t, ylab="t", pch=". ", col="green") Generamos un gráfico de t. A. points(A*interseccion, ordinary. t*interseccion , pch="*") Generamos puntos * en los datos correspondientes a la intersección. 97

points(tablatop 10 tpumilio$A[1], tablatop 10 tpumilio$ordinary. t[1], pch= points(tablatop 10 tpumilio$A[1], tablatop 10 tpumilio$ordinary. t[1], pch="*", col="blue", cex=2); text(12. 267057, 4. 738344, "CG 2229") Generamos un punto * azul en el dato correspondiente al primer gen de tablatop 10 tpumilio. 98

Gracias por venir!!!!!!!!! Está claro que ahora no van a tener más dudas……jajaja………. pero Gracias por venir!!!!!!!!! Está claro que ahora no van a tener más dudas……jajaja………. pero en caso de que las tengan……………. http: //www. dm. uba. ar/materias/analisis_expl_y_conf_de_datos_de_exp_de _marrays_Mae/2006/1/ Página de la materia Análisis Exploratorio y Confirmatorio de Datos de Experimentos de Microarrays dictada por la Dra. Diana Kelmansky. 99