Looping over a matrix passed to .C

m <- matrix(1:25,5)
indx <- combn(5, 3)
indx

mxy <- m[cbind(indx[1,], indx[2,])]


mxy

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

I've written a function in R which takes a symmetrical matrix as input and
processes all triplicate combinations of values from the matrix. The function
looks something like:

my_fun <- function(m) {
 if( nrow(mat) != ncol(mat) ) {
   stop("'m' must be a square matrix")
 }

 size <- nrow(m)

 for(x in 1:(size -2)) {
   for(y in (x+1):(size -1)) {
     xy <- m[x,y]

     for(z in (y+1):size ) {
       xz <- m[x,z]
       yz <- m[y,z]

       # do something with xy, xz, yz
     }
   }
 }
}

I'd like to speed this up since when size gets > a few thousand, I estimate it
would take 3yrs to complete the "do something" task. I could implement the "do
something" in C and have it called from within the nested for loops in R, but I
think I should get better performance if I implement the for loops in C as well.
As such, I'm trying to get my head around looping through matrix values when
passed to .C()

As I understand it, once a matrix (n x m in size) is passed to .C() it is seen
as an unwrapped array of length n*m. Could someone help/guide me in implementing
this?

Cheers,
Nathan


- --
- --------------------------------------------------------
Dr. Nathan S. Watson-Haigh
OCE Post Doctoral Fellow
CSIRO Livestock Industries
Queensland Bioscience Precinct
St Lucia, QLD 4067
Australia

Tel: +61 (0)7 3214 2922
Fax: +61 (0)7 3214 2900
Web: http://www.csiro.au/people/Nathan.Watson-Haigh.html
- --------------------------------------------------------

-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.9 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iEYEARECAAYFAkmSPM0ACgkQ9gTv6QYzVL5HbwCfRdA+7madbF5zuoKJRbuh4/tE
sLwAn39h1JxNF6MjaF+AgxCq3XVe0X7T
=KEDP
-----END PGP SIGNATURE-----

Why don't you generate all the possible combinations using 'combn' and
then use that matrix to address your data:

dim(m)

m <- matrix(1:25,5)
indx <- combn(5, 3)
indx

     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1    1    1    1    1    1    2    2    2     3
[2,]    2    2    2    3    3    4    3    3    4     4
[3,]    3    4    5    4    5    5    4    5    5     5

mxy <- m[cbind(indx[1,], indx[2,])]


mxy

 [1]  6  6  6 11 11 16 12 12 17 18

You can use each of the columns as x,y,z values and compute everything at once.

On Tue, Feb 10, 2009 at 9:49 PM, Nathan S. Watson-Haigh
<nathan.watson-haigh at csiro.au> wrote:
I've written a function in R which takes a symmetrical matrix as input and
processes all triplicate combinations of values from the matrix. The function
looks something like:

my_fun <- function(m) {
 if( nrow(mat) != ncol(mat) ) {
   stop("'m' must be a square matrix")
 }

 size <- nrow(m)

 for(x in 1:(size -2)) {
   for(y in (x+1):(size -1)) {
     xy <- m[x,y]

     for(z in (y+1):size ) {
       xz <- m[x,z]
       yz <- m[y,z]

       # do something with xy, xz, yz
     }
   }
 }
}

I'd like to speed this up since when size gets > a few thousand, I estimate it
would take 3yrs to complete the "do something" task. I could implement the "do
something" in C and have it called from within the nested for loops in R, but I
think I should get better performance if I implement the for loops in C as well.
As such, I'm trying to get my head around looping through matrix values when
passed to .C()

As I understand it, once a matrix (n x m in size) is passed to .C() it is seen
as an unwrapped array of length n*m. Could someone help/guide me in implementing
this?

Cheers,
Nathan

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

jim holtman wrote:

Why don't you generate all the possible combinations using 'combn' and
then use that matrix to address your data:

I had thought of using this, however the size of data I'm using will more likely
be orders of magnitude larger, such as:

dim(m)

[1] 10000 10000

Which takes looooong time to compute using:
indx <- combn(10000, 3)

m <- matrix(1:25,5)
indx <- combn(5, 3)
indx

     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1    1    1    1    1    1    2    2    2     3
[2,]    2    2    2    3    3    4    3    3    4     4
[3,]    3    4    5    4    5    5    4    5    5     5

mxy <- m[cbind(indx[1,], indx[2,])]


mxy

 [1]  6  6  6 11 11 16 12 12 17 18

You can use each of the columns as x,y,z values and compute everything at once.

On Tue, Feb 10, 2009 at 9:49 PM, Nathan S. Watson-Haigh
<nathan.watson-haigh at csiro.au> wrote:
I've written a function in R which takes a symmetrical matrix as input and
processes all triplicate combinations of values from the matrix. The function
looks something like:

my_fun <- function(m) {
 if( nrow(mat) != ncol(mat) ) {
   stop("'m' must be a square matrix")
 }

 size <- nrow(m)

 for(x in 1:(size -2)) {
   for(y in (x+1):(size -1)) {
     xy <- m[x,y]

     for(z in (y+1):size ) {
       xz <- m[x,z]
       yz <- m[y,z]

       # do something with xy, xz, yz
     }
   }
 }
}

I'd like to speed this up since when size gets > a few thousand, I estimate it
would take 3yrs to complete the "do something" task. I could implement the "do
something" in C and have it called from within the nested for loops in R, but I
think I should get better performance if I implement the for loops in C as well.
As such, I'm trying to get my head around looping through matrix values when
passed to .C()

As I understand it, once a matrix (n x m in size) is passed to .C() it is seen
as an unwrapped array of length n*m. Could someone help/guide me in implementing
this?

Cheers,
Nathan