You're right I had N and M defined outside of the function and rook and rsub 
were picking up on that.  The following is a bit better and more cleaned up 
version with i-th order option:

rook = function(Y,i) {
        N = nrow(Y);
        M = ncol(Y);
        rsub = function(Z,i) {
                X = matrix(0,N,M);
                X[1:(N-i),] = X[1:(N-i),] + Z[(1+i):N,];
                X[(1+i):N,] = X[(1+i):N,] + Z[1:(N-i),];
                X[,1:(M-i)] = X[,1:(M-i)] + Z[,(1+i):M];
                X[,(1+i):M] = X[,(1+i):M] + Z[,1:(M-i)];
                return(X);
        }
        return(rsub(Y,i)/rsub(matrix(1,N,M),i));
}

Notice that the variable "i" can be passed any value even one that causes an 
error.

--Brett


-----Original Message-----
From: Mills, Jason [mailto:[EMAIL PROTECTED] 
Sent: Friday, February 17, 2006 4:11 PM
To: Pontarelli, Brett
Subject: RE: [R] Matrix indexing in a loop

Hi Brett, thanks for the tip.

I tried this function on my sample matrix and got the error message "Error in 
rsub(fun) : object "N" not found".  Your code looks like it should work, so I 
must be doing some wrong.  I will continue to experiment.

As for the neighbor pattern, the convention follows the rules of chess.
I would consider a 2nd order rook case to include only 4 elements:

a[i-2,j]
a[i+2,j]
a[i,j-2]
a[i,j+2]

Even 3rd order Rook would still only include 4 elements.  As another example, a 
Queen neighborhood includes 8 elements, independent of the lag order.


I haven't started using 3D arrays yet.  I am attempting spatio-temporal 
analysis and I have thought about representing my landscape (two
dimensions) over time (the third dimension).  For now, I'm trying to get a 
handle on working in two dimensions.

Thanks.  

Jason






-----Original Message-----
From: Pontarelli, Brett [mailto:[EMAIL PROTECTED]
Sent: Friday, February 17, 2006 4:04 PM
To: Mills, Jason; [email protected]
Subject: RE: [R] Matrix indexing in a loop

Do you have to use a loop?  The following function should do what you want for 
the 1st order:

rook = function(Y) {
        rsub = function(Z) {
                X = matrix(0,nrow(Z),ncol(Z));
                X[1:(N-1),1:M] = X[1:(N-1),1:M] + Z[2:N,1:M];
                X[2:N,1:M] = X[2:N,1:M] + Z[1:(N-1),1:M];
                X[1:N,1:(M-1)] = X[1:N,1:(M-1)] + Z[1:N,2:M];
                X[1:N,2:M] = X[1:N,2:M] + Z[1:N,1:(M-1)];
                return(X);
        }
        return(rsub(Y)/rsub(matrix(1,nrow(Y),ncol(Y))));
}

I'm not sure I understand how the higher orders work.  For example, an interior 
element for the 1st order is always divided by 4.  Is an interior element for a 
3rd order divided by 4 or 8 or something else?
Also, how are you implementing your 3D matrices?

--Brett


-----Original Message-----
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Mills, Jason
Sent: Friday, February 17, 2006 1:36 PM
To: [email protected]
Subject: [R] Matrix indexing in a loop


How do you specify matrix location a[i,j] (or a[i-1,j], etc.) in a "for"
loop?  

I am looking for a flexible method of indexing neighbors over a series of lags 
(1,2,3...) and I may wish to extend this method to 3D arrays.


Example:

Data matrix
> fun
     [,1] [,2] [,3]
[1,]    1    5    9
[2,]    2    6   10
[3,]    3    7   11
[4,]    4    8   12


For each element a[i,j] in "fun", sum the 1st order (Rook's) neighbors:

a[i-1,j]

a[i+1,j]

a[i,j-1]

a[i,j+1]

Then divide by the number of elements included as neighbors-- this number 
depends on the location of a[i,j] in the matrix.


Insert the product of the neighbor calculation for each a[i,j] into the 
corresponding position b[i,j] in an empty matrix with the same dimensions as 
"fun".


For example, element [2,2] in "fun" should yield element [2,2] in a new matrix 
equal to 24/4=6.  Of course, element [1,1] in the new matrix should be the 
product of only two numbers.


Thanks

J. Mills

        [[alternative HTML version deleted]]

______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html

______________________________________________
[email protected] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

Reply via email to