El problema

Planteamiento

Estudiamos si hay relación entre conjuntos de genes previamente definidos y fenotipo.
Estos grupos vendrán definidos según distintos criterios.
Por ejemplo, corresponder a una ruta metabólica,
localizados en un mismo cromosoma o bien definidos utilizando términos de Gene Ontology o KEGG.
De un modo genérico: un conjunto de genes que represente un un hecho biológico relevante.

Ejemplos

Un ejemplo simulado de Efron y Tibshirani: et1

Se consideran 1000 genes y 50 muestras.
Se supone que las 25 primeras muestras son controles y las últimas 25 es el grupo de tratamiento.
Definimos los grupos de genes como:
- las 20 primeras filas (de la matriz de expresión) corresponden al primer grupo, de la 21 a la 40 el segundo y así sucesivamente.
Los niveles de expresión los generamos aleatoriamente independientemente y con la misma distribución.
La distribución común es una normal estándar, $Y_{ij} \sim N(0,1)$ .
En el primer grupo añadimos un valor constante de 2.5 a los primeros 10 genes en las muestras correspondientes a tratamiento (últimas 25 columnas de la matriz de expresión).
En consecuencia, lo que hacemos es que solamente el primer grupo tiene la mitad de los genes asociados con el fenotipo y la otra mitad sin ningún tipo de asociación. -El resto de grupos no tiene ninguna asociación con fenotipo.

 N = 1000; n = 50
 set.seed(280562) ## Para obtener los mismos valores generados
 et1  = matrix(rnorm(N*n),nrow = N,ncol = n)
 et1[1:10,26:50] = et1[1:10,26:50] + 2.5
 y.et = factor(rep(1:2,rep(25,2)),levels = 1:2, 
               labels = c("Normal","Tratamiento"))

Cuantificamos una medida de asociación fenotipo-expresión.

et1.tt = genefilter::rowttests(et1,y.et)$statistic

En figura tenemos una estimación de la densidad. Podemos ver cómo se aprecia, por debajo de 8, el valor que corresponde al primer grupo.

pacman::p_load("ggplot2")
df = data.frame(et1.tt)
(p = ggplot(df,aes(x=et1.tt))+geom_density())

Enriqueciendo el conjunto de genes

El vector $t = (t_1,\ldots,t_N)$ constituye una descripción de la asociación marginal gen-fenotipo.
Enriquecer el estadístico $t$ para el conjunto consiste en considerar un resumen de $t$ sobre $S$ .
Si denotamos el estadístico de enriquecimiento por
1. $t(S)$ algunas opciones son $t(S) = \sum_{i \in S} \frac{t_i}{n_S}.$
2. Definimos:
  - $t^{+} = \max\{t,0\}$
  - $t^{-} = - \min\{t,0\}$
  - $\bar{t}^{+}_S = \sum_{i \in S} \frac{t^{+}_i}{n_S}$
  - $\bar{t}^{-}_S = \sum_{i \in S} \frac{t^{-}_i}{n_{S}}$ La medida de enriquecimiento maxmean es $t(S) = \max\{\bar{t}^{+}_S, \bar{t}^{-}_S\}$ Estamos tomando la media de las partes positivas $t_i^+$ , la media de las partes negativas $t_i^-$ y nos quedamos con el máximo de ambos valores.

Aleatorizando columnas (muestras)

Para un gen tenemos su perfil de expresión $u = (u_1,\ldots,u_n)$ y tenemos el vector $y$ asociado a las muestras.
Sea $\pi$ denota una permutación aleatoria de $(1,\ldots,n)$ de forma que los $\pi(i)$ será el índice que ocupa la posición $i$ -ésima en la permutación $\pi$ .
Denotaremos por $\pi_0$ la permutación que nos devuelve el orden original: $\pi_0(j) = j$ .
Si el vector asociado a las muestras es $y$ entonces tendremos el vector (permutado de $y$ ) $y_{\pi} = (y_{\pi(1)},\ldots,y_{\pi(n)})$ .
Obviamente, $y = y_{\pi_0}$ .
El valor observado de la asociación gen-fenotipo para $u$ e $y$ será $t_{\pi_0}$ .
Si consideramos $B$ permutaciones aleatorias entonces tendremos los valores $t_{\pi_1},\ldots,t_{\pi_B}$ para las $B$ permutaciones aleatorias.
Si no hay asociación gen-fenotipo entonces el valor de $t_{\pi_0}$ debe de ser como los valores $t_{\pi_1},\ldots,t_{\pi_B}$ .
De hecho, cualquier ordenación de $t_{\pi_0},t_{\pi_1},\ldots,t_{\pi_B}$
tiene la misma probabilidad.
El p-valor sería la proporción de $t_{\pi_b}$ s mayores que $t_{\pi_0}$ , es decir, $p_r = \frac{|\{b: t_{\pi_b} > t_{\pi_0}\}|}{B+1}.$
En el caso en que una mayor asociación se exprese como un valor muy grande (positivo) o muy pequeño (negativo) de $t$ entonces el p valor vendría dado por $p_r = \frac{|\{b: |t_{\pi_b}| > |t_{\pi_0}|\}|}{B+1},$ ya que tendríamos un test bilateral.

Aleatorizando filas (genes o características)

Denotemos por $S_0$ el conjunto de genes original en el que tenemos interés.
Elegimos al azar un grupo $S$ de genes (filas) del mismo cardinal que $S_0$ .
Lo aleatorio ahora no es el vector $y$ sino el conjunto $S$ de genes.
Si $S$ es aleatorio entonces también lo es $t(S)$ (ahora las muestras tienen su ordenación original).
A la distribución de probabilidad de $t(S)$ se le llama distribución de aleatorización del estadístico de enriquecimiento.
Tomamos $B$ selecciones aleatorias de $n_{S_0}$ genes.
Tenemos los conjuntos $S_b$ con $b= 1, \ldots, B$ y los valores del estadístico de enriquecimiento observados son $t(S_0)$ (para el grupo original) y $t(S_b)$ con $b= 1, \ldots, B$ para los seleccionados al azar.
El p-valor se define análogamente como $p_c = \frac{|\{b: t(S_b) > t(S_0)\}|}{B+1}.$ si solamente rechazamos para valores grandes (positivos).
En el caso bilateral tendremos $p_{ae} = \frac{|\{b: |t(S_b)| > |t(S_0)|\}|}{B+1}.$

Análisis de conjuntos de genes

El problema

Planteamiento

Diferencias con la aproximación gen-a-gen o marginal

Conjunto(s) de genes

¿Hay asociación entre un conjunto de genes y el fenotipo?

Ejemplos

Un ejemplo simulado de Efron y Tibshirani: et1

Un ejemplo simulado de Efron y Tibshirani: et2

Grupos

Enriquecimiento

Cuantificando asociación gen-fenotipo

Enriqueciendo el conjunto de genes

Distribuciones condicionadas a los datos

Aleatorizando columnas (muestras)

Aleatorizando filas (genes o características)

¿Cómo hacerlo?

Realizando un análisis de grupos de genes

Limma::genSetTest

Limma::wilcoxGST

Colecciones de grupos de genes

Gene Ontology

Explorando las colecciones de genes

Using `lapply` and `sapply`

Análisis de gse21942 con paquete tami

Cargamos paquetes y datos

Un primer análisis de gse21942

Un segundo análisis

Informes con los resultados

`PRJNA297664`