HiveBrain v1.2.0
Get Started
← Back to all entries
patternpythonMinor

Speeding up nested for-loop analysis code

Submitted by: @import:stackexchange-codereview··
0
Viewed 0 times
analysisloopnestedcodeforspeeding

Problem

I have two matrices, the first one is mat 2000500, the second one is 6000500 tep
I did this code for some analysis to get another matrix that has a specific elements.

Problem is :
This code is very very slow, I need to make faster.Any suggestion to do so is welcome

mat=matrix(sample(c(0,1,2),200,replace =T),nrow=500,ncol=2000)
tep =matrix(rnorm(200,5,2),nrow=500,ncol=6000)`


The code is trying to create one matrix (500 row by 2000 col). Its elements depend on the following :
Let the new matrix is A .

if mat[1,1] = 0 then A[1,1]=tep[1,1]/(sum(tep[1,1:3]))
if mat[1,1] = 1 then A[1,1]=tep[1,2]/(sum(tep[1,1:3]))
if mat[1,1] = 2 then A[1,1]=tep[1,3]/(sum(tep[1,1:3]))
if mat[1,2] = 0 then A[1,2]=tep[1,4]/(sum(tep[1,4:6]))
if mat[1,2] = 1 then A[1,2]=tep[1,5]/(sum(tep[1,4:6])) 
if mat[1,2] = 2 then A[1,2]=tep[1,6]/(sum(tep[1,4:6]))`
......


and so on..
The one column in mat is corresponding to three columns in tep.
Every three values in three columns in matrix tep at a j row is a group we calculate from each group there will be one value according to the example above, so we will end up with a 500 row by 2000 column matrix.

allele=c(0,1,2)
A=matrix(nrow=500,ncol=2000)
 start=seq(1,6000,3)
for(k in 1:500){
for(i in 1:2000){
for(j in start){
  temp.z= mat[k,i]
  temp.pl=tep[k,j:(j+2)]
  loc=which(allele==temp.z)
  temp.pl=temp.pl[loc]/(sum(temp.pl))
  A[k,i]=temp.pl
  rm(temp.pl,temp.z,loc)

      }
   }
}

Solution

Try the following:

mat0 <- tep[, c(TRUE, FALSE, FALSE)]
mat1 <- tep[, c(FALSE, TRUE, FALSE)]
mat2 <- tep[, c(FALSE, FALSE, TRUE)]

denominator <- mat0 + mat1 + mat2
numerator   <- mat0 * (mat - 1) * (mat - 2) / (0 - 1) / (0 - 2) +
               mat1 * (mat - 0) * (mat - 2) / (1 - 0) / (1 - 2) +
               mat2 * (mat - 0) * (mat - 1) / (2 - 0) / (2 - 1)

A <- numerator / denominator

Code Snippets

mat0 <- tep[, c(TRUE, FALSE, FALSE)]
mat1 <- tep[, c(FALSE, TRUE, FALSE)]
mat2 <- tep[, c(FALSE, FALSE, TRUE)]

denominator <- mat0 + mat1 + mat2
numerator   <- mat0 * (mat - 1) * (mat - 2) / (0 - 1) / (0 - 2) +
               mat1 * (mat - 0) * (mat - 2) / (1 - 0) / (1 - 2) +
               mat2 * (mat - 0) * (mat - 1) / (2 - 0) / (2 - 1)

A <- numerator / denominator

Context

StackExchange Code Review Q#67460, answer score: 3

Revisions (0)

No revisions yet.