[R] array problem

Gabor Grothendieck ggrothendieck at myway.com
Sun Jan 4 17:51:47 CET 2004



Replace your line that updates A with this:

p <- unique(c(i,j,1:5))
f <- function(x) diag( matrix(x,4,4) )
AA <- apply( outer(A,G/B), p[-(1:2)], f )   # form product
A <- aperm( array( AA, dim(A) ), order(p) ) # reshape

Here are a couple of tests.  You might want to do some
more tests yourself as well since these are the only
ones I did:

> 
> i <- 3; j <- 5
> G <- 100 * matrix(1:4,2)
> B <- 1+0*G  # all ones
> A <- array(1:32,c(2,2,2,2,2))
> A2 <- A
> A2[,,1,,1] <- A2[,,1,,1] * G[1,1]
> A2[,,1,,2] <- A2[,,1,,2] * G[1,2]
> A2[,,2,,1] <- A2[,,2,,1] * G[2,1]
> A2[,,2,,2] <- A2[,,2,,2] * G[2,2]
> 
> test <- function(A,i,j) {
+ p <- unique(c(i,j,1:5))
+ f <- function(x) diag( matrix(x,4,4) )
+ AA <- apply( outer(A,G/B), p[-(1:2)], f )   # form product
+ A <- aperm( array( AA, dim(A) ), order(p) ) # reshape
+ A
+ }
> 
> identical(test(A,i,j),A2)
[1] TRUE
> 
> 
> i <- 2; j <- 4
> A3 <- A
> A3[,1,,1,] <- A3[,1,,1,] * G[1,1]
> A3[,1,,2,] <- A3[,1,,2,] * G[1,2]
> A3[,2,,1,] <- A3[,2,,1,] * G[2,1]
> A3[,2,,2,] <- A3[,2,,2,] * G[2,2]
> 
> identical(test(A,i,j),A3)
[1] TRUE
> 




--- 
Date: Sun, 04 Jan 2004 21:42:55 +0800 
From: Z P <nusbj at hotmail.com>
To: <r-help at stat.math.ethz.ch> 
Subject: [R] array problem 

 
 
Dear all,

I define , for n=5 or any integer greater than 0.

A<-array((1/2)^n , c(rep(2,n)))

then for any i not equal to j, and 1<=i,j<=n,

B<-apply(a,c(i,j),sum)

now B is a 2 by 2 matrix, I also define another costant 2 by 2 matrix G,

How can I change the values of each elements of array A, according the rule 
that,

for example, i=3,j=5,

A[i1,i2,m,i4,l]<-A[i1,i2,m,i4,l]*G[m,l]/B[m,l] , where m,l=1,2 and 
i1,i2,i4=1,2

I can control this given any i and j, however, I must do the iteration

for i in 1:(n-1) {
for j in (i,n)
{B<-apply(a,c(i,j),sum)
here change the value of every elements of A according to my rule
}
}

Is there any easy way to change the value of A in the iteration? Thank you.




More information about the R-help mailing list