Hi all!

I need to design a tool to help understand web analytics data in a kind of A/B testing scenario. But it's not a straight A/B test of one piece of content vs another. Instead it's one process for generating a page vs another, and there's many pages on which users can convert/not convert, and some pages may have inherently higher or lower conversion rates. Also while there may be a lot of overall data, when you just look at the page level, many pages have just a few visits, and usually under just condition A or B for a page.

So I researched and found that I should be using a hierarchical bayesian model. I've tried to implement the model in both Python and JAGS, but first it was failing and now I'm not sure it's right. Can you guys look at my model (see below for XXXXXX) and let me know what I'm doing wrong?

Model (Try.bug)

var alpha_A, beta_A, alpha_B, beta_B, p_A[N], p_B[N], delta[N], n_A[N], n_B[N], c_A[N], c_B[N]; model { for (i in 1:N) { c_A[i] ~ dbin(p_A[i],n_A[i]) c_B[i] ~ dbin(p_B[i],n_B[i]) p_A[i] ~ dbeta(alpha_A,beta_A) p_B[i] ~ dbeta(alpha_B,beta_B) delta[i] <- p_A[i]-p_B[i] } alpha_A ~ XXXXXX alpha_B ~ XXXXXX beta_A ~ XXXXXX beta_B ~ XXXXXX }

Data (Try.r)

"N" <- 60 "c_A" <- structure(c(0,6,0,3,0,8,0,4,0,6,1,5,0,5,0,7,0,3,0,7,0,4,0,5,0,4,0,5,0,4,0,2,0,4,0,5,0,8,2,7,0,6,0,3,0,3,0,8,0,4,0,4,2,6,0,7,0,3,0,1)) "n_A" <- structure(c(0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3)) "c_B" <- structure(c(5,0,2,2,2,0,2,0,2,0,0,0,5,0,4,0,3,1,2,0,2,0,2,0,0,0,3,0,6,0,4,1,5,0,2,0,6,0,1,0,2,0,4,0,4,1,1,0,3,0,5,0,0,0,5,0,2,0,7,1)) "n_B" <- structure(c(9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3,9,0,9,0,9,0,3,0,9,0,9,0,9,3))

These sample data generated using p_A~Beta(30,20), p_B~Beta(20,30), 0 correlation between p_A and p_B, and n=[0,9,0,9,0,9,3] or [9,0,9,0,9,0,3] repeating.

Control (Try.cmd)

model in Try.bug data in Try.r compile, nchains(2) initialize update 400 monitor set p_A, thin(3) monitor set p_B, thin(3) monitor set delta, thin(3) update 1000 coda *, stem(Try)

Originally my model used dgamma(1,20) for the priors on alpha and beta (I've marked this as XXXXXX), but I was always getting `Error in node beta_B, Slicer stuck at value with infinite density`

or `Error in node beta_A, Slicer stuck at value with infinite density`

. Interestingly I never got it for either of the alphas.

I wasn't sure what the infinite densities might be, but I decided to try to remove the peak near zero for alphas/betas by changing the XXXXXX to dgamma(2,20). Unfortunately, I was still getting the same error, but by sheer chance, after running it two dozen times, one time it finished and after inspecting the output, I realized I was seeing a lot of p values near 0 or 1, which I was NOT expecting since the underlying distributions have p_A closely clustered around 0.6 and p_B around 0.4.

To try to force the sampler away from the alpha or beta 0 values corresponding to those 0 or 1 p values, I changed XXXXXX to dpar(0.5,1), and that got the sampler to finish on each try. First output here

Even though it's now completing, I suspect there's something wrong with my model since I had to force the priors away from certain (incorrect) values to get the model to complete.

One hunch is that it has something to do with all the 0 measurements. But since n is 0 for those, c=0=Binom(0 trials, p success rate) should be equally probable given any p, and it shouldn't be affecting the posterior of p (much less alpha or beta)

Anyway, can someone help me understand what's going on?

Last edit: Fabio 2013-09-18