---
title: "Applying bi-clustering algorithms for specifying partial rankings based on U21 data"
author: 'Zsolt T. Kosztyan, Zsuzsanna Banasz, Vivien Valeria Csanyi, Andras Telcs'
output:
  html_document:
    fig_caption: yes
    fig_height: 8
    fig_width: 11
    highlight: haddock
    theme: lumen
    toc: yes
    toc_float:
      collapsed: no
  html_notebook:
    fig_caption: yes
    highlight: haddock
    theme: lumen
    toc: yes
    toc_float:
      collapsed: no
  pdf_document:
    keep_tex: yes
    latex_engine: lualatex
    toc: yes
  word_document:
    toc: yes
---

Workflow

1. Loading libraries and data
  + Loading necessary libraries
  + Loading necessary data
2. Estimation of the number of bi-clusters
  + Seriation of data
  + Seriation of normalized data
  + Estimation of the number of bi-clusters
3. Applying iBBiG on normalized data
  + Calculation of top league bi-cluster
  + Calculation of F-tests of the top league bi-cluster
  + Heat maps of the top league
  + Checking stability of the top league bi-cluster
  + Calculation of partial rankings
  + Comparison of partial rankings to the adequate part of the original ranking
4. Applying iBBiG on reversed normalized data
  + Calculation of lower league bi-cluster
  + Calculation of F-tests of the lower league bi-cluster
  + Heat maps of the lower league
  + Checking stability of the lower league bi-cluster
  + Calculation of partial rankings
  + Comparison of partial rankings to the adequate part of the original ranking
5. Applying BicARE
  + Calculation of middle league bi-cluster
  + Calculation of F-tests of the middle league bi-cluster
  + Heat maps of the middle league
  + Checking stability of the middle league bi-cluster
  + Calculation of partial rankings
  + Comparison of partial rankings to the adequate part of the original ranking
6. Specifying overlaps
  + Venn diagram of indicators
  + Venn diagram of countries
7. Summing up the conclusions

#Loading libraries and data

##Loading necessary libraries
```{r,warning=FALSE,message=FALSE}
#Seriation and biclustering
library(seriation) #For seriate a matrix
library(biclust) #General biclustering functions
library(iBBiG) #iBBiG biclustering
library(BicARE) #BicARE biclustering
library(BcDiag) #Diagnostics of biclusters

#Additional statistic packages
library(factoextra) #Additional distance functions
library(corrgram) #Correlation diagram
library(VennDiagram) #VennDiagram
library(venneuler) #VennDiagram

#For formatted and interactive tables and figures
library(knitr) #Formatted tables, figures
```

Citations

```{r,warning=FALSE,message=FALSE}
#Seriation and biclustering
print(citation("seriation"),style="text") #For seriate a matrix
print(citation("biclust"),style="text") #General biclustering functions
print(citation("iBBiG"),style="text") #iBBiG biclustering
print(citation("BicARE"),style="text") #BicARE biclustering
print(citation("BcDiag"),style="text") #Diagnostics of biclusters

#Additional statistic packages
print(citation("factoextra"),style="text") #Additional distance functions
print(citation("VennDiagram"),style="text") #VennDiagram
print(citation("venneuler"),style="text") #VennDiagram

#For formatted and interactive tables and figures
print(citation("knitr"),style="text") #Formatted tables, figures
```

##Loading necessary data
```{r,warning=FALSE,message=FALSE}
load("U21_BIC.RData")
#ALL: 50 by 24 dataframe of normalized data, with two letters indicator names
#U21: 50 by 24 dataframe of normalized data, with original indicator names
#mtx: 50 by 24 matrix of normalized data, with two letters indicator names
#orig_mtx: 50 by 24 matrix of unnormalized data, with two letters indicator names
#U21_rank: 50 by 1 matrix (vector) of U21'S country ranks
#weights: 1 by 24 matrix (vector) of weights 
```

#Estimation of the number of bi-clusters
##Seriation of data

Firstly, to obtain a preliminary picture of the possible bi-clusters, the seriation of the original data is done.
```{r,warning=FALSE,message=FALSE}
ORIG_SER <- c(seriate(get_dist(orig_mtx,"euclidean"),method="HC_COMPLETE"),seriate(get_dist(t(orig_mtx),"spearman"),method="HC_COMPLETE"))
ORIG_ORDERED<-orig_mtx[get_order(ORIG_SER,dim=1),get_order(ORIG_SER,dim=2)]
rownames(ORIG_ORDERED) <- rownames(mtx[get_order(ORIG_SER,dim=1),get_order(ORIG_SER,dim=2)])
```

After seriation, plotting the heat map of the original data (see **Fig. 1**) can help to detect the possible bi-clusters. On **Fig. 1** blue cells indicate a lower value of the indicators, and red cells mean a higher value of the indicators.

```{r fig_orig_seriated,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 7,fig.cap="**Fig. 1** Heat map of the original matrix after seriation"}

hmap(as.matrix(ORIG_ORDERED),col=bluered(100),showdist="both")
```

##Seriation of normalized data

After the normalization of the original data, seriation is done on the normalized data. 
```{r,warning=FALSE,message=FALSE}
ALL_SER <- c(seriate(get_dist(ALL,"euclidean"),method="HC_COMPLETE"),seriate(get_dist(t(ALL),"spearman"),method="HC_COMPLETE"))
ALL_ORDERED<-ALL[get_order(ALL_SER,dim=1),get_order(ALL_SER,dim=2)]
rownames(ALL_ORDERED) <- rownames(mtx[get_order(ALL_SER,dim=1),get_order(ALL_SER,dim=2)])
```

##Estimating the number of bi-clusters
After seriation, plotting the heat map of the normalized original data (see **Fig. 2**) can help to detect the possible bi-clusters. 
```{r fig_seriated,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 7,fig.cap="**Fig. 2** Heat map of the normalized and seriated matrix"}

hmap(as.matrix(ALL_ORDERED),col=bluered(100),showdist="both")
```

After seriation, two bigger homogeneous blocks can be identified based on **Fig. 2**. The block of the red cells on the top left corner on **Fig. 2** indicates the top league, while the bigger blue block, which indicates the remaining (lower) league, can be discovered in the bottom of the figure. The dendrogram of two-way clustering also shows that regarding rows and columns two main blocks can be specified.

#Applying iBBiG on normalized data

Applying the iBBiG method on the original normalized data can specify the top league bi-clusters. First of all, the number of bi-clusters needed to be determined.
```{r,warning=FALSE,message=FALSE}
res <- iBBiG(binaryMatrix=binarize(mtx,threshold = 0.50),nModules = 2,alpha=0.05,pop_size = 100,mutation = 0.08,stagnation = 50,selection_pressure = 1.2,max_sp = 15,success_ratio = 0.8)
summary(res)
```
The results above show that the iBBiG algorithm found 2 bi-clusters. In the first bi-cluster there are 23 rows (countries) and 19 columns (indicators), in the second bi-cluster there are 22 rows (countries) and 5 columns (indicators). The next step is to calculate F-tests to determine which bi-cluster is significant.

```{r,message=FALSE,warning=FALSE}
Obs.FStat <- NULL
for (i in c(1:res@Number)){
  Obs.FStat[[i]] <- computeObservedFstat(x=mtx,bicResult=res,number=i)
}
kable(Obs.FStat,caption = "**Table 1** The results of row and column effects for two bi-clusters")
```
In **Table 1** the results of the F-tests can be seen. In the case of the first bi-cluster, both the row and column effect are significant. On the other hand, in the case of the second bi-cluster the row effect is not significant; thus the second bi-cluster is not significant. Only the first bi-cluster is significant which is the top league (League A).

##Calculating top league bi-cluster

The next step is to identify the top league bi-cluster.
```{r,warning=FALSE,message=FALSE}
res <- iBBiG(binaryMatrix=binarize(mtx,threshold = 0.50),nModules = 1,alpha=0.05,pop_size = 100,mutation = 0.08,stagnation = 50,selection_pressure = 1.2,max_sp = 15,success_ratio = 0.8)
summary(res)
```
The results of the iBBiG algorithm can be seen above. The algorithm found one bi-cluster in which there are 23 rows (countries) and 19 columns (indicators). The following step is to calculate the F-tests to determine whether the bi-cluster is significant or not.

##Calculating F-tests of the top league bi-cluster

In **Table 2** the results of the row and column effects (F-tests) can be read. Both the row and column effects are significant, so the top league bi-cluster is significant.
```{r,message=FALSE,warning=FALSE}
Obs.FStat <- NULL
for (i in c(1:res@Number)){
  Obs.FStat[[i]] <- computeObservedFstat(x=mtx,bicResult=res,number=i)
}
kable(Obs.FStat,caption = "**Table 2** The results of row and column effects for the top league")
```

To see the differences between those countries that are included in the top league and those that are not in the top league, one can plot the means, medians, variances, and median absolute deviations (MAD) for rows (countries) (see **Fig. 3**). (The horizontal axis of **Fig. 3** shows the countries, the vertical axis shows the mentioned parameters on normalized (0-1) data of indicators.) On **Fig. 3** the red line indicates those countries that are included in the top league, while the black line shows the remaining countries. The countries in the top league have higher means and medians, and have lower variances.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 3** Mean, median, variance and MAD for the top league countries"}
exploreBic(dset=mtx,bres=res,mname='biclust',pfor='all',gby='genes',bnum=1)
```

One can do the same for the columns (indicators). **Fig. 4** shows the differences between the means, medians, variances, and median absolute deviations (MAD) for the indicators. (The horizontal axis of **Fig. 4** shows the indicators, the vertical axis shows the mentioned parameters on normalized (0-1) data of countries.) The red line indicates those indicators that are in the top league, while the black line shows the remaining indicators.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 4** Mean, median, variance and MAD for the top league indicators"}
exploreBic(dset=mtx,bres=res,mname='biclust',pfor='all',gby='conditions',bnum=1)
```


##Heat maps of the top league

**Fig. 5** depicts the heat map of the top league. The top league contains 23 countries and 19 indicators. The red cells mean higher values in the case of indicators, while blue cells symbol lower values.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 5** Heat map of the top league"}
drawHeatmap(x=mtx,bicResult=res,number=1,beamercolor = TRUE,paleta=bluered(100))
```

**Fig. 6** shows the heat map of the results of iBBiG on normalized data. The top league can be seen in the upper left corner of **Fig. 6** which is enlarged on **Fig.5**.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9, fig.cap="**Fig. 6** Heat map of the output after applying iBBiG on normalized data"}
drawHeatmap(x=mtx,bicResult=res,local=FALSE,number=1,beamercolor = TRUE,paleta=bluered(100))
```


##Checking stability of the top league bi-cluster

The next step is to check the stability of the top league bi-cluster. The results of the bootstrap can be seen below.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
Bootstrap <- diagnoseColRow(x=as.matrix(ALL),bicResult=res,number=1,
  nResamplings=100,replace=TRUE)
Bootstrap
```

Since the *bootstrapPvalueRow* and the *bootstrapPvalueCol* are lower than 0.01, the bi-cluster is stable.

##Calculating partial rankings

The next step is to calculate partial rankings based on the indicators and countries in the top league. The calculation of each score of each indicator of the countries are done according to the method used in the original U21 ranking. The original weights of U21's indicators are used, and the total scores for the countries are calculated using the selected indicators in the given bi-cluster. **Table 3** shows the result of partial ranking in the case of the top league.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
C<-biclust::bicluster(orig_mtx,res,number=1)
B<-as.data.frame(C[[1]])
selectedweight<-weights[,colnames(B)]
colsR <- grep("R", names(B), value=T)
colsC <- grep("C", names(B), value=T)
colsO <- grep("O", names(B), value=T)
colsE <- grep("E", names(B), value=T)
swR<-weights[,colsR]
swC<-weights[,colsC]
swO<-weights[,colsO]
swE<-weights[,colsE]
BR<-rowSums(B[,colsR]*swR)*100/max(rowSums(B[,colsR]*swR))
BC<-rowSums(B[,colsC]*swC)*100/max(rowSums(B[,colsC]*swC))
BO<-rowSums(B[,colsO]*swO)*100/max(rowSums(B[,colsO]*swO))
BE<-rowSums(B[,colsE]*swE)*100/max(rowSums(B[,colsE]*swE))
B$R_Score<-BR
B$E_Score<-BE
B$C_Score<-BC
B$O_Score<-BO
B$Overall_Score<-rowSums(cbind(BR,BE,BC,BO)*c(0.2,0.2,0.2,0.4))*100/max(rowSums(cbind(BR,BE,BC,BO)*c(0.2,0.2,0.2,0.4)))
B$Rank<-rank(-B$Overall_Score)
kable(B,caption = "**Table 3** League A by iBBiG",digits = 2)
```


##Comparing partial rankings to the adequate part of the original ranking

The final step is to compare the results of the partial rankings to the corresponding part of the original rankings. The results of the Spearman's rank correlation shows a strong positive relationship between the two rankings.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
cor(B$Rank,U21_rank[rownames(B),],method="spearman")
cor.test(B$Rank,U21_rank[rownames(B),],method="spearman")
```

#Applying iBBiG on reversed normalized data

To identify the lower league(s), the reverse data matrix is needed to be calculated firstly. Applying the iBBiG method on the reversed normalized data can specify the lower league bi-clusters. 

```{r,warning=FALSE,message=FALSE} 
rmtx <- 1-mtx 
rres <- iBBiG(binaryMatrix=binarize(rmtx,threshold = NA),nModules = 2,alpha=0.08,pop_size = 100,mutation = 0.08,stagnation = 50,selection_pressure = 1.2,max_sp = 15,success_ratio = 0.8)
summary(rres) 
``` 

The results above show that the iBBiG algorithm found 2 bi-clusters on the reversed data. In the first bi-cluster there are 38 rows (countries) and 19 columns (indicators), in the second bi-cluster there are 11 rows (countries) and 7 columns (indicators). The next step is to calculate F-tests to determine which bi-cluster is significant.

```{r,message=FALSE,warning=FALSE} 
Obs.FStat <- NULL 
for (i in c(1:rres@Number)){
  Obs.FStat[[i]] <- computeObservedFstat(x=rmtx,bicResult=rres,number=i) 
  } 
kable(Obs.FStat,caption = "**Table 4** The results of row and column effects of the reversed data for the two bi-clusters")
```

In **Table 4** the results of the F-tests can be seen. In the case of the first bi-cluster, both the row and column effect are significant. On the other hand, in the case of the second bi-cluster the row effect is not significant; thus the second bi-cluster is not significant. Only the first bi-cluster is significant which is the lower league (League C).

##Calculating lower league bi-cluster

The next step is to identify the lower league bi-cluster.

```{r,warning=FALSE,message=FALSE,echo=FALSE} 
  rres <- iBBiG(binaryMatrix=binarize(rmtx,threshold = NA),nModules = 1,alpha=0.3,pop_size = 100,mutation = 0.08,stagnation = 50,selection_pressure = 1.2,max_sp = 15,success_ratio = 0.8)
  score<-rres@Clusterscores
summary(rres) 
```

The results of the iBBiG algorithm on the reversed data can be seen above. The algorithm found one bi-cluster in which there are 38 rows (countries) and 19 columns (indicators). The following step is to calculate the F-tests to determine whether the bi-cluster is significant or not.

##Calculating F-tests of the lower league bi-cluster

In **Table 5** the results of the row and column effects (F-tests) can be read. Both the row and column effects are significant, so the lower league bi-cluster is significant.

```{r,message=FALSE,warning=FALSE} 
Obs.FStat <- NULL 
for (i in c(1:rres@Number)){
  Obs.FStat[[i]] <- computeObservedFstat(x=rmtx,bicResult=rres,number=i) 
  } 
kable(Obs.FStat,caption = "**Table 5** The results of row and column effects of the lower league")
```

To see the differences between those countries that are included in the lower league and those that are not in the lower league, one can plot the means, medians, variances, and median absolute deviations (MAD) for rows (countries) (see **Fig. 7**). On **Fig. 7** the red line indicates those countries that are included in the lower league, while the black line shows the remaining countries. The countries in the lower league have lower means, medians, and variances.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 7** Mean, median, variance and MAD for lower leauge countries"}
exploreBic(dset=mtx,bres=rres,mname='biclust',pfor='all',gby='genes',bnum=1)
```

One can do the same for the columns (indicators). **Fig. 8** shows the differences between the means, medians, variances, and median absolute deviations (MAD) for the indicators. The red line indicates those indicators that are in the lower league, while the black line shows the remaining indicators.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 8** Mean, median, variance and MAD for lower leauge indicators"}
exploreBic(dset=mtx,bres=rres,mname='biclust',pfor='all',gby='conditions',bnum=1)
```


##Heat maps of the lower league

**Fig. 9** depicts the heat map of the lower league. The lower league contains 38 countries and 19 indicators. The red cells mean higher values in the case of indicators, while blue cells symbol lower values.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 9** Heat map of the lower leauge"}
drawHeatmap(x=mtx,bicResult=rres,number=1,paleta=bluered(100),beamercolor = TRUE)
``` 

**Fig. 10** shows results of iBBiG on reversed normalized data. The lower league can be seen in the upper left corner of **Fig. 10** which is enlarged on **Fig.9**.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9, fig.cap="**Fig. 10** Heat map of the output after applying iBBiG on reversed normalized data"}
drawHeatmap(x=mtx,bicResult=rres,local=FALSE,number=1,paleta=bluered(100),beamercolor = TRUE) 
```

##Checking stability of the lower league bi-cluster

The next step is to check the stability of the lower league bi-cluster. The results of the bootstrap can be seen below.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
Bootstrap <- diagnoseColRow(x=as.matrix(ALL),bicResult=rres,number=1,
  nResamplings=100,replace=TRUE)
Bootstrap
```

Since the *bootstrapPvalueRow* and the *bootstrapPvalueCol* are lower than 0.01, the bi-cluster is stable.


##Calculating partial rankings

The next step is to calculate partial rankings based on the indicators and countries in the lower league. The calculation of each score of each indicator of the countries are done according to the method used in the original U21 ranking. The original weights of U21�s indicators are used, and the total scores for the countries are calculated using the selected indicators in the given bi-cluster. **Table 6** shows the result of partial ranking in the case of the lower league.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
C<-biclust::bicluster(orig_mtx,rres,number=1)
B<-as.data.frame(C[[1]])
selectedweight<-weights[,colnames(B)]
colsR <- grep("R", names(B), value=T)
colsC <- grep("C", names(B), value=T)
colsO <- grep("O", names(B), value=T)
swR<-weights[,colsR]
swC<-weights[,colsC]
swO<-weights[,colsO]
BR<-rowSums(B[,colsR]*swR)*100/max(rowSums(B[,colsR]*swR))
BC<-rowSums(B[,colsC]*swC)*100/max(rowSums(B[,colsC]*swC))
BO<-rowSums(B[,colsO]*swO)*100/max(rowSums(B[,colsO]*swO))
B$R_Score<-BR
B$C_Score<-BC
B$O_Score<-BO
B$Overall_Score<-rowSums(cbind(BR,BC,BO)*c(0.2,0.2,0.4))*100/max(rowSums(cbind(BR,BC,BO)*c(0.2,0.2,0.4)))
B$Rank<-rank(-B$Overall_Score)
kable(B,caption = "**Table 6** League C by iBBiG",digits = 2)
```

##Comparing partial rankings to the adequate part of the original ranking

The final step is to compare the results of the partial rankings to the corresponding part of the original rankings. The results of the Spearman�s rank correlation shows a strong positive relationship between the two ranking.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
cor(B$Rank,U21_rank[rownames(B),],method="spearman")
cor.test(B$Rank,U21_rank[rownames(B),],method="spearman")
```

#Applying BicARE on normalized data

With the BicARE method the middle league(s) can be specified. In order to use the BicARE method, the first step is to convert the BicARE object to a Biclust object. 
```{r,message=FALSE,warning=FALSE}
bicare2biclust <- function(x){
	if(class(x)=="Biclust"){ 
		return(x)
	} else if(class(x)=="biclustering"){
		Parameters <- list(numberofbicluster=x$param[1,2],residuthreshold=x$param[2,2],genesinitialprobability=x$param[3,2],samplesinitialprobability=x$param[4,2],numberofiterations=x$param[5,2],date=x$param[6,2])
		RowxNumber <- t(x$bicRow==1)
		NumberxCol <- x$bicCol==1
		Number <- as.numeric(dim(RowxNumber)[2])
		info <- list()
		return(new("Biclust",Parameters=Parameters,RowxNumber=RowxNumber,NumberxCol=NumberxCol,Number=Number,info=info))
	}
}
```

In the next step the BicARE algorithm determines the number of the bi-clusters.

```{r,warning=F,message=FALSE}
BICARE_res <- FLOC(Data=as.matrix(ALL),k=2,pGene=1,pSample=1,r=1e-16,N=17,M=6,t=10000,blocGene=NULL,blocSample=NULL)
bres <- bicare2biclust(BICARE_res)
summary(bres)
bres@Parameters$residuthreshold
```

The results above show that the BiCARE algorithm found 2 bi-clusters. In both bi-clusters there are 17 rows (countries) and 6 columns (indicators) which means that the first and the second bi-cluster are same. Because of this, in the following, only the first bi-cluster will be considered. This bi-cluster contains the middle league countries and indicators. 

##Calculating middle league bi-cluster

The next step is to identify the middle league bi-cluster.

```{r,warning=F,message=FALSE}
BICARE_res <- FLOC(Data=as.matrix(ALL),k=1,pGene=1,pSample=1,r=1e-16,N=17,M=6,t=10000,blocGene=NULL,blocSample=NULL)
bres <- bicare2biclust(BICARE_res)
summary(bres)
bres@Parameters$residuthreshold
```

The results of the BicARE algorithm can be seen above. The algorithm found one bi-cluster in which there are 17 rows (countries) and 6 columns (indicators). The following step is to calculate the F-tests to determine whether the bi-cluster is significant or not.

##Calculating F-tests of the middle league bi-cluster

In **Table 7** the results of the row and column effects (F-tests) can be read. Both the row and column effects are significant, so the middle league bi-cluster is significant.

```{r,message=FALSE,message=FALSE}
RObs.FStat <- NULL
for (i in c(1:bres@Number)){
  RObs.FStat[[i]] <- computeObservedFstat(x=mtx,bicResult=bres,number=i)
}
kable(RObs.FStat,caption = "**Table 7** The results of row and column effects of the middle league")
```

To see the differences between those countries that are included in the middle league and those that are not in the middle league, one can plot the means, medians, variances, and median absolute deviations (MAD) for rows (countries) (see **Fig. 11**). On **Fig. 11** the red line indicates those countries that are included in the middle league, while the black line shows the remaining countries. **Fig. 11.** shows that in the case of the countries in the middle league the variances differ the most.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 11** Mean, median, variance and MAD for middle league countries"}
exploreBic(dset=mtx,bres=bres,mname='biclust',pfor='all',gby='genes',bnum=1)
```

One can do the same for the columns (indicators). **Fig. 12** shows the differences between the means, medians, variances, and median absolute deviations (MAD) for the indicators. The red line indicates those indicators that are in the middle league, while the black line shows the remaining indicators. In the case of the indicators in the middle league, all parameters are different.

```{r,message=FALSE,warning=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 12** Mean, median, variance and MAD for middle league indicators"}
exploreBic(dset=mtx,bres=bres,mname='biclust',pfor='all',gby='conditions',bnum=1)
```


##Heat maps of the middle league

**Fig. 13** depicts the heat map of the middle league. The middle league contains 17 countries and 6 indicators. The red cells mean higher values in the case of indicators, while blue cells symbol lower values.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 7, fig.cap="**Fig. 13** Heat map of the middle league"}
drawHeatmap(x=mtx,bicResult=bres,number=1,paleta=bluered(100),beamercolor = TRUE) 
``` 

**Fig. 14** shows the results of BicARE on normalized data. The middle league can be seen in the upper left corner of **Fig. 14** which is enlarged on **Fig.13**.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9, fig.cap="**Fig. 14** Heat map of the output after applying BicARE on normalized data"}
drawHeatmap(x=mtx,bicResult=bres,local=FALSE,number=1,paleta=bluered(100),beamercolor = TRUE) 
```

##Checking stability of the middle league bi-cluster

The next step is to check the stability of the middle league bi-cluster. The results of the bootstrap can be seen below.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
Bootstrap <- diagnoseColRow(x=as.matrix(ALL),bicResult=bres,number=1,
  nResamplings=100,replace=TRUE)
Bootstrap
```

Since the *bootstrapPvalueRow* and the *bootstrapPvalueCol* are lower than 0.01, the bi-cluster is stable.

##Calculating partial rankings

The next step is to calculate partial rankings based on the indicators and countries in the middle league. The calculation of each score of each indicator of the countries are done according to the method used in the original U21 ranking. The original weights of U21�s indicators are used, and the total scores for the countries are calculated using the selected indicators in the given bi-cluster. **Table 8** shows the result of partial ranking in the case of the middle league.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
C<-biclust::bicluster(orig_mtx,bres,number=1)
B<-as.data.frame(C[[1]])
selectedweight<-weights[,colnames(B)]
colsR <- grep("R", names(B), value=T)
colsO <- grep("O", names(B), value=T)
swR<-weights[,colsR]
swO<-weights[,colsO]
BR<-rowSums(B[,colsR]*swR)*100/max(rowSums(B[,colsR]*swR))
BO<-rowSums(B[,colsO]*swO)*100/max(rowSums(B[,colsO]*swO))
B$R_Score<-BR
B$O_Score<-BO
B$Overall_Score<-rowSums(cbind(BR,BO)*c(0.2,0.4))*100/max(rowSums(cbind(BR,BO)*c(0.2,0.4)))
B$Rank<-rank(-B$Overall_Score)
kable(B,caption = "**Table 8** League B by BicARE",digits = 2)
```

##Comparing partial rankings to the adequate part of the original ranking

The final step is to compare the results of the partial rankings to the corresponding part of the original rankings. The results of the Spearman�s rank correlation shows a strong positive relationship between the two ranking.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9}
cor(B$Rank,U21_rank[rownames(B),],method="spearman")
cor.test(B$Rank,U21_rank[rownames(B),],method="spearman")
```


#Specifying overlaps

Bi-cluster might (or might not) overlap. Analyzing the overlaps can give a subtler picture of countries and indicators.

##Venn diagram of indicators

**Fig. 15** shows the Venn diagram of the indicators. The red circle contains those indicators that are in the upper league (League A), in the the blue circle those indicators can be found that are in the lower league (League C). The green circle indicates the middle league's (League B) indicators.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9, fig.cap="**Fig. 15** Venn diagram of the indicators"}
T<-biclust::bicluster(orig_mtx,res,number = 1)
League_A<-T[[1]]
R<-biclust::bicluster(orig_mtx,rres,number = 1)
League_C<-R[[1]]
B<-biclust::bicluster(orig_mtx,bres,number = 1)
League_B<-B[[1]]
draw.triple.venn(area1=nrow(as.matrix(colnames(League_A))),area2=nrow(as.matrix(colnames(League_B))),area3=nrow(as.matrix(colnames(League_C))),n12=nrow(as.matrix(intersect(colnames(League_A),colnames(League_B)))),n13=nrow(as.matrix(intersect(colnames(League_A),colnames(League_C)))),n23=nrow(as.matrix(intersect(colnames(League_B),colnames(League_C)))),n123=nrow(as.matrix(intersect(intersect(colnames(League_A),colnames(League_B)),colnames(League_C)))),fill=c("red","green","blue"),category=c("League A","League B","League C"),main="Venn Diagram of indicators")
```

One interesting finding is, that all of the indicators in the middle league are part of the lower league's indicators.


##Venn diagram of countries

**Fig. 16** shows the Venn diagram of the countries. The red circle contains those countries that are in the upper league (League A), in the the blue circle those countries can be found that are in the lower league (League C). The green circle indicates the middle league's (League B) countries.

```{r,warning=FALSE,message=FALSE,fig.width = 7, fig.height = 9, fig.cap="**Fig. 16** Venn diagram of the countries"}
draw.triple.venn(area1=nrow(as.matrix(rownames(League_A))),area2=nrow(as.matrix(rownames(League_B))),area3=nrow(as.matrix(rownames(League_C))),n12=nrow(as.matrix(intersect(rownames(League_A),rownames(League_B)))),n13=nrow(as.matrix(intersect(rownames(League_A),rownames(League_C)))),n23=nrow(as.matrix(intersect(rownames(League_B),rownames(League_C)))),n123=nrow(as.matrix(intersect(intersect(rownames(League_A),rownames(League_B)),rownames(League_C)))),fill=c("red","green","blue"),category=c("League A","League B","League C"),main="Venn Diagram of countries")
```

Seven countries are in all of the three leagues. There are 11 countries which belong only to the upper league, and there are 18 countries which are just in the lower league.

#Summing up the conclusions

The above proposed bi-clustering method can help to identify which countries can be compared based on which indicators in university rankings. With the iBBiG method the upper and lower league(s) (*League(s) A and League(s) C*) can be determined and the BicARE method can specify the middle league(s) (*League(s) B*). 
