Quantile in Python

quantile.py

#!/usr/bin/env python

import numpy as np


def Quantile(data, q, precision=1.0):
    """
    Returns the q'th percentile of the distribution given in the argument
    'data'. Uses the 'precision' parameter to control the noise level.
    """
    N, bins = np.histogram(data, bins=precision*np.sqrt(len(data)))
    norm_cumul = 1.0*N.cumsum() / len(data)

    return bins[norm_cumul > q][0]


def RunTheQuantileTest():
    """
    Run a test on normally distributed data points to see if we get reasonable
    quantiles.
    """
    from matplotlib import pyplot as plt

    data = np.random.normal(size=2000000)

    q01 = Quantile(data, 0.01)
    q99 = Quantile(data, 0.99)

    print "error in  1% quantile", ((1.0*(data < q01).sum() / len(data)) - 0.01)
    print "error in 99% quantile", ((1.0*(data < q99).sum() / len(data)) - 0.99)
    print q01, np.percentile(data,  1.0)
    print q99, np.percentile(data, 99.0)

    plt.axvline(q01)
    plt.axvline(q99)

    plt.hist(data, bins=200)
    plt.show()


if __name__ == "__main__":
    RunTheQuantileTest()