I think that the code you propose does exactly the same thing as the OP's code for z normalization...
I have generated 2 heatmaps using heatmap.2 function. Following is my command and it's output. First heatmap contains z-score which is been calculated by heatmap.2 function with parameter "scale=row". While second heatmap I have generated is with values calculated by my own function name "normalisation" which also calculate the Z-score. Surprisingly output is different in both case. I wonder how heatmap.2 calculate the Z-Score. The way I have calculated Z-score is right ? Please help.
Command For First Heatmap --> h1
heatmap.2 (
mat ,
col=colorRampPalette(c("blue","white","red"))(dim(mat)[1]*dim(mat)[2]),
density ="none",
Colv = T,
Rowv=T,
trace ="none",
scale="row",
dendrogram = "both",
main = paste("Cg_Pol_2_Top_2k_Highly_Variable"),#,dim(mat)[1],sep =" "),
adjCol = c(1,0),
margins = c(10,4),
lmat = rbind(c(4,3),c(2,1)),
lwid = c(1.5,4),
lhei = c(0.5,5),
labRow=c(""),
keysize=1,
)
Command for heatmap 2
heatmap.2 (
normalisation(mat),
col=colorRampPalette(c("blue","white","red"))(dim(mat)[1]*dim(mat)[2]),
density ="none",
Colv = T,
Rowv=T,
trace ="none",
scale="none",
dendrogram = "both",
main = paste("Cg_Pol_2_Top_2k_Highly_Variable"),#,dim(mat)[1],sep =" "),
adjCol = c(1,0),
margins = c(10,4),
lmat = rbind(c(4,3),c(2,1)),
lwid = c(1.5,4),
lhei = c(0.5,5),
labRow=c(""),
keysize=1,
)
normalization<-function(x){
dimm=dim(x)
for(i in 1:dimm[1]){
x[i,]=(x[i,]-mean(x[i,]))/sd(x[i,])
}
return(x)
}
2 answers
Better late than never...
The actual code snippet from heatmap.2 is:
if (scale == "row") {
retval$rowMeans <- rm <- rowMeans(x, na.rm = na.rm)
x <- sweep(x, 1, rm)
retval$rowSDs <- sx <- apply(x, 1, sd, na.rm = na.rm)
x <- sweep(x, 1, sx, "/")
}
else if (scale == "column") {
retval$colMeans <- rm <- colMeans(x, na.rm = na.rm)
x <- sweep(x, 2, rm)
retval$colSDs <- sx <- apply(x, 2, sd, na.rm = na.rm)
x <- sweep(x, 2, sx, "/")
}
The standard way to scale (by column) to the Z-scale in R is with the scale function, i.e., subtract the mean and then divide by the standard deviation [of each column] - to do this by row, which is typically related to genes, we would transpose the data, apply scale(), and then transpose back (t(scale(t(x)))), which is ultimately the same as what is happening in the heatmap.2 code. The proof:
generate random matrix of numbers
mat <- matrix(rexp(200, rate=.1), ncol=20)
heatmap.2 functions
na.rm <- TRUE
scaleRow <- function(x) {
rm <- rowMeans(x, na.rm = na.rm)
x <- sweep(x, 1, rm)
sx <- apply(x, 1, sd, na.rm = na.rm)
x <- sweep(x, 1, sx, "/")
return(round(x, 6))
}
scaleCol <- function(x) {
colMeans <- rm <- colMeans(x, na.rm = na.rm)
x <- sweep(x, 2, rm)
colSDs <- sx <- apply(x, 2, sd, na.rm = na.rm)
x <- sweep(x, 2, sx, "/")
return(round(x, 6))
}
check that heatmap.2 functions produce same data as manually using scale()
all((scaleRow(mat) == round(t(scale(t(mat))), 6)) == TRUE)
[1] TRUE
all((scaleCol(mat) == round(scale(mat), 6)) == TRUE)
[1] TRUE
Kevin
If I am correct
z score = (x-mean) / sd
where x is a specific value. In your code x is also the row.
You may try:
with new is a new matrix to avoid using already proccessed data
normalization<-function(x){
new = x
dimm=dim(x)
for(row in 1:dimm[1]){
for( col in 1:dimm[2]){
new [row,col] = x[row,col] - mean(x[row,]) / sd(x[row,])
}
}
return(new)
}
Hope this help.
could you details a bit more please? I may be wrong but it seems to me that the OP's function have two mistakes that I try to correct:
First modifying the matrix in the same time it compute z-score (mean and sd will vary for the same row)
Secondly in
x[i,] = (X[i,] - mean(x[i,]) ) / sd(x[i,] )the X[i,] correspond to the row and not to a specific value? This is why I proposed to interate through collumn too.
In fact you are totally right I miss the fact than R make operation directly on vector and matrix. I apologyze for the mistake Thanks to point it out.
Log in to answer this question.
You use
normalisation(mat)and notnormaliZation(mat)in your code. Is it a typo or could "normalisation" be a different function than the one written below ?