---
title: "Test de Fisher unilateral"  
author: "Guillermo Ayala Gallego"  
date: "`r Sys.Date()`"  
address: Departamento de Estadística e I.O. Universidad de Valencia  
format: 
    revealjs:  
      scrollable: true
      echo: true
      embed-resources: true  
    pdf:
       toc: true
       number-sections: false
       colorlinks: true
---


## Análisis de sobre representación  

- Hasta este momento hemos obtenido una ordenación de los genes.
- Esta lista la hemos estudiado pretendiendo que cuando mayor sea la 
expresión diferencial del gen este aparezca antes en la lista.
- De este modo el primer gen en la lista es el **marginalmente**
tiene una mayor expresión diferencial y así sucesivamente.
- De hecho, el p-valor no es más que una medida de esa diferenciación.
- Una expresión muy utilizada en la literatura es que que hay una asociación
entre el gen (su expresión) y el fenotipo.
- Cuando fijamos una tasa de error (FDR) lo que hacemos es fijar un 
punto de corte.
- Los genes que están antes del punto de
corte se consideran significativos y los que siguen no. 
- ¿Cómo interpretamos este grupo de genes?

## Test de Fisher  

-  Sea $G$ es el conjunto de genes considerado (posiblemente los contemplados
en nuestro estudio).  
- $S_0$ ($\subset G$) el conjunto de genes que nuestro estudio ha indicado 
como significativo.
- $S_1$ un conjunto de genes predefinido (misma función, misma localización).


+-------------+------------------+-----------------+-----------------+
|             |   $S_1$          |   $S_1^c$       |   Total         |
+=============+==================+=================+=================+
|   $S_0$     |   $n_{11}$       |   $n_{12}$      |   $n_{1 \cdot}$ |
+-------------+------------------+-----------------+-----------------+    
|    $S_0^c$  |    $n_{21}$      |    $n_{22}$     |    $n_{2 \cdot}$|  
+-------------+------------------+-----------------+-----------------+    
|   Total     |   $n_{\cdot 1}$  |   $n_{\cdot 2}$ |   $N$           |
+-------------+------------------+-----------------+-----------------+

- Bajo la hipótesis de que no hay ningún tipo de asociación entre 
fila y columna,  entonces la probabilidad de **lo observado** es
$$
P(N_{11} = n_{11}) = \frac{\binom{n_{1 \cdot}}{n_{11}}
\binom{n_{2 \cdot}}{n_{\cdot 1}- n_{11}}}{\binom{n}{n_{\cdot 1}}}.
$$


- Supongamos que estamos contrastando la posible sobrerepresentación.
- **Bajo la hipótesis de independencia** (condicionada a las marginales),
**rechazaríamos** la hipótesis nula de independencia para un valor mayor 
o igual al observado.
- El p-valor sería 
$$
p = P(N_{11} \geq n_{11}) = \sum_{t= n_{11}}^{\min\{n_{1 \cdot},n_{\cdot 1}\}} =
\frac{\binom{n_{1 \cdot}}{t}
\binom{n_{2 \cdot}}{n_{\cdot 1}- t}}{\binom{N}{n_{\cdot 1}}}.
$$

- Un ejemplo 

+-------------+------------------+-----------------+-----------------+
|             |   $S_1$          |   $S_1^c$       |   Total         |
+=============+==================+=================+=================+
|   $S_0$     |       30         |      40         |       70        |
+-------------+------------------+-----------------+-----------------+    
|    $S_0^c$  |     120          |     156         |      276        |  
+-------------+------------------+-----------------+-----------------+    
|   Total     |      150         |     196         |     346         |
+-------------+------------------+-----------------+-----------------+

Podemos utilizar el **test exacto de Fisher direccional** (o unilateral
o de una cola).

```{r }
#| label: fisher1
conteos = matrix(c(30,120,40,156),ncol=2)
fisher.test(conteos,alternative = "greater")
```


## Colecciones de grupos de genes    

1. Cargamos paquete y datos  

```{r }
#| label: fisher2
pacman::p_load(EnrichmentBrowser,Biobase,SummarizedExperiment)
data(gse21942,package="tamidata")
se21942 = makeSummarizedExperimentFromExpressionSet(gse21942)
```

Podemos ver que los genes vienen identificados por los
correspondientes al fabricante AffyID o PROBEID.

```{r  }
#| label: fisher3
head(rowData(se21942))
```
Y los metadatos o variables fenotípicas.

```{r  }
#| label: fisher4
head(colData(se21942))
```

2. Preparamos el **ExpressionSet**.

```{r  }
#| label: fisher5
se21942=probe2gene(se21942)
```
¿Qué hemos hecho?

Modificar los identificadores de los genes sustituyéndolos por los
ENTREZID.

```{r  }
#| label: fisher6
head(rowData(se21942))
```

Introducir una variable fenotípica **GROUP** con valores
0 y 1.

3. Análisis de expresión diferencial utilizando el procedimiento **Limma**.

```{r }
#| label: fisher7
se21942 = deAna(expr = se21942) ## t-test moderados
```

Podemos ver los resultados en **fData**.

```{r }
#| label: fisher8
head(rowData(se21942))
```

4. Grupos de genes con Gene Ontology.

```{r }
#| label: fisher9
#| eval: false
hsaGO =  getGenesets("hsa",onto="BP")
```

```{r}
#| label: fisher10
#| echo: false
dirTamiData = "/home/gag/ownCloud/alltami/tami-data/"
load(paste(dirTamiData,"hsaGO.rda",sep=""))
```

¿Cómo son estos grupos? ¿A qué se refieren?

```{r }
#| label: fisher11
head(names(hsaGO))
```

¿Qué genes, con código ENTREZID, forman el grupo 3?

```{r }
#| label: fisher12
hsaGO[[3]]
names(hsaGO[3])
```
O bien con
```{r }
#| label: fisher13
hsaGO$"GO:0000012_single_strand_break_repair"
```

5. Sobre representación con GO

Realizamos un análisis de sobre representación con el test de Fisher 
unilateral. 

```{r }
#| label: fisher14
se21942.oraGO= sbea(method="ora", se=se21942, gs=hsaGO,
                      perm=0, alpha=0.05)
```

```{r }
#| label: fisher15
gsRanking(se21942.oraGO)
```

## Análisis de sobre-representación de tamidata::gse21942


1. Análisis de expresión diferencial  

- Utiliza el procedimiento del test de la t moderado de **Limma**.

```{r}
#| label: fisher16
pacman::p_load(tami)
gse21942_deo = dema(x=gse21942,
                    y="FactorValue..DISEASE.STATE.",
                    test = rowtmod,correction = "BH",
                    fdr = .0001,foutput = "gse21942")
```

```{r}
#| label: fisher17
browseURL(glimpse(gse21942_deo))
```


2. Utilizamos un test de Fisher unilateral para realizar un análisis
   de sobre-representación. 

```{r}
#| label: fisher18
gse21942.tidy = tidy(gse21942_deo)
sel0 = which(gse21942.tidy[,"adjp"] < .05)
set0 = gse21942.tidy[sel0,"ENTREZID"]
gse21942_gsao = ora(set0,gsc=hsaGO)
head(gse21942_gsao)
```


