Cálculo de la probabilidad de superposición de la lista de genes entre una secuencia de ARN y un conjunto de datos de chip ChIP

Espero que alguien en estos foros pueda ayudarme con este problema básico en los estudios de expresión génica.

Hice una secuenciación profunda de un tejido experimental y uno de control. Luego obtuve valores de enriquecimiento doble de genes en la muestra experimental sobre el control. El genoma de referencia tiene ~ 15,000 genes. 3.000 de los 15.000 genes están enriquecidos por encima de un cierto límite en mi muestra de interés en comparación con el control.

Entonces: A = población total de genes = 15,000 B = subpoblación enriquecida en RNA-Seq = 3,000.

En un experimento anterior de ChIP-chip, encontré 400 genes enriquecidos por ChIP-chip. De los 400 genes del chip ChIP, 100 genes están en el grupo de 3.000 transcripciones de RNA-Seq enriquecidas.

Entonces: C = número total de genes enriquecidos en chips ChIP = 400.

¿Cuál es la probabilidad de que mis 100 genes ChIP-chip se enriquezcan con RNA-Seq solo por casualidad? En otras palabras, ¿cuál es la forma más prudente de calcular si mi superposición observada entre B y C (100 genes) es mejor que la obtenida por casualidad? Por lo que he leído hasta ahora, la mejor manera de probar esto es mediante el uso de distribución hipergeométrica.

Utilicé una calculadora en línea (stattrek.com) para configurar una prueba de distribución hipergeométrica con los siguientes parámetros: - tamaño de pop = 15,000 - # de éxitos en población = 3,000 - tamaño de muestra = 400, - # de éxitos en muestra = 100. Obtengo lo siguiente para Probabilidad hipergeométrica P (x = 100) = 0.00224050636447747

El número real de genes superpuestos entre B y C = 100. ¿Es esto mejor que por casualidad? No parece que sea si la posibilidad de que un gen sea enriquecido es 1: 5 (3,000 de 15,000). Es por eso que no entiendo cómo es que mi P (x = 100) que calculé anteriormente es 0.0022. Eso equivale a una probabilidad del 0.2% de que la superposición ocurra por casualidad. ¿No debería ser esto mucho más alto?

Si tomé muestras de 400 genes aleatorios de la gran lista de 15,000, entonces se esperaría que 80 de estos genes se enriquecieran por casualidad (1: 5). El número de genes que en realidad se superponen es de 100, por lo que es un poco mejor que por casualidad.

También intenté encontrar una solución usando las funciones dhyper o phyper en R (usando lo que vi en otra publicación): A = todos los genes en el genoma (15,000) B = genes enriquecidos en RNA-Seq (3,000) C = ChIP genes enriquecidos en chip (400) Aquí está la entrada / salida R (adaptada de una publicación anterior de intercambio de pila):

> totalpop <- 15000    
> sample1 <- 3000    
> sample2 <- 400    
> dhyper(0:2, sample1, totalpop-sample1, sample2)    
[1] 4.431784e-40 4.584209e-38 2.364018e-36    
> phyper(-1:2, sample1, totalpop-sample1, sample2)    
[1] 0.000000e+00 4.431784e-40 4.628526e-38 2.410304e-36    

No estoy seguro de cómo interpretar estos números. Creo que 2.36e-36 es la probabilidad de obtener una superposición completa entre B y C solo por casualidad. Pero esto no tiene sentido, ya que esa probabilidad está mucho más cerca de 1: 5. Si empiezo con 15,000 genes, 3,000 se enriquecerán. Del mismo modo, si empiezo con 400 genes de chips ChIP, 80 de ellos deberían enriquecerse solo en la secuencia de RNA debido a las posibilidades 1: 5 de enriquecimiento en ese conjunto de datos.

¿Cuál es la forma correcta de calcular el valor p, de acuerdo con la distribución hipergeométrica, para la superposición de B y C?

Estás cerca, con tu uso de dhypery phyper, pero no entiendo de dónde 0:2y de dónde -1:2vienes.

El valor p que desea es la probabilidad de obtener 100 o más bolas blancas en una muestra de tamaño 400 de una urna con 3000 bolas blancas y 12000 bolas negras. Aquí hay cuatro formas de calcularlo.

sum(dhyper(100:400, 3000, 12000, 400))
1 - sum(dhyper(0:99, 3000, 12000, 400))
phyper(99, 3000, 12000, 400, lower.tail=FALSE)
1-phyper(99, 3000, 12000, 400)

Estos dan 0.0078.

dhyper(x, m, n, k)da la probabilidad de dibujar exactamente x. En la primera línea, resumimos las probabilidades para 100 - 400; en la segunda línea, tomamos 1 menos la suma de las probabilidades de 0 a 99.

phyper(x, m, n, k)da la probabilidad de obtener xo menos, por lo que phyper(x, m, n, k)es lo mismo que sum(dhyper(0:x, m, n, k)).

El lower.tail=FALSEes un poco confuso. phyper(x, m, n, k, lower.tail=FALSE)es lo mismo que 1-phyper(x, m, n, k), y también lo es la probabilidad de x+1o más. [Nunca recuerdo esto y siempre tengo que verificarlo dos veces].

$\ge$

max(dhyper(0:400, 3000, 12000, 400))$\sim$

Aquí hay una imagen de la distribución hipergeométrica en este caso. Puedes ver que está centrado en 80 (20% de 400) y que 100 está bastante lejos en la cola derecha.

Míralo de esta manera ... Si lo asumiste como un binomio, lo que puede no ser correcto, pero debería ser bastante aproximado ... tu sigma ^ 2 es .8 * .2 * 400 = 64, entonces sigma = 8. Entonces de 80 a 100 has pasado 2.5 desviaciones estándar ... Esto es bastante significativo ... Debería tener un pequeño valor p.