[R] *Building* a covariance matrix efficiently
rex.dwyer at syngenta.com
rex.dwyer at syngenta.com
Mon Mar 14 17:43:16 CET 2011
Tjerk,
This is just a pseudo code outline of what you need to do:
M = matrix(0, number of variables, number of variables)
V = rep(0, number of variables)
N = 0
While (more observations to read) {
X <- next observation
V <- V + X
M <- M + outer(X,X)
N <- N+1
}
Compute covariance matrix from elements of V,M, and N
You just need to refer to the formula defining covariance.
Outlook seems to think all my variables should be upper case.
HTH
Rex
-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On Behalf Of Tsjerk Wassenaar
Sent: Monday, March 14, 2011 10:14 AM
To: R-help
Subject: [R] *Building* a covariance matrix efficiently
deaRs,
I want to build a covariance matrix out of the data from a binary
file, that I can read in chunk by chunk, with each chunk containing a
single observation vector X. I wonder how to do that most efficiently,
avoiding the calculation of the full symmetric matrices XX'. The
trivial non-optimal approach boils down to something like:
Q <- matrix(rnorm(100000),ncol=200)
M <- matrix(0,ncol=200,nrow=200)
for (i in 1:nrow(Q))
M <- M + tcrossprod(Q[i,])
I would appreciate pointers to help me fill this lacuna in my R skills :)
Cheers,
Tsjerk
--
Tsjerk A. Wassenaar, Ph.D.
post-doctoral researcher
Molecular Dynamics Group
* Groningen Institute for Biomolecular Research and Biotechnology
* Zernike Institute for Advanced Materials
University of Groningen
The Netherlands
______________________________________________
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.
message may contain confidential information. If you are not the designated recipient, please notify the sender immediately, and delete the original and any copies. Any use of the message by you is prohibited.
More information about the R-help
mailing list