|
Dear All,
I would like to thank Nick Ellis,
John Fennick, and Bert Gunter who kindly replied to my query about forming
iteratively three-way/four-way tables by combining adjacent rows. I greatly
appreciate their help and suggestions. All strongly suggest to avoid explicit
loops.
In response to 11 requests to post a summary, I am posting the
most helpful and appropriate suggestion. Dr. Ellis' reply has given me
an interesting idea for improving speed and fully take advantage of the
power of S. Here is Dr. Ellis' vectorized program:
You can avoid all
the looping for a two-way table as
follows.
mat<-rbind(c(0,0,1,2),c(1,1,3,0),c(2,0,0,3)) n <-
nrow(mat) adjmat <-
mat[-1,]+mat[-n,] lapply(seq(n-1),function(i) rbind(adjmat[i,],mat[-c(i,i+1),])[order(c(i,seq(1,len=i-1),seq(i+1,len=n-i-1))),]) [[1]]:
[,1] [,2] [,3] [,4] [1,] 1
1 4 2 [2,]
2 0 0
3
[[2]]: [,1] [,2] [,3] [,4]
[1,] 0 0
1 2 [2,] 3
1 3 3
adjmat is a matrix of adjacent
sums. The lapply loops over 1 to n-1 binding the relevant row of adjmat with the
relevant rows of the original matrix, and ordering them
appropriately.
I'm not sure what you mean by extending to 3- or 4-way. If
you mean summing adjacent rows again, then it's quite simple. You use the above
method after converting the 3 or 4 dimensional array to a matrix. That's easily
done by changing the dim attribute.
# 3D case array3 <-
array(1:(4*3*2),dim=c(4,3,2)) # 4x3x2 array dim(array3) <- c(4,3*2) n
<- nrow(array3) adjmat <-
array3[-1,]+array3[-n,] lapply(seq(n-1),function(i) { a
<-rbind(adjmat[i,],array3[-c(i,i+1),])[order(c(i,seq(1,len
=i-1),seq(i+1,len=n-i-1))),] dim(a) <-
c(3,3,2) a }) dim(array3) <- c(4,3,2)
# 4D case array4
<- array(1:(5*4*3*2),dim=c(5,4,3,2)) # 5x4x3x2 array dim(array4) <-
c(5,4*3*2) n <- nrow(array4) adjmat <-
array4[-1,]+array4[-n,] lapply(seq(n-1),function(i) { a
<- rbind(adjmat[i,],array4[-c(i,i+1),])[order(c(i,seq(1,len=i-1),seq(i+1,len=n-i-1))),] dim(a)
<- c(4,4,3,2) a }) dim(array4) <- c(5,4,3,2)
As a general
rule, avoid using for loops; S is very slow. Try to find a way to vectorize your
calculations by treating data as a whole. For example, your code for (j in
1:cols) newmat[index,j] <- mat[index,j] + mat[index+1,j] is the same
as newmat[index,] <- mat[index,] + mat[index+1,]
Ray
Haraf.
|