---
title: "ClusterJudge"
author: "Adrian Pasculescu"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
%\VignetteIndexEntry{Vignette Title}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
references:
- author:
- family: Gibbons
given: F.D
- family: Roth
given: F.P.
container-title: Genome Research
id: GibRoth2002
issued:
month: 10
volume: 12
year: 2002
publisher: Genome Research
title: 'Judging the Quality of Gene Expression-Based Clustering Methods Using Gene Annotation'
- author:
- family: Cho
given: R.C. et all
container-title: Molecular Cell
id: Cho1998
issued:
month: 7
volume: 2
year: 1998
publisher: Science Direct
title: 'A genome-wide transcriptional analysis of the mitotic cell cycle'
- author:
- family: Tavazoie
given: S. et all
container-title: Nature Genetics
id: Taivasoie1999
issued:
month: 7
volume: 22
year: 1999
publisher: Nature Publishing Group
title: 'Systematic determination of genetic network architecture'
- author:
- family: Cover
given: T.M. et all
container-title: book
id: Cover1991
issued:
month: 1
volume: (ed. D.L. Schilling)
year: 1991
publisher: Wiley-Interscience
title: 'Elements of information theory'
---
```{r echo=FALSE, warning=FALSE, message=FALSE}
devtools::load_all('.')
```
## Introduction
The quality of a cluster analysis can be judged based on external information.
For example when the clustered entities possess attributes, the total mutual information between
the clustering and each of the important and non-redundant attributes of entities might be helpful.
Gene Ontology is an example of a system of attributes that can help in assessing gene clustering resulting
from genomics or proteomics experiments.
A clustering that `agrees well` with a system of attributes should have a pronounced decrease of
the total mutual information when the number of random swaps of entities between clusters
is increased [@GibRoth2002].
The following simple steps must be performed to evaluate a clustering using ClusterJudge:
1. Obtain the clusters as a structure containing entities and their cluster IDs;
2. Obtain the entity-attribute associations and verify the entity-attribute associations for consistency;
3. Often, the number of attributes can be consolidated (reduced) by using a set of _uncertainty_ levels;
4. 'Judge' the original clustering versus clustering of random changes that are applied to the original data.
## Obtaining the clustering
In most of the cases the clustering is obtained from experimental data after
applying some statistical methods such as k-means or cutting a hierarchical clustering at a predefined 'height'.
The following example extracts the Yeast cell cycle clustering of genes related to the Yeast cell
cycle [@Cho1998] [@Taivasoie1999]:
```{r warning=FALSE, message=FALSE}
library('yeastExpData')
```
```{r echo=TRUE}
data(ccyclered)
head(ccyclered)
```
The clustering should not contain any ambiguities ( i.e. one entity should be assigned to just one cluster)
and no NA should be present.
For this package, the clusters must be organized as a named vector:
```{r}
clusters <- ccyclered$Cluster
### convert from Gene names to the new standard of Saccharomyces Genome Database (SGD) gene ids
ccyclered$SGDID <- sub('^S','S00',ccyclered$SGDID)
names(clusters) <- ccyclered$SGDID
str(clusters)
```
Note that, in this example, we modified the names of the genes to match the names of entities in
the Gene Ontology data (see next section).
## Obtaining Entity Attribute associations
For genomics or proteomics data, one source of attributes is the Gene Ontology (http://www.geneontology.org/).
For Yeast Gene Ontology attributes can also be downloaded from the [Saccharomyces Genome Database (SGD)]
(http://www.yeastgenome.org/)
In the ClusterJudge R package there is already a downloaded mapping of Yeast genes to attributes:
```{r, fig.show='hold', echo=TRUE}
data(Yeast.GO.assocs);
str(Yeast.GO.assocs);
head(Yeast.GO.assocs);
validate_association(Yeast.GO.assocs)
```
The *download_Yeast_GO_mapping()* function is also available to download a fresh copy of the mapping.
For other species one can use specialized R and Bioconductor packages such as biomaRt
as in the following `toy` example.
```{r, eval=FALSE, echo=TRUE}
library(biomaRt)
rn <- useDataset("rnorvegicus_gene_ensembl", mart=useMart("ensembl"))
rgd.symbol=c("As3mt", "Borcs7", "Cyp17a1", "Wbp1l", "Sfxn2", "Arl3") ### exemplify for a limited set of genes
entity.attr <- getBM(attributes=c('rgd_symbol','go_id'), filters='rgd_symbol', values=rgd.symbol, mart=rn)
```
The _validate_association()_ function quickly verifies if the input structure has two columns,
if there are no NAs or NULLs and if there are no duplicated associations.
## Consolidating Entity-Attribute associations
In some cases, the attributes may not bring additional useful information. For example, an attribute
may be too specific as it is only assigned to one or two genes (e.g. it may simply be a new name of a gene).
In other cases, two or more attributes may be very correlated and appear on almost the same set of genes.
For the first case a consolidated association will ignore those attributes that belong to very few genes
(*min.entities.per.attr*)
```{r, fig.width=6, fig.height=4}
entities_attribute_stats(Yeast.GO.assocs) ### shows number of entities per attribute distribution
Yeast.GO.assocs.cons1 <- consolidate_entity_attribute(
entity.attribute = Yeast.GO.assocs
, min.entities.per.attr =3 ### keep only attributes associated to 3 or more entities
, mut.inf=FALSE
)
dim(Yeast.GO.assocs)
dim(Yeast.GO.assocs.cons1) ### shows reduction in the number of associations
```
For the second case we use the pre-calculated mutual information [@Cover1991] between attributes and impose
a maximum _uncertainty_ level where uncertainty is defined as mutual information divided by the maximum
of mutual information.
```{r, fig.width=6, fig.height=4}
data(mi.GO.Yeast)
Yeast.GO.assocs.cons <- consolidate_entity_attribute(
entity.attribute = Yeast.GO.assocs
, min.entities.per.attr =3
, mut.inf=mi.GO.Yeast ### use precalculated mutual information
, U.limit = c(0.1, 0.001) ### calculate consolidated association for these uncertainty levels
) ### shows distribution of the number of pairs of attributes by Uncertainty
str(Yeast.GO.assocs.cons)
```
The lower the Uncertainty limit the more attributes are ignored and the faster the judging
of clusters will be performed.
## Calculating the mutual information
In the case where the entity-attributes associations are downloaded from other organisms or other sources,
mutual information can be calculated using this function: *attribute_mut_inf* :
```{r}
data(Yeast.GO.assocs)
### because it takes time, we use a small sampled subset of associations
entity.attribute.sampled <- Yeast.GO.assocs[sample(1:nrow(Yeast.GO.assocs),100),]
mi.GO.Yeast.sampled <- attribute_mut_inf(
entity.attribute = entity.attribute.sampled
, show.progress = FALSE ## for this small example do not print progress
)
str(mi.GO.Yeast.sampled)
```
The mutual information of the pair of attributes A,B is defined as
mi(A,b) = H(A)+H(B) - H(cbind(A,B))
where H is the entropy.
## Judging the clustering for selected `uncertainty` levels
The _clusters_ named vector prepared above and one of the consolidated entity-attribute
associations (for example the one at 0.001 level of Uncertainty) are needed at this step.
The _ClusterJudge_ will swap at random more and more entities between clusters and will calculate
and plot the new total mutual information obtained at each iteration.
If the total mutual information decreases by increasing the number of random swaps we can conclude
that there is a good agreement between the clustering of experimental data and the ensemble of information
provided by the entity-attribute associations.
```{r, fig.width=6, fig.height=4}
mi.by.swaps<-clusterJudge(
clusters = clusters
, entity.attribute=Yeast.GO.assocs.cons[["0.001"]]
, plot.notes='Yeast clusters judged at uncertainty level 0.001 - Ref: Tavazoie S,& all
`Systematic determination of genetic network architecture. Nat Genet. 1999`'
, plot.saveRDS.file= 'cj.rds') ### save the plot for later use
p <- readRDS('cj.rds') ### retrieve the previous plot
pdf('cj.pdf'); plot(p); dev.off() ### plot on another device
```
In most of the cases mutual information should have a decreasing trend when the number of shuffles increases.
However, in rare cases, a consistent increasing trend might suggest novel discoveries; only when both experimental data
and clustering are *very reliable*. In these rare cases a *high confident* clustering might complement or even contradict
information accumulated in time into entity-attribute associations.
*ClusterJudge* can also be seen as a complement to many of the Bioconductor packages
related to evaluation and comparison of Clustering (such as [clustComp](http://www.bioconductor.org/packages/devel/bioc/html/clustComp.html ),
[clusterExperiment](http://www.bioconductor.org/packages/devel/bioc/html/clusterExperiment.html ),
[clusterSignificance](http://www.bioconductor.org/packages/devel/bioc/html/ClusterSignificance.html ),
[GOpro](http://www.bioconductor.org/packages/devel/bioc/html/GOpro.html ),
[multiClust](http://www.bioconductor.org/packages/devel/bioc/html/multiClust.html ) ).
## References