Kaspar Pflugshaupt
2000-Feb-25 16:01 UTC
[R] Summary: Partial correlation coefficients in R. Thanks everybody!
Hello all, here's a collection of answers I got on my question concerning partial correlation coefficients: Some people gave a simple formula for the three-variable-case, as did Dave Lucy: pcor <- function(v1, v2, v3) { c12 <- cor(v1, v2) c23 <- cor(v2, v3) c13 <- cor(v1, v3) partial <- (c12-(c13*c23))/(sqrt(1-(c13^2)) * sqrt(1-(c23^2))) return(partial) } The general algorithm was given by John Logsdon, among others: The partial correlation coefficients are the negative scaled off-diagonal inverse variance. So compute the variance-covariance matrix, invert, scale to the diagonal and negate and you have it. And an implementation by Martyn Plummer (here, too, I received several): pcor2 <- function(x){ conc <- solve(var(x)) resid.sd <- 1/sqrt(diag(conc)) pcc <- - sweep(sweep(conc, 1, resid.sd, "*"), 2, resid.sd, "*") return(pcc) } This is the version I'm using now, together with a test for significance of each coefficient (H0: coeff=0): f.parcor <- function (x, test = F, p = 0.05) { nvar <- ncol(x) ndata <- nrow(x) conc <- solve(cor(x)) resid.sd <- 1/sqrt(diag(conc)) pcc <- -sweep(sweep(conc, 1, resid.sd, "*"), 2, resid.sd, "*") colnames(pcc) <- rownames(pcc) <- colnames(x) if (test) { t.df <- ndata - nvar t <- pcc/sqrt((1 - pcc^2)/t.df) pcc <- list(coefs = pcc, significant = t > qt(1 - (p/2), df = t.df)) } return(pcc) } Thanks to everybody for your help! Kaspar -- Kaspar Pflugshaupt Geobotanisches Institut Zuerichbergstr. 38 CH-8044 Zuerich Tel. ++41 1 632 43 19 Fax ++41 1 632 12 15 mailto:pflugshaupt at geobot.umnw.ethz.ch privat:pflugshaupt at mails.ch http://www.geobot.umnw.ethz.ch -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Possibly Parallel Threads
- partial correlation coefficients in R?
- Re: [R] R 1.1 congrat; undocumented behaviour of recordPlot (PR#578)
- eqscplot() in library MASS: fails when given only one point (PR#1164)
- eqscplot() in library MASS: fails when given only one point (PR#1162)
- list of all objects - just being curious