Rejection sampler for the epistemic trust model of Shafto, Eaves, et al (2012)

etrs.py

"""
Copyright (C) 2015 Baxter Eaves
License : WTFPL V2

Rejection sampler for the epistemic trust model of Shafto, Eaves, et al (2012)
"""


import random

N_HYPOTHESIS = 4  # Arbitrary choice from Shafto, Eaves et al (2012)


def get_theta_k(alpha_k, beta_k):
    theta_k = random.betavariate(alpha_k, beta_k)
    return theta_k


def get_knowledgeable(theta_k):
    is_knowledgeable = random.random() < theta_k
    return is_knowledgeable


def get_theta_h(alpha_h, beta_h):
    theta_h = random.betavariate(alpha_h, beta_h)
    return theta_h


def get_helpful(theta_h):
    is_knowledgeable = random.random() < theta_h
    return is_knowledgeable


def get_belief(world, is_knowledgeable):
    if is_knowledgeable:
        belief = world
    else:
        belief = random.randrange(N_HYPOTHESIS)

    return belief


def get_action(belief, is_helpful):
    if is_helpful:
        action = belief
    else:
        actions = [a for a in range(N_HYPOTHESIS) if a != belief]
        action = random.choice(actions)

    return action


def sample_thetas(worlds, actions, params_k=(.5, .5,), params_h=(.5, .5),
                  n_samples=10000):
    """ Draw from the distribution f(\theta_k, \theta_h | observations, prior)

    Parameters
    ----------
    worlds : list<int>
        list of the hypotheses (true states of the world) observed
    actions : list <int>
        list of the observed actions
    params_k : tuple(alpha_k, beta_k,), optional
        Prior prameters on theta_k
    params_h : tuple(alpha_h, beta_h,), optional
        Prior prameters on theta_h
    n_samples : int, optional
        The number of samples to collect

    Returns
    -------
    thetas : list<tuple(theta_k, theta_h)>
        list of samples
    """
    thetas = []
    # keep drawing samples until we have enough
    while len(thetas) < n_samples:

        # thetas persist across trials (demonstrations)
        theta_k = get_theta_k(*params_k)
        theta_h = get_theta_h(*params_h)

        match = True
        # for each observation...
        for world_obs, action_obs in zip(worlds, actions):
            is_knowledgeable = get_knowledgeable(theta_k)
            is_helpful = get_helpful(theta_h)

            # we can fix the state of the world, becuase P(w) is uniform
            belief = get_belief(world_obs, is_knowledgeable)
            action = get_action(belief, is_helpful)

            # If the sampled action does not match the observed action,
            # start over.
            if action != action_obs:
                match = False
                break

        if match:
            thetas.append((theta_k, theta_h,))

    return thetas


def prob_correct_label(thetas, action_obs):
    count = 0.0
    norm = 0.0

    # for each value in our 'function approximation'
    for theta_k, theta_h in thetas:
        is_knowledgeable = get_knowledgeable(theta_k)
        is_helpful = get_helpful(theta_h)

        world = random.randrange(N_HYPOTHESIS)
        belief = get_belief(world, is_knowledgeable)
        action = get_action(belief, is_helpful)

        if action == action_obs:
            norm += 1.0
            if action == world:
                count += 1.0
    return count/norm


def prob_endorse_informant_a(action_a, action_b, thetas_a, thetas_b):
    p_a = prob_correct_label(thetas_a, action_a)
    p_b = prob_correct_label(thetas_b, action_b)

    return p_a / (p_a + p_b)


if __name__ == "__main__":
    N_TRIALS = 3
    worlds = [0]*N_TRIALS
    actions_a = [0]*N_TRIALS
    actions_b = [1]*N_TRIALS

    thetas_a = sample_thetas(worlds, actions_a)
    thetas_b = sample_thetas(worlds, actions_b)

    action_a = 1
    action_b = 0

    pchoosea = prob_endorse_informant_a(action_a, action_b, thetas_a, thetas_b)

    print("P(endorse | \Xi) = %f" % (pchoosea,))