GENIE3.Rmd

---
title: "GENIE3 vignette"
author: "Van Anh Huynh-Thu (vahuynh@ulg.ac.be)"
date: "`r Sys.Date()`"
package: "`r pkg_ver('GENIE3')`"
vignette: >
  %\VignetteIndexEntry{GENIE3}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
output: 
  BiocStyle::html_document:
    toc: yes
    number_sections: false
    css: corrected.css
  pdf_document:
    toc: yes
---
This is the documentation for the R implementation of GENIE3.

The GENIE3 method is described in:

> Huynh-Thu V. A., Irrthum A., Wehenkel L., and Geurts P. (2010) Inferring regulatory networks from expression data using tree-based methods. *PLoS ONE*, 5(9):e12776.

## Format of expression data

### Format of steady-state expression data

The `GENIE3()` function takes as input argument a gene expression matrix `exprMatr`. Each row of that matrix must correspond to a gene and each column must correspond to a sample. The gene names must be specified in `rownames(exprMatr)`. The sample names can be specified in `colnames(exprMatr)`, but this is not mandatory. For example, the following command lines generate a fake expression matrix (for the purpose of this tutorial only):

```{r}
exprMatr <- matrix(sample(1:10, 100, replace=TRUE), nrow=20)
rownames(exprMatr) <- paste("Gene", 1:20, sep="")
colnames(exprMatr) <- paste("Sample", 1:5, sep="")
head(exprMatr)
```

This matrix contains the expression data of 20 genes from 5 samples. The expression data does not need to be normalised in any particular way (but whether it is normalized/filtered/log-transformed WILL affect the results!). 

## How to run GENIE3

### Run GENIE3 with the default parameters

The following command runs GENIE3 on the expression data `exprMatr` with the default parameters:

```{r}
library(GENIE3)
set.seed(123) # For reproducibility of results
weightMat <- GENIE3(exprMatr)
```
```{r}
dim(weightMat)
weightMat[1:5,1:5]
```

The algorithm outputs a matrix containing the weights of the putative regulatory links, with higher weights corresponding to more likely regulatory links. `weightMat[i,j]` is the weight of the link directed from the $i$-th gene to $j$-th gene.

### Restrict the candidate regulators to a subset of genes

By default, all the genes in `exprMatr` are used as candidate regulators. The list of candidate regulators can however be restricted to a subset of genes. This can be useful when you know which genes are transcription factors.

```{r}
# Genes that are used as candidate regulators
regulators <- c(2, 4, 7)
# Or alternatively:
regulators <- c("Gene2", "Gene4", "Gene7")
weightMat <- GENIE3(exprMatr, regulators=regulators)
```

Here, only `Gene2`, `Gene4` and `Gene7` (respectively corresponding to rows 2, 4 and 7 in `exprMatr`) were used as candidate regulators. In the resulting `weightMat`, the links that are directed from genes that are not candidate regulators have a weight equal to 0.

To request different regulators for each gene & return as list:
```{r}
regulatorsList <- list("Gene1"=rownames(exprMatr)[1:10],
                       "Gene2"=rownames(exprMatr)[10:20],
                       "Gene20"=rownames(exprMatr)[15:20])
set.seed(123)
weightList <- GENIE3(exprMatr, nCores=1, targets=names(regulatorsList), regulators=regulatorsList, returnMatrix=FALSE)
```


### Change the tree-based method and its settings

GENIE3 is based on regression trees. These trees can be learned using either the Random Forest method ^[Breiman L. (2001) Random forests. *Machine learning*, 45(1):5-32.] or the Extra-Trees method ^[Geurts P., Ernst D. and Wehenkel L. (2006) Extremely randomized trees. *Machine learning*, 36(1):3-42.]. The tree-based method can be specified using the `tree.method` parameter (`tree.method="RF"` for Random Forests, which is the default choice, or `tree.method="ET"` for Extra-Trees).

Each tree-based method has two parameters: `K`  and `ntrees`. `K` is the number of candidate regulators that are randomly selected at each tree node for the best split determination. Let $p$ be the number of candidate regulators. `K` must be either:

* `"sqrt"`, which sets $K=\sqrt{p}$. This is the default value.
* `"all"`, which sets $K=p$.
* Or any integer between $1$ and $p$.

The parameter `ntrees` specifies the number of trees that are grown per ensemble. It can be set to any strictly positive integer (the default value is 1000).

An example is shown below:


```{r eval=FALSE}
# Use Extra-Trees (ET) method
# 7 randomly chosen candidate regulators at each node of a tree
# 5 trees per ensemble
weightMat <- GENIE3(exprMatr, treeMethod="ET", K=7, nTrees=50)
```


### Parallel GENIE3

To decrease the computing times, GENIE3 can be run on multiple cores. The parameter `ncores` specifies the number of cores you want to use. For example:

```{r eval=FALSE}
set.seed(123) # For reproducibility of results
weightMat <- GENIE3(exprMatr, nCores=4, verbose=TRUE)
```

Note that `seet.seed` allows to get the same results across different runs, but only within `nCores==1` or `nCores>1`. e.g. A run with `set.seed(123)` and `nCores=1` and another with the same seed but `nCores>1` may provide different results.


### Obtain more information

```{r}
?GENIE3
```


## Get the list of the regulatory links

### Get all the regulatory links

You can obtain the list of all the regulatory links (from most likely to least likely) with this command:

```{r}
linkList <- getLinkList(weightMat)
dim(linkList)
head(linkList)
```

The resulting `linkList` matrix contains the ranking of links. Each row corresponds to a regulatory link. The first column shows the regulator, the second column shows the target gene, and the last column indicates the weight of the link.

(Note that the ranking that is obtained will be slightly different from one run to another. This is due to the intrinsic randomness of the Random Forest and Extra-Trees methods. The variance of the ranking can be decreased by increasing the number of trees per ensemble.)


### Get only the top-ranked links

Usually, one is only interested in extracting the most likely regulatory links. The optional parameter `report.max` sets the number of top-ranked links to report:

```{r eval=FALSE}
linkList <- getLinkList(weightMat, reportMax=5)
```


### Get only the links with a weight higher than some threshold

Alternatively, a threshold can be set on the weights of the links:

```{r eval=FALSE}
linkList <- getLinkList(weightMat, threshold=0.1)
```

### *Important note* on the interpretation of the weights
The weights of the links returned by `GENIE3()` **do not have any statistical meaning** and only provide a way to rank the regulatory links. There is therefore no standard threshold value, and caution must be taken when choosing one.

### Obtain more information

```
?getLinkList
```