speed up the deviance_residual_transform

[miderxi]

mydeviance_residual_transform <- function(count_mat){
resid_mat <- matrix(nrow = nrow(count_mat), ncol = ncol(count_mat))
#row is cell,coll is gene
row_sum <- rowSums(count_mat) #每个细胞总read
col_sum <- colSums(count_mat) #每个基因总read
total_sum <- sum(count_mat)

for (i in 1:nrow(count_mat)){
U_i <- row_sum[i]col_sum/total_sum
cross_prod <- 2count_mat[i,]log((count_mat[i,]+1e-6)/(U_i+1e-6)) + 2(row_sum[i] - count_mat[i,])*log((row_sum[i]-count_mat[i,])/(row_sum[i]-U_i))
l <- sign(count_mat[i,]-U_i) * sqrt(cross_prod)
resid_mat[i, ] <- l
}
return(resid_mat)
}

dev_res <- mydeviance_residual_transform(t(as.matrix(combined_obj@assays$RNA@counts[1:1000,1:1000])))
dev_res2 <- deviance_residual_transform(t(as.matrix(combined_obj@assays$RNA@counts[1:1000,1:1000])))

max(dev_res -dev_res2)
[1] 1.856293e-13