####loading previous saved FPKM data

setwd('~/Documents/work/paper_write/genome_training/R_training/')
fpkm <- read.table('vehicle_drug_feature_counts.fpkm.txt',
                   header = T, sep = '\t', row.names = 1)
head(fpkm)
##                      c1        c2        t1        t2
## PF3D7_0621350  0.000000  7.214287  3.381795  4.057177
## PF3D7_1314600 36.160450 28.121355  7.404854 12.813007
## PF3D7_0209400 90.954041 83.522938 66.440631 59.494751
## PF3D7_1142200  4.206178 10.297819  6.194962  1.520215
## PF3D7_0507200 12.145795 12.447656 15.430355 11.706097
## PF3D7_0806000 13.944896 10.452966 10.511226 10.951166



####biological correlation calculation and heatmap

corF <- cor(fpkm)
corF
##           c1        c2        t1        t2
## c1 1.0000000 0.9636821 0.8379175 0.8543578
## c2 0.9636821 1.0000000 0.8391606 0.8452566
## t1 0.8379175 0.8391606 1.0000000 0.9936924
## t2 0.8543578 0.8452566 0.9936924 1.0000000
heatmap(corF)

##using heatmap.2 to draw heatmap
library(gplots)
## 
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
## 
##     lowess
heatmap.2(corF, Colv= F, dendrogram = 'row',
          trace = 'n',      ## hist plot line inside heatmap
          denscol = 'navy'  ## hist plot line in color key
          )

##change the direction of column lable
heatmap.2(corF, Colv= F, dendrogram = 'row',
          trace = 'n', 
          denscol = 'navy', srtCol = 0 ## the angel of y column lable
          )

## change the color of heatmap
colPlate <- colorRampPalette(c('red', 'white', 'navy'))(100)
heatmap.2(corF, Colv= F, dendrogram = 'row',
          trace = 'n',  
          denscol = 'navy', srtCol = 0, ## the angel of y column lable
          col = colPlate
          )

## removing the lable in the heatmap
## change the size of the key
heatmap.2(corF, Colv= F, dendrogram = 'row',
          trace = 'n',  
          denscol = NA, srtCol = 0, ## the angel of y column lable
          col = colPlate, 
          #------------#
          key.title = '',
          key.ylab = NA, key.xlab = NA,
          keysize=0.75, 
          #remove the y-axis#
          #key.par = list(mar = c(1,0.3,1,1))
          )



####Drawing heatmap by ggplot

corF <- cor(fpkm)
corF
##           c1        c2        t1        t2
## c1 1.0000000 0.9636821 0.8379175 0.8543578
## c2 0.9636821 1.0000000 0.8391606 0.8452566
## t1 0.8379175 0.8391606 1.0000000 0.9936924
## t2 0.8543578 0.8452566 0.9936924 1.0000000
##load package
library(reshape2)
melted_corF <- melt(corF) # change the correlation matrix to the pair 
length(corF)
## [1] 16
nrow(melted_corF)
## [1] 16
#heatmap drawing
library(ggplot2)
ggplot(data = melted_corF, aes(x=Var1, y=Var2, fill=value)) + 
  geom_tile()     

# Get lower triangular of the correlation matrix
get_upper_tri<-function(cormat){
  cormat[lower.tri(cormat)] <- NA
  return(cormat)
}

# Only keep lower triangular
upper_tri        <- get_upper_tri(corF)
melted_upper_tri <- melt( upper_tri )


ggplot(data = melted_upper_tri, aes(Var2, Var1, fill = value))+
  geom_tile(color = "white")+
  scale_fill_gradient2(low = "navy", high = "red", mid = "white", 
                       midpoint = 0.9, limit = c(0.8,1), 
                       #space = "Lab", 
                       name="Pearson\nCorrelation") +
  theme_minimal()+ 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, 
                                   size = 12, hjust = 1))+
  coord_fixed()