[R] Partial Correlation

[Moumita Das]

Why do we get Partial correlation values greater than 1?

 I have used the default function pcor.mat :--

I have manipulated the default pcor.mat function a bit so ignore tha
variables corr_type,element1_in_no,element2_in_no,P.Please ignore the
“pairwise” section and have a look at athe “listwise ” part i.e else part.

*pcor.mat <-
function(x,y,z,method="p",na.rm=T,corr_type,element1_in_no,element2_in_no,P){
*

*                *

*                print("pcor.mat")*

*                x <- c(x)*

*                y <- c(y)*

*                z <- as.data.frame(z)*

*                print(z)*

*                #print("element1_in_no")*

*                #print(element1_in_no)*

*                *

* *

*                if(dim(z)[2] == 0){           *

*                                stop("There should be given data\n")*

*                }*

* *

*                data <- data.frame(x,y,z)*

*                *

*                *

*                if(corr_type=="pairwise")*

*                {*

*                                *

*                                print("inside pairwise")*

*                                rxx.z
<-P[as.numeric(element1_in_no),as.numeric(element2_in_no)]*

*                                #print("rxx.z")*

*                                #print(rxx.z)*

*                                *

*                                return(rxx.z)*

* *

*                }*

*                else*

*                {*

*                                print("inside listwise")*

*                                if(na.rm == T){*

*                                                data = na.omit(data)*

*                                }*

*                                *

*                                xdata <- na.omit(data.frame(data[,c(1,2)]))
#i1,C1*

*                                print("printing
xdata.......................................")*

*                                print(xdata)*

*                                Sxx <- cov(xdata,xdata,m=method)*

*
print("Sxx...................................................")*

*                                print(Sxx )*

* *

*                                xzdata <- na.omit(data)*

*                                xdata <- data.frame(xzdata[,c(1,2)])*

*                                zdata <- data.frame(xzdata[,-c(1,2)])*

*                                print("zdata..........................")*

*                                print(zdata)*

*                                Sxz <- cov(xdata,zdata,m=method)*

*
print("Sxz.....................................
")*

*                                print(Sxz)*

* *

*                                zdata <-
na.omit(data.frame(data[,-c(1,2)]))*

*                                Szz <- cov(zdata,zdata,m=method)*

*                                print("Szz................................
")*

*                                print(Szz)*

*                                #print("new type par corr")*

*                                #P<-partialCorr_matrix(data)*

*                                *

*                }*

*                *

*                *

*                # is Szz positive definite?*

*                zz.ev <- eigen(Szz)$values*

*                if(min(zz.ev)[1]<0){*

*                                *

*                                stop("\'Szz\' is not positive definite!\n")
*

*                }*

*                *

*                *

*                # partial correlation*

*                Sxx.z <- Sxx - Sxz %*% solve(Szz) %*% t(Sxz)*

*                print("Sxx.z")*

*                print(qr(Sxx.z))*

*                rxx.z <- cov2cor(Sxx.z)[1,2]*

*                *

*                return(rxx.z)*

*}*

*Cov2cor function:-----------------------------------------*

*cov2cor<-function (V) *

*{*

*   print("inside cov2cor")*

*  *

*   *

*   p <- (d <- dim(V))[1]*

*    if (!is.numeric(V) || length(d) != 2L || p != d[2L]) *

*        stop("'V' is not a square numeric matrix")*

*    Is <- sqrt(1/diag(V))*

*    print("Is")*

*    print(Is)*

*    if (any(!is.finite(Is))) *

*        warning("diag(.) had 0 or NA entries; non-finite result is
doubtful")*

*    r <- V*

*    r[] <- Is * V * rep(Is, each = p)*

*     print("r")*

*     print(r[])*

*     *

*    r[cbind(1L:p, 1L:p)] <- 1*

*    r*

*    *

*}*

* *

Sxx , Sxz , Szz  all these three values I have calculated and they match
with SPSS result.

After that I am not understanding why my partial correlation result doesn’t
match with the SPSS result.

Sxx.z <- Sxx - Sxz %*% solve(Szz) %*% t(Sxz)  # how to calculate in
SPSS.After this as you can see in pcor.mat function rxx.z <-
cov2cor(Sxx.z)[1,2] is computed.



Don’t know where things are going wrong.**

*I will attach the datas of categories and items which form the categories.*

*Category 1 has items -1,5 *

*Category 2 has items – 2,4,6,7,8,9,11*

*Category 3 has items -3,12,14*

*Category 4 has items -10,13,15*

* Category values are actually mean of the items which form the category.*

Here in pcor.mat function i passed a item-category pair dataset and rets all
items and categories as predictors.

For Partial pairwise i can understand this function won't work ,because some
items which form a particular category ,cannot be partailly correlated to
that catgory.

but listwise shuld have any problem.
**

* *Is it due to LINPACK or LAPACK package.I don't have LAPACK package,which
is used in teh solve.default package.

*
*

*
*

Thanks

Moumita*
*

* *

* *