Re: [Rd] Bug in cor function (PR#7689)

From: Peter Dalgaard <p.dalgaard_at_biostat.ku.dk>
Date: Sun 13 Feb 2005 - 20:12:36 EST

cig69410@syd.odn.ne.jp writes:

> I can't hardly accept the result of cor function with
> pairwize.colplete.obs or complete.obs
>
> insert print statements in cor function,
>
>
> + if (method != "pearson") {
> + Rank <- function(u) if (is.matrix(u))
> + apply(u, 2, rank, na.last = "keep")
> + else rank(u, na.last = "keep")
> + x <- Rank(x)
> + print(x) # add
> + if (!is.null(y)) {
> + y <- Rank(y)
> + print(y) # add
> + }
> + }
> + .Internal(cor(x, y, na.method, method == "kendall"))
>
> and, data is
> > x <- c(7, 9, 8, 0, NA, NA)
> > y <- c(2, 3, 4, NA, 4, 3)
>
> and, call cor function
> > cor(x, y, use="pair", method="sp")
>
> order of x, and y are
> [1] 2 4 3 1 NA NA
> [1] 1.0 2.5 4.5 NA 4.5 2.5
>
> alas!! and the result is
> [1] 0.4271211
>
> oh! no!!
>
> the result must be 0.5

And which part of the following did you fail to understand?

      For 'cov()', a non-Pearson method is unusual but available for
     the sake of completeness.  Note that '"spearman"' basically
     computes 'cor(R(x), R(y))' (or 'cov(.,.)') where 'R(u) := rank(u,
     na.last="keep")'. Notice also that the ranking is (currently) done
     removing only cases that are missing on the variable itself,  which
     may not be what you expect if you let 'use' be '"complete.obs"' or
     '"pairwise.complete.obs"'.

If you have improved code to contribute, you're welcome (notice that this requires reranking for every pair of variables in an n x n correlation matrix in the pairwise case), but there's really not much point in reporting issues that are already known and documented.

-- 
   O__  ---- Peter Dalgaard             Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics     2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark      Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard@biostat.ku.dk)             FAX: (+45) 35327907

______________________________________________
R-devel@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Received on Sun Feb 13 19:24:15 2005

This archive was generated by hypermail 2.1.8 : Sun 13 Feb 2005 - 19:27:29 EST