From: Charles C. Berry <cberry_at_tajo.ucsd.edu>

Date: Mon, 21 Jul 2008 13:51:48 -0700

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Mon 21 Jul 2008 - 20:58:52 GMT

Date: Mon, 21 Jul 2008 13:51:48 -0700

On Mon, 21 Jul 2008, Michela Cameletti wrote:

> Dear R user,

*> I'm trying to find a solution for optimizing my code. I have to run a 50.000
**> iteration long simulation and it is absolutely necessary to have an
**> optimized code.
**>
*

What you have is just a quadratic form, so t( y ) %*% A %*% y will do it.

You need to rearrange your subscripts, this is one approach:

*> X2 <- matrix( aperm(X, c(1,3,2) ),nr=d ) # X is as you defined it.
**> nother.result <- crossprod( matrix(A%*%X2,nc=k), matrix( X2, nc=k) )
*

Using dim( X2 ) <- ... rather than matrix() will probably speed this further.

**HTH,
**
Chuck

p.s. This is the road to perdition. You can optimize the heck out of some piece of linear algebra in R only to find that it needs to be written in C, but the R code that you have written is so hard to read that you have to start from scratch. IIRC, writing block-trace algorithms for mixed models takes one down this path.

> I have to do this operation

*> *sum_t ( t(X_t) %*% A %*% X_t )*
**> where X_t is a (d*k) matrix which changes in time and A is a constant in
**> time (d*d) matrix.
**> I have put all my X_t in a three dimensional array X of dimension (d,k,T).
**>
**> At the moment for computing the sum over time I'm doing a for loop and
**> saving the resulting (k*k) matrix in a list and at the end I sum the T
**> matrices in this list. I'm wondering if there is a better way to do this.
**>
**> Here an example of what I have to do:
**>
**> *d=3
**> k=2
**> T=4
**>
**> X = array(rnorm(d*k*T),dim=c(d,k,T))
**> A = matrix(rnorm(d*d),d,d)
**>
**> e1 = list()
**> for (t in 1:T){ #I would like to avoid this
**> e1[[t]] = t(X[,,t])%*%A%*%X[,,t]
**> }
**>
**> ##############################
**> #Function for doing the sum of matrices in a list
**> ##############################
**>
**> sumMatrices <- function(matrices){
**> if (length(matrices) > 2) matrices[[1]] + Recall(matrices[-1])
**> else matrices[[1]] + matrices[[2]]
**> }
**> ##############################
**>
**> result = sumMatrices(e1)
**> *
**>
**>
**> Thank you in advance for all your help,
**> best
**> Michela
**>
**> [[alternative HTML version deleted]]
**>
**> ______________________________________________
**> R-help_at_r-project.org mailing list
**> https://stat.ethz.ch/mailman/listinfo/r-help
**> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
**> and provide commented, minimal, self-contained, reproducible code.
**>
*

Charles C. Berry (858) 534-2098 Dept of Family/Preventive Medicine E mailto:cberry_at_tajo.ucsd.edu UC San Diegohttp://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Mon 21 Jul 2008 - 20:58:52 GMT

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.2.0, at Mon 21 Jul 2008 - 21:32:49 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*