Как получить заднюю часть разницы с помощью MCMCpack?

Question

Я пытаюсь получить распределение заднего использования MCMCpack разницы между двумя показателями конверсии, сродни A и B вместе участок this PyMC tutorial.

я могу получить апостериорные два выбранных скоростей просто отлично, но я изо всех сил пытаюсь реализовать пробоотборную дельта. Любые идеи?

Edit Истинная дельта (что было бы неизвестно, если бы мы не сфабриковал данные и то, что мы хотим, чтобы оценить с помощью MCMC) разница между этими двумя показателями true_p_a и true_p_b т.е. 0,01.

# define true success rates 
true_p_a = 0.05 
true_p_b = 0.04 

# set sample sizes 
n_samples_a = 1000 
n_samples_b = 1000 

# fabricate some data 
set.seed(10); 
obs_a = rbinom(n=n_samples_a, size=1, prob=true_p_a) 
set.seed(1); 
obs_b = rbinom(n=n_samples_b, size=1, prob=true_p_b) 

# what are the observed conversion rates? 
mean(obs_a) #0.056 
mean(obs_b) #0.042 

# convert to number of successes 
successes_a = sum(obs_a) #56 
successes_b = sum(obs_b) #42 

# calculate the posterior 
require(MCMCpack) 

simulations = 20000 

posterior_a = MCbinomialbeta(successes_a ,n_samples_a, alpha=1, beta=1,mc=simulations) 
posterior_b = MCbinomialbeta(successes_b ,n_samples_b, alpha=1, beta=1,mc=simulations) 

posterior_delta = ????           

posterior_density_a = density(posterior_a) 
posterior_density_b = density(posterior_b) 


# plot the posteriors 
require(ggplot2) 
ggplot() + 
    geom_area(aes(posterior_density_a$x, posterior_density_a$y), fill="#7ad2f6", alpha=.5) + 
    geom_vline(aes(xintercept=.05), color="#7ad2f6", linetype="dotted", size=2) + 
    geom_area(aes(posterior_density_b$x, posterior_density_b$y), fill="#014d64", alpha=.5) + 
    geom_vline(aes(xintercept=.04), color="#014d64", linetype="dotted", size=2) + 
    scale_x_continuous(labels=percent_format(), breaks=seq(0,0.1, 0.01)) 

Answer 1

Вы боретесь только потому, что вы не полностью приняли байесовского мышления. Все в порядке, у меня было много одинаковых концептуальных проблем, когда я начинал. (Этот вопрос устарел, поэтому вы, вероятно, уже поняли это).

Байесовская задняя плотность включает в себя всю доступную информацию о распределении параметров модели. Таким образом, чтобы вычислить функцию любого из параметров модели, вы просто вычисляете эту функцию для каждой ничьей из заднего распределения. Вам не нужно беспокоиться о стандартных ошибках и асимптотических выводах, потому что у вас уже есть вся необходимая информация.

В этом случае, поскольку разница между параметрами является константой, и у вас много данных, неопределенности относительно дельта мало. Он оценивается в среднем 0,014 с SD (не SE) 0,009.

Ваш код с готовым анализом:

# define true success rates 
    true_p_a = 0.05 
    true_p_b = 0.04 

    # set sample sizes 
    n_samples_a = 1000 
    n_samples_b = 1000 

    # fabricate some data 
    set.seed(10); 
    obs_a = rbinom(n=n_samples_a, size=1, prob=true_p_a) 
    set.seed(1); 
    obs_b = rbinom(n=n_samples_b, size=1, prob=true_p_b) 

    # what are the observed conversion rates? 
    mean(obs_a) #0.056 
    mean(obs_b) #0.042 

    # convert to number of successes 
    successes_a = sum(obs_a) #56 
    successes_b = sum(obs_b) #42 

    # calculate the posterior 
    require(MCMCpack) 

    simulations = 20000 

    posterior_a = MCbinomialbeta(successes_a ,n_samples_a, alpha=1, beta=1,mc=simulations) 
    posterior_b = MCbinomialbeta(successes_b ,n_samples_b, alpha=1, beta=1,mc=simulations) 

    # Subtract the posterior deltas, look at empirical summaries and plot the empirical density function 

    posterior_delta = posterior_a - posterior_b 

    summary(posterior_delta) 

    require(ggplot2) 

    ggplot(data.frame(delta=as.numeric(posterior_delta)),aes(x=delta)) + geom_density() + theme_minimal()