Compare Fisher transformation to rdist

rdist_test.py

#!/usr/bin/python
#emacs: -*- mode: python-mode; py-indent-offset: 4; tab-width: 4; indent-tabs-mode: nil -*- 
#ex: set sts=4 ts=4 sw=4 noet:

"""The target is to compare Fisher transformation to simply taking std values from rdist
"""

import numpy as np
from scipy.stats import rdist, norm, t

for p in [0.5, 0.1, 0.05, 0.01]:
  print "original p=%.2f" % p
  for ss_log2 in xrange(1, 10):           # log2 of sample size (ss)
    ss = 2**ss_log2                     # sample size
    rd = rdist(ss-1, 0, 1)              # rdist distribution sample

    r_p = rd.isf(p)

    t_p = r_p / np.sqrt((1-r_p**2)/(ss-2))

    # Lets convert to Fisher score and figure out what p-value that
    # would be

    Fr = 1/2.*np.log((1+r_p)/(1.-r_p))
    z = np.sqrt(ss - 3) * Fr
    print " ss=%3i r(p<=%.3g rdist)=%5.3f rdist.sf(r)=%.3g z(Fisher tr)=%5.3f sf(z)=%5.3f t=%5.3f sf(t)=%5.3f" \
          % (ss, p, r_p, rd.sf(r_p), z, norm().sf(z), t_p, t(ss).sf(t_p))