# Getting weird values for a statistical test in python

I'm trying to perform a Mann-Kendall Test of some data. I'm using the code in the following link (https://github.com/mps9506/Mann-Kendall-Trend/blob/master/mk_test.py) modified a little bit bellow for the result to be in an array and it only gives back a p-value (p) and a tau value (z).

def mk_test(x, alpha=0.05):
"""
This function is derived from code originally posted by Sat Kumar Tomer
([email protected])
The purpose of the Mann-Kendall (MK) test (Mann 1945, Kendall 1975, Gilbert
1987) is to statistically assess if there is a monotonic upward or downward
trend of the variable of interest over time. A monotonic upward (downward)
trend means that the variable consistently increases (decreases) through
time, but the trend may or may not be linear. The MK test can be used in
place of a parametric linear regression analysis, which can be used to test
if the slope of the estimated linear regression line is different from
zero. The regression analysis requires that the residuals from the fitted
regression line be normally distributed; an assumption not required by the
MK test, that is, the MK test is a non-parametric (distribution-free) test.
Hirsch, Slack and Smith (1982, page 107) indicate that the MK test is best
viewed as an exploratory analysis and is most appropriately used to
identify stations where changes are significant or of large magnitude and
to quantify these findings.
Input:
x:   a vector of data
alpha: significance level (0.05 default)
Output:
trend: tells the trend (increasing, decreasing or no trend)
h: True (if trend is present) or False (if trend is absence)
p: p value of the significance test
z: normalized test statistics

Examples
--------
>>> x = np.random.rand(100)
>>> trend,h,p,z = mk_test(x,0.05)
"""
n = len(x)

# calculate S
s = 0
for k in range(n-1):
for j in range(k+1, n):
s += np.sign(x[j] - x[k])

# calculate the unique data
unique_x = np.unique(x)
g = len(unique_x)

# calculate the var(s)
if n == g:  # there is no tie
var_s = (n*(n-1)*(2*n+5))/18
else:  # there are some ties in data
tp = np.zeros(unique_x.shape)
for i in range(len(unique_x)):
tp[i] = sum(x == unique_x[i])
var_s = (n*(n-1)*(2*n+5) - np.sum(tp*(tp-1)*(2*tp+5)))/18

if s > 0:
z = (s - 1)/np.sqrt(var_s)
#result = (s - 1)/np.sqrt(var_s)
elif s == 0:
z = 0
#result = 0
elif s < 0:
z = (s + 1)/np.sqrt(var_s)
#result = (s + 1)/np.sqrt(var_s)

# calculate the p_value
p = 2*(1-norm.cdf(abs(z)))  # two tail test
result= np.append(p,z)
h = abs(z) > norm.ppf(1-alpha/2)

return np.array(result)


I am then using the following code to perform the test.

out = np.empty((0))
for i in range(145):
for j in range(192):
out1 = mk_test(yrmax[:,i,j], alpha=0.05)
out = np.append(out, out1, axis=0)


I think something is going wrong somewhere while performing the test because I would expect to get z values between -1 and 1, but I'm getting back some values that are greater than 1. Is the coding going wrong or am I misunderstanding what z is and it's not actually tau and hence why I'm getting values I'm not expecting?

• $z$ is not $\tau$. Nov 1, 2017 at 7:15
• Attn close voters: Reopened as the underlying issue is plainly statistical - the answer is an explanation of the difference between $\tau$ and $z(\tau)$, for goodness sake. Nov 1, 2017 at 22:40

The $\tau$ statistic is (in effect) the number of pairs of points showing an increase minus the number of pairs of points showing a decrease all divided by the number of pairs of points showing an increase or decrease (which if there's no ties, is just the number of pairs of points).

Incorporating adjustment for ties in the series (but assuming one value per time instant, so no ties there), $\tau=S/D$, where $D = \sqrt{[{n\choose 2}-\sum_j {t_j\choose 2}\}]{n\choose 2}}$ where $t_j$ is the number of ties at the $j$-th value which is non-unique (and if there are no ties, $D={n\choose 2}$ as mentioned above)

So (out of all pairs of points) it's the proportion of increasing pairs minus the proportion of decreasing pairs.

It's a number between -1 and 1.

A z score would divide by the standard deviation of the statistic (say the difference in counts divide by the standard deviation of the difference in counts or $\tau$ divided by its standard deviation, given the population Kendall correlation is 0), and as a result it can be much larger than 1 in absolute value -- in large samples the upper bound approaches $1.5\sqrt{n}$.

The larger $n$ is the smaller the standard deviation in $\tau$ under the null, so a given $\tau$ may correspond to a small-sized $Z$ at a small sample size and a large-sized $Z$ at a large sample size. That is, the relationship between $Z$ and $\tau$ is changes with sample size.

If the null hypotheses is true the $Z$ value will have mean 0 and standard deviation $1$ and in large samples will be approximately standard normal. In those circumstances it will nearly always be between say -2.5 and 2.5 or so, but when the null is false it may be quite far from 0.

With the Mann-Kendall test it seems to be standard to incorporate the continuity correction directly into the z-score (instead of using the correction when calculating the p-value); this just moves the difference in counts one* closer toward 0 before standardizing.

* why 1, rather than 1/2? Because the counts two sets of counts are perfectly negatively correlated so the difference in counts moves by two (one goes up the other goes down); the continuity correction is then half the smallest change.

• Thank you for the clarification! I've never taken a stats class, so I'm still learning some stuff. Do you have a recommendation coding wise what to do to get tau? Nov 1, 2017 at 16:28
• Okay, let's assume a maximum of one data point at any given instant of time. See in the code where it has calculated S and var_s and then calculates z? In there where it's computing var_s, compute $D$ as in my edited answer -- where my $t_j$ is tp` in the code. Then $\tau=S/D$. You should get the same value from anything that calculates Kendall's tau though, and I bet there are already a bunch of ways to just call something for that. Nov 1, 2017 at 22:31