-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy path02-notas.Rmd
205 lines (143 loc) · 6.62 KB
/
02-notas.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
Universidad Nacional Autónoma - Licenciatura en ciencias genómicas
Bioinformática y estadística 2 - Análisis de datos de secuenciación masiva (Dr. Leonardo Collado)
**Salvador González Juárez**
# Notas 24 de marzo
## Objetos de Bioconductor para datos de expresion
### Paquete rnaseqGene
- URL: [https://bioconductor.org/packages/release/workflows/html/rnaseqGene.html](https://bioconductor.org/packages/release/workflows/html/rnaseqGene.html)
Bioconductor tiene muchos paquetes que admiten el análisis de datos de secuencia de alto rendimiento, incluida la secuenciación de ARN (RNA-seq). Los paquetes que usaremos en este flujo de trabajo incluyen paquetes centrales mantenidos por el equipo central de Bioconductor para trabajar con anotaciones de genes (ubicaciones de genes y transcripciones en el genoma, así como búsqueda de ID de genes). También utilizaremos paquetes contribuidos para el análisis estadístico y la visualización de datos de secuenciación.
- SummarizedExperiment: Matriz donde se puede visualizar las muestras (columnas) de multiples genes (filas). Tambien puede tener multiples tablas e información del experimento.
- GRanges: Estructura para almacenar de forma eficiente la información de los genes. Es obtenido de una dependencia con el mismo nombre: GenomicRanges ([https://bioconductor.org/packages/release/bioc/html/GenomicRanges.html](https://bioconductor.org/packages/release/bioc/html/GenomicRanges.html)).
- tracklayer: Útil para leer archivos de texto de diferentes formatos, con enfasis bioinformático.
---
# SummarizedExperiment
```{r}
## Lets build our first SummarizedExperiment object
library("SummarizedExperiment")
## ?SummarizedExperiment
## De los ejemplos en la ayuda oficial
## Creamos los datos para nuestro objeto de tipo SummarizedExperiment
## para 200 genes a lo largo de 6 muestras
nrows <- 200
ncols <- 6
## Números al azar de cuentas
set.seed(20210223)
counts <- matrix(runif(nrows * ncols, 1, 1e4), nrows)
## Información de nuestros genes
rowRanges <- GRanges(
rep(c("chr1", "chr2"), c(50, 150)),
IRanges(floor(runif(200, 1e5, 1e6)), width = 100),
strand = sample(c("+", "-"), 200, TRUE),
feature_id = sprintf("ID%03d", 1:200)
)
names(rowRanges) <- paste0("gene_", seq_len(length(rowRanges)))
## Información de nuestras muestras
colData <- DataFrame(
Treatment = rep(c("ChIP", "Input"), 3),
row.names = LETTERS[1:6]
)
## Juntamos ahora toda la información en un solo objeto de R
rse <- SummarizedExperiment(
assays = SimpleList(counts = counts),
rowRanges = rowRanges,
colData = colData
)
## Exploremos el objeto resultante
rse
## Número de genes y muestras
dim(rse)
## IDs de nuestros genes y muestras
dimnames(rse)
## Nombres de tablas de cuentas que tenemos (RPKM, CPM, counts, logcounts, etc)
assayNames(rse)
## El inicio de nuestra tabla de cuentas
head(assay(rse))
## Información de los genes en un objeto de Bioconductor
rowRanges(rse)
## Tabla con información de los genes
rowData(rse) # es idéntico a 'mcols(rowRanges(rse))'
## Tabla con información de las muestras
colData(rse)
## Explora el objeto rse de forma interactiva
library("iSEE")
#iSEE::iSEE(rse)
```
## ¿Como utilizar el objeto SummarizedExperiment?
**class** → Nombre del objeto.
**dim** → Dimesiones de la matriz.
**assay** → Son los datos de un experimento, pueden haber varias tablas con distintos datos.
**rownames** → el nombre de los genes.
**rowData names** → nombre general de las filas.
**colnames** → el nombre de las muestras.
**colData names** → nombre general de las columnas (muestras)
---
## Ejercicio: Explica que sucede en las siguientes líneas de código de R.
**Comando 1**
```{R}
rse[1:2, ]
```
- Obtenemos los primeros dos genes de todas las muestras.
**Comando 2**
```{R}
rse[, c("A", "D", "F")]
```
- Obtenemos todos los genes de las muestras A, D y F.
---
## Ejercicio con spatialLIBD
**Descarga un PDF que reproduzca la imagen de la diapositiva. Incluye ese PDF en tu repositorio de notas del curso.**
```{r}
## Descarguemos unos datos de spatialLIBD
sce_layer <- spatialLIBD::fetch_data("sce_layer")
## Revisemos el tamaño de este objeto
pryr::object_size(sce_layer)
#iSEE::iSEE(sce_layer)
```
- La imagen se encuentra en el archivo *ReducedDimensionPlot1.pdf*
**Explora en con un heatmap la expresión de los genes MOBP, MBP y PCP4. Si hacemos un clustering (agrupamos los genes), ¿cúales genes se parecen más?**
- Son mas parecidos MOBP y MBP, como podemos observar en la imagen de *ComplexHeatmapPlot2.pdf*
**¿En qué capas se expresan más los genes MOBP y MBP?**
- Ambos genes tienen más presencia en las capas WM, L6 y un poco en L1, como podemos observar en la imagen de *ComplexHeatmapPlot1.pdf*
# Datos de RNA-seq a través de recount3
**Proyecto recunt**: precesar datos de RNA-seq publicos para facilitar el acceso y uso de la información en análisis trancriptómicos. Actualmente solo hay datos de humano y ratón.
```{r}
## Load recount3 R package
library("recount3")
## Revisemos todos los proyectos con datos de humano en recount3
human_projects <- available_projects()
## Encuentra tu proyecto de interés. Aquí usaremos
## SRP009615 de ejemplo
proj_info <- subset(
human_projects,
project == "SRP009615" & project_type == "data_sources"
)
## Crea un objetio de tipo RangedSummarizedExperiment (RSE)
## con la información a nivel de genes
rse_gene_SRP009615 <- create_rse(proj_info)
## Explora el objeto RSE
rse_gene_SRP009615
## Explora los proyectos disponibles de forma interactiva
#proj_info_interactive <- interactiveDisplayBase::display(human_projects)
## Selecciona un solo renglón en la tabla y da click en "send".
## Aquí verificamos que solo seleccionaste un solo renglón.
#stopifnot(nrow(proj_info_interactive) == 1)
## Crea el objeto RSE
#rse_gene_interactive <- create_rse(proj_info_interactive)
## Convirtamos las cuentas por nucleotido a cuentas por lectura
## usando compute_read_counts().
## Para otras transformaciones como RPKM y TPM, revisa transform_counts().
assay(rse_gene_SRP009615, "counts") <- compute_read_counts(rse_gene_SRP009615)
## Para este estudio en específico, hagamos más fácil de usar la
## información del experimento
rse_gene_SRP009615 <- expand_sra_attributes(rse_gene_SRP009615)
colData(rse_gene_SRP009615)[
,
grepl("^sra_attribute", colnames(colData(rse_gene_SRP009615)))
]
#iSEE::iSEE(rse_gene_SRP009615)
```
# Ejercicio: Utiliza iSEE para reproducir la imagen
**Pistas**:
- Utiliza el dynamic feature selection
- Utiliza información de las columnas para el eje X
- Utiliza información de las columnas para los colores
- La imagen se encuentra en el archivo *FeatureAssayPlot1.pdf*