for equal weights post

equal_weights_monte_carlo.py

from pylab import *
from numpy.random import dirichlet, rand

def _unit_weight(dim):
    return ones(dim) / float(dim)

def _feature_vec(dim, storage = None):
    result = rand(dim)
    result[where(result > 0.5)] = 1.0
    result[where(result <= 0.5)] = 0.0
    return result

def test_ranking(dim, nsamples=10000):
    u = _unit_weight(dim)
    h = dirichlet(ones(dim))

    diff_count = 0
    for i in range(nsamples):
        v = _feature_vec(dim)
        w = _feature_vec(dim)

        u_delta = dot(u, v-w)
        h_delta = dot(h, v-w)
        if sign(u_delta * h_delta) < 0:
            diff_count += 1
    return float(diff_count) / nsamples


def ranking(dim, vec_samples=1000, h_samples=1000):
    data = zeros(h_samples)
    for i in range(h_samples):
        data[i] = test_ranking(dim, vec_samples)
    return (mean(data), std(data))

if __name__=="__main__":
    n_dim = 250
    d = arange(n_dim)
    m = zeros(n_dim)
    s = zeros(n_dim)
    for n in range(n_dim):
        m[n], s[n] = ranking(n)
        print "Up to " + str(n) + " dimensions"
    errorbar(d, m, yerr=s)
    xlabel("number of dimensions")
    ylabel("Error fraction")
    show()