Hello,
(I realize this mail is a bit lengthy, but I would appreciate it if someone could comment on it).
I believe that I found a bug in your implementation of Kendall's Tau. I
have evaluated the implementation (to verify a self-written
implementation). When the results turned out to be different I
investigated the current SciPy implementation at the following URL:
http://svn.scipy.org/svn/scipy/trunk/scipy/stats/stats.py
(I am aware of the fact that there is also a Kendall's Tau implementation in mstats.py, but
have not evaluated that implementation yet).
I will give some explanation of my interpretation of Kendall's tau, an
example showing the differences between SciPy's and my implementation and a
possible fix for SciPy's implementation.
Your implementation is Kendall's tau-b with tie correction (same as
mine). I take as my reference definition, the one in the following
'poster' paper:
http://portal.acm.org/citation.cfm?id=1277935
(this same definition appears in other places as well, this is the
shortest resource I could find)
Recall that Kendall's tau calculates a score t given two rankings R1 and
R2. Variables P, Q and T are all characteristics of the pairs in those
rankings.
The definition given in the reference is:
t = (P - Q) / SQRT((P + Q + T) * (P + Q + U))
where P is the number of concordant pairs, Q the number of discordant
pairs, T the number ties in R1 and U the number of ties in R2.
An example:
-----------
Let's use two identical rankings with a tie:
A B C
R1 = [1, 1, 2]
R2 = [1, 1, 2]
There are three pair combinations in these lists, namely: (A, B), (A, C)
and (B, C). It is obvious that _one_ of these combinations has a tie for
both lists (the (A,B) combination which is (1,1) for both R1 and R2).
So, since there is one tie in both list we have T = U = 1
We find that there are two concordant pairs in both lists (A, C) and
(B,C) so P = 2. There are no discordant pairs, so Q = 0. With all
variables given, we can now calculate Kendall's tau for R1 and R2:
t = (2 - 0) / SQRT((2 + 0 + 1)*(2 + 0 + 1))
t = 2 / SQRT(3*3)
t = 2 / 3
t = 0.6666666
However, using scipy (svn HEAD) as follows:
import scipy.stats.stats as s
s.kendalltau([1,1,2], [1,1,2])
Yields t = 1.0:
(1.0, 0.11718509694604401)
Which I believe is wrong (or at least: has no correction for ties, as is
claimed in the source code). If there are three combinations and one of
these is a tie, and the other two combinations are concordant, it makes
sense that Kendall's tau-b should yield 2 / 3.
The cause and fix
-----------------
Playing around with SciPy's code (and comparing it with my own) I believe I
discovered a probable cause for this difference in SciPy's code. Again, I used the
implementation at the following URL:
http://svn.scipy.org/svn/scipy/trunk/scipy/stats/stats.py
(please take look at the implementation first, otherwise you will not
understand my explanation)
In the 'kendalltau(x,y)' function we see a test for ties and an 'else'
branch. In the 'else' branch the values of 'n1' and 'n2' are incremented
if there is a tie (conforming to +T and +U in the formula given above).
However, I believe that the 'if' conditions here are wrong:
1) Consider that if 'a1' has value '0' it is tied (the same goes for
'a2'). In the else branch I see:
if a1:
n1 = n1 + 1
if a2:
n2 = n2 + 1
So, here the addition takes places on the variables (n1, n2) if there is
NO tie, instead of if there is a tie. Hence, this explains the different
outcome. Translating this back to the formula gives me T = U = 0, which
would yield:
t = (2 - 0) / SQRT((2 + 0 + 0)*(2 + 0 + 0))
t = 2 / SQRT(2*2)
t = 2 / 2
t = 1.0
Which is indeed consistent with the SciPy outcome. Henceforth, I believe
the solution to this is to correct the condition in the if statements in
the Kendall's tau function:
if not a1:
n1 = n1 + 1
if not a2:
n2 = n2 + 1
Closing
-------
Of course, my interpretation of Kendall's Tau could be wrong. Since I
can not exclude that possibility I would appreciate it if one of you could
check and see if you reach the same conclusion. Maybe the base formula that
SciPy uses is different.
I have compared your implementation also to that implemented in the R
project, however their source code suggests that they do not adjust for
ties (effectively implementing Kendall's tau-a).
--
With kind regards,
Almer S. Tigelaar
University of Twente