[R] partial correlations
Ko-Kang Wang
Ko-Kang at xtra.co.nz
Wed Aug 1 12:45:42 CEST 2001
I got this from Associate Professor Brian McArdle's multivariate statistics course. However I have not tested it yet...
parcor.test <- function( x, y, q=1, alternative = "two.sided" ) {
if( is.matrix( x ) && ncol( x ) >= 2 ) {
y <- x[, 2]
x <- x[, 1]
}
if( length( x ) != length( y ) )
stop( "x and y should be the same length " )
n <- length( x )
if( n <= 2 )
stop( "x and y should effectively be longer than 2" )
# Pearson's product moment correlation:
coef <- cor( x, y )
if( is.na( coef ) )
stop( "too small variance" )
if( abs( coef ) > 0.999999 )
stop( paste( "correlation is", coef ) )
attr( coef, "names" ) <- "cor"
stat <- coef * sqrt( ( n - 2-q )/( 1 - coef^2 ) ) # t-statistic
attr( stat, "names" ) <- "t"
p.value <- if( names( stat ) == "t" ) {
switch( alternative,
greater = 1 - pt( stat, n - 2 - q ),
less = pt( stat, n - 2 - q ),
two.sided = 2 * pt( - abs( stat ), n - 2 - q ),
)
}
attr( p.value, "Names" ) <- NULL #a kludge for nicer output
if( names( stat ) == "t" ) {
pars <- n - 2 - q
attr( pars, "names" ) <- "df"
}
else pars <- NULL
null.value <- 0
attr( null.value, "names" ) <- c( "coef" )
z <- list( statistic = stat, parameters = pars, p.value =
p.value, estimate = coef, null.value = null.value,
alternative = alternative )
z$data.name <- paste( deparse( substitute( x ) ), "and", deparse(
substitute( y ) ) )
attr( z, "class" ) <- "htest"
return( z )
}
--------------------------------------------------------------------------------------------
Ko-Kang Kevin Wang
Head of Statistical Analysis Division
Software Developers' Klub (SDK)
University of Auckland
New Zealand
-----Original Message-----
From: owner-r-help at stat.math.ethz.ch [mailto:owner-r-help at stat.math.ethz.ch]On Behalf Of Ragnar Beer
Sent: Wednesday, August 01, 2001 9:46 PM
To: r-help at stat.math.ethz.ch
Subject: [R] partial correlations
Howdy!
I need to calculate partial correlations and I just can't find out how to
do that with R. Can anybody help?
Ragnar
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
More information about the R-help
mailing list