NB2_alina.Rmd

---
title: "Generating Log2FC and Zscore tables and Heatmaps"
author: "Keshav Prasada Gubbi"
date: "3/29/2022"
output: html_document
---


```{r}
library(tidyverse)
library(EnhancedVolcano)
library(data.table)
library(gplots)
library(RColorBrewer)
```

## log2FC values table - per Gene :

```{r}
binded_frame <- read.csv("D:/R/Rtuts/Data/WG__Discussion_Alina's_EPEC_project/csv/Data_Alina_MasterTable_13strains.csv")
```


```{r}
by_gene_forlog2FC <- binded_frame %>% group_by(symbol) %>% 
  summarise(log2FoldChange, Strain_name) %>% ungroup()
head(by_gene_forlog2FC)
```

```{r}
log2FC_476 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S476") %>%
           summarise(log2FC_476 = paste0(log2FoldChange)) 
log2FC_754 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S754") %>%
           summarise(log2FC_754 = paste0(log2FoldChange))
log2FC_755 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S755") %>%
           summarise(log2FC_755 = paste0(log2FoldChange))
log2FC_757 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S757") %>%
           summarise(log2FC_757 = paste0(log2FoldChange))
log2FC_758 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S758") %>%
           summarise(log2FC_758 = paste0(log2FoldChange))
log2FC_760 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S760") %>%
           summarise(log2FC_760 = paste0(log2FoldChange))
log2FC_761 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S761") %>%
           summarise(log2FC_761 = paste0(log2FoldChange))
log2FC_762 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S762") %>%
           summarise(log2FC_762 = paste0(log2FoldChange))
log2FC_763 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S763") %>%
           summarise(log2FC_763 = paste0(log2FoldChange))
log2FC_764 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S764") %>%
           summarise(log2FC_764 = paste0(log2FoldChange))
log2FC_765 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S765") %>%
           summarise(log2FC_765 = paste0(log2FoldChange))
log2FC_766 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S766") %>%
           summarise(log2FC_766 = paste0(log2FoldChange))
log2FC_768 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S768") %>%
           summarise(log2FC_768 = paste0(log2FoldChange))
log2FC_769 <- by_gene_forlog2FC %>% 
           group_by(symbol) %>% 
           filter(Strain_name == "S769") %>%
           summarise(log2FC_769 = paste0(log2FoldChange))

```


```{r}
merged_L2fc1 <- merge(data.frame(log2FC_476, row.names=NULL), data.frame(log2FC_754, row.names=NULL), by = 1, all = TRUE)
merged_L2fc2 <- merge(data.frame(log2FC_755, row.names=NULL), data.frame(log2FC_757, row.names=NULL), by = 1, all = TRUE)
merged_L2fc3 <- merge(data.frame(log2FC_758, row.names=NULL), data.frame(log2FC_760, row.names=NULL), by = 1, all = TRUE)
merged_L2fc4 <- merge(data.frame(log2FC_761, row.names=NULL), data.frame(log2FC_762, row.names=NULL), by = 1, all = TRUE)
merged_L2fc5 <- merge(data.frame(log2FC_763, row.names=NULL), data.frame(log2FC_764, row.names=NULL), by = 1, all = TRUE)
merged_L2fc6 <- merge(data.frame(log2FC_765, row.names=NULL), data.frame(log2FC_766, row.names=NULL), by = 1, all = TRUE)
merged_L2fc7 <- merge(data.frame(log2FC_768, row.names=NULL), data.frame(log2FC_769, row.names=NULL), by = 1, all = TRUE)

```



```{r}
merged_L2fc_A1 <- merge(data.frame(merged_L2fc1, row.names=NULL), data.frame(merged_L2fc2, row.names=NULL), by = 1, all = TRUE)
merged_L2fc_A2 <- merge(data.frame(merged_L2fc3, row.names=NULL), data.frame(merged_L2fc4, row.names=NULL), by = 1, all = TRUE)
merged_L2fc_A3 <- merge(data.frame(merged_L2fc5, row.names=NULL), data.frame(merged_L2fc6, row.names=NULL), by = 1, all = TRUE)
```

```{r}
merged_L2fc_B1 <- merge(data.frame(merged_L2fc_A1, row.names=NULL), data.frame(merged_L2fc_A2, row.names=NULL), by = 1, all = TRUE)
```

```{r}
merged_L2fc_B2 <- merge(data.frame(merged_L2fc_A3, row.names=NULL), data.frame(merged_L2fc7, row.names=NULL), by = 1, all = TRUE)
```

```{r}
merged_L2fc <- merge(data.frame(merged_L2fc_B1, row.names=NULL), data.frame(merged_L2fc_B2, row.names=NULL), by = 1, all = TRUE)
```

```{r}
write.csv(merged_L2fc,"D:/R/Rtuts/Data/WG__Discussion_Alina's_EPEC_project/csv/L2fc_MergedTable_ByGene.csv")

```

```{r}
head(merged_L2fc)
```


# Determining Zscores and Heatmaps - per gene and per strain!

```{r}
#merged_L2fc_a <- as.data.frame(apply(merged_L2fc, 2, as.numeric))
merged_L2fc <- merged_L2fc %>% distinct(symbol,.keep_all=TRUE)
str(merged_L2fc)
```
```{r}
merged_L2fc <- merged_L2fc %>% mutate_at(c(2:15), as.double)
str(merged_L2fc)
```

this is the crucial part and this is where the problem lies - on how to calculate the mean and also zscore


```{r}
zs_476 = scale(merged_L2fc$log2FC_476, center = TRUE, scale = TRUE)
zs_754 = scale(merged_L2fc$log2FC_754, center = TRUE, scale = TRUE)
zs_755 = scale(merged_L2fc$log2FC_755, center = TRUE, scale = TRUE)
zs_757 = scale(merged_L2fc$log2FC_757, center = TRUE, scale = TRUE)
zs_758 = scale(merged_L2fc$log2FC_758, center = TRUE, scale = TRUE)
zs_760 = scale(merged_L2fc$log2FC_760, center = TRUE, scale = TRUE)
zs_761 = scale(merged_L2fc$log2FC_761, center = TRUE, scale = TRUE)
zs_762 = scale(merged_L2fc$log2FC_762, center = TRUE, scale = TRUE)
zs_763 = scale(merged_L2fc$log2FC_763, center = TRUE, scale = TRUE)
zs_764 = scale(merged_L2fc$log2FC_764, center = TRUE, scale = TRUE)
zs_765 = scale(merged_L2fc$log2FC_765, center = TRUE, scale = TRUE)
zs_766 = scale(merged_L2fc$log2FC_766, center = TRUE, scale = TRUE)
zs_768 = scale(merged_L2fc$log2FC_768, center = TRUE, scale = TRUE)
zs_769 = scale(merged_L2fc$log2FC_769, center = TRUE, scale = TRUE)
```

```{r}
GeneName <- merged_L2fc$symbol #obtaining genenames from earlier df
ZS_table <- data.frame(GeneName, zs_476,zs_754,zs_755,zs_757,zs_758,zs_760,
                       zs_761,zs_762,zs_763,zs_764,zs_765,zs_766,zs_768,zs_769)
```

```{r}
ZS_table <- ZS_table %>% distinct(GeneName,.keep_all=TRUE) #To remove any duplicate or repeated entries, in terms of genes.
```

```{r}
ZS_table <- data.frame(ZS_table, row.names = 1) # convert the genenames column into rownames as its needed for heatmaps of geneexpressions.
```

```{r}
write.csv(ZS_table,"D:/R/Rtuts/Data/WG__Discussion_Alina's_EPEC_project/csv/zscore_table.csv")
```

### Fixing the NA values

The NA values are making it very hard to work with this tables. Replace the NA with the mean of minimum of all strains : -25.56672

```{r}
#mean(
min(ZS_table$zs_476, na.rm = TRUE)
min(ZS_table$zs_754, na.rm = TRUE)
min(ZS_table$zs_755, na.rm = TRUE)
min(ZS_table$zs_757, na.rm = TRUE)
min(ZS_table$zs_758, na.rm = TRUE)
min(ZS_table$zs_760, na.rm = TRUE)
min(ZS_table$zs_762, na.rm = TRUE)
min(ZS_table$zs_761, na.rm = TRUE)
min(ZS_table$zs_763, na.rm = TRUE)
min(ZS_table$zs_764, na.rm = TRUE)
min(ZS_table$zs_766, na.rm = TRUE)
min(ZS_table$zs_768, na.rm = TRUE)
min(ZS_table$zs_769, na.rm = TRUE)
#)
```

```{r}
head(ZS_table)
```

```{r}
#replace all NA values with mean of minimum of values : -25.56672
ZS_table <- ZS_table %>% replace(is.na(.), -25.56672)
head(ZS_table)
```


# Generating a heatmap for Gene expression

### METHOD1: Scaling the table column wise : this is what I found on internet - Calculating the zscore per column.

```{r}
# Extract just the numeric data into a matrix with named rows by gene
geneExp_matrix <- as.matrix(ZS_table)
head(geneExp_matrix)
```

```{r}
heatmap(geneExp_matrix)
tiff("D:/R/Rtuts/Data/WG__Discussion_Alina's_EPEC_project/heatmaps/h1.tif", compression = "lzw")
heatmap(geneExp_matrix)
dev.off()
```


```{r}
heatmap(geneExp_matrix, Rowv=NA,  Colv=NA)
tiff("D:/R/Rtuts/Data/WG__Discussion_Alina's_EPEC_project/heatmaps/h2.tif", compression = "lzw")
heatmap(geneExp_matrix, Rowv=NA,  Colv=NA)
dev.off()
```

```{r}
# Scaled within rows, no column clustering
heatmap.2(geneExp_matrix, col=rev(brewer.pal(9,"RdBu")), scale="row")
```

```{r}

heatmap(geneExp_matrix, col=rev(brewer.pal(9,"RdBu")))
```



```{r}
heatmap(geneExp_matrix,col=brewer.pal(9,"RdGy"))

```



## Method2 : Row-wise

```{r}
transposed_matrix<- t(geneExp_matrix)
```

```{r}
heatmap(transposed_matrix)
```

```{r}
heatmap(transposed_matrix, Rowv=NA,  Colv=NA)
```

```{r}
heatmap(transposed_matrix, col=rev(brewer.pal(9,"RdBu")))
```

```{r}
heatmap.2(transposed_matrix)
```