/* A/B testing in Picat. From http://rpubs.com/rasmusab/exercise_2_bayesian_ab_testing """ Exercise 2: Bayesian A/B testing for Swedish Fish Incorporated with Stan Rasmus Bååth Swedish Fish Incorporated is the largest Swedish company delivering fish by mail order, but you probably already knew that. The marketing department have done a pilot study and tried two different marketing methods: A: Sending a mail with a colorful brochure that invites people to sign up for a one year salmon subscription. B: Sending a colorful brochure that invites people to sign up for a one year salmon subscription nd that includes a free salmon. The marketing department sent out 16 mails of type A and 16 mails of type B. Six Danes that received a mail of type A signed up for one year of salmon, and ten Danes that received a mail of type B signed up! The marketing department now wants to know, which method should we use, A or B? At the bottom of this document you’ll find a solution. But try yourself first! Question I: Build a Bayesian model in Stan that answers the question: What is the probability that method B is better than method A? """ This is a port of my Racket/Gamble model gamble_ab_testing.rkt This program was created by Hakan Kjellerstrand, hakank@gmail.com See also my Picat page: http://www.hakank.org/picat/ */ import ppl_distributions,ppl_utils. main => go. /* var : rateA Probabilities (truncated): 0.679781617447632: 0.0027322404371585 0.66409340950851: 0.0027322404371585 0.651084620343095: 0.0027322404371585 0.648681496572067: 0.0027322404371585 ......... 0.155947619190415: 0.0027322404371585 0.15560814512689: 0.0027322404371585 0.13884041278569: 0.0027322404371585 0.123382013348575: 0.0027322404371585 mean = 0.389533 Histogram: 0.123: 1 ### (0.003 / 0.003) 0.137: 1 ### (0.003 / 0.005) 0.151: 2 ##### (0.005 / 0.011) 0.165: 3 ######## (0.008 / 0.019) 0.179: 2 ##### (0.005 / 0.025) 0.193: 5 ############# (0.014 / 0.038) 0.207: 5 ############# (0.014 / 0.052) 0.221: 11 ############################# (0.030 / 0.082) 0.235: 3 ######## (0.008 / 0.090) 0.249: 9 ####################### (0.025 / 0.115) 0.262: 11 ############################# (0.030 / 0.145) 0.276: 14 ##################################### (0.038 / 0.183) 0.290: 9 ####################### (0.025 / 0.208) 0.304: 12 ############################### (0.033 / 0.240) 0.318: 16 ########################################## (0.044 / 0.284) 0.332: 21 ####################################################### (0.057 / 0.342) 0.346: 11 ############################# (0.030 / 0.372) 0.360: 19 ################################################## (0.052 / 0.423) 0.374: 20 #################################################### (0.055 / 0.478) 0.388: 19 ################################################## (0.052 / 0.530) 0.402: 23 ############################################################ (0.063 / 0.593) 0.415: 19 ################################################## (0.052 / 0.645) 0.429: 16 ########################################## (0.044 / 0.689) 0.443: 11 ############################# (0.030 / 0.719) 0.457: 11 ############################# (0.030 / 0.749) 0.471: 12 ############################### (0.033 / 0.781) 0.485: 7 ################## (0.019 / 0.801) 0.499: 12 ############################### (0.033 / 0.833) 0.513: 9 ####################### (0.025 / 0.858) 0.527: 9 ####################### (0.025 / 0.883) 0.541: 8 ##################### (0.022 / 0.904) 0.555: 6 ################ (0.016 / 0.921) 0.569: 8 ##################### (0.022 / 0.943) 0.582: 9 ####################### (0.025 / 0.967) 0.596: 4 ########## (0.011 / 0.978) 0.610: 1 ### (0.003 / 0.981) 0.624: 1 ### (0.003 / 0.984) 0.638: 2 ##### (0.005 / 0.989) 0.652: 2 ##### (0.005 / 0.995) 0.666: 2 ##### (0.005 / 1.000) HPD intervals: HPD interval (0.5): 0.29859935413050..0.44048093559243 HPD interval (0.84): 0.22551471750509..0.54798259937576 HPD interval (0.9): 0.21775911292888..0.57901682219422 HPD interval (0.95): 0.18677282388637..0.59081055111755 HPD interval (0.99): 0.15560814512689..0.66409340950851 HPD interval (0.999): 0.12338201334857..0.67978161744763 HPD interval (0.9999): 0.12338201334857..0.67978161744763 HPD interval (0.99999): 0.12338201334857..0.67978161744763 var : rateB Probabilities (truncated): 0.938697629114006: 0.0027322404371585 0.921840754301213: 0.0027322404371585 0.828973054806224: 0.0027322404371585 0.815111212346289: 0.0027322404371585 ......... 0.375526512216556: 0.0027322404371585 0.375144928868462: 0.0027322404371585 0.366619092117352: 0.0027322404371585 0.341878234567995: 0.0027322404371585 mean = 0.622827 Histogram: 0.342: 1 ## (0.003 / 0.003) 0.357: 0 (0.000 / 0.003) 0.372: 3 ####### (0.008 / 0.011) 0.387: 0 (0.000 / 0.011) 0.402: 1 ## (0.003 / 0.014) 0.416: 4 ########## (0.011 / 0.025) 0.431: 5 ############ (0.014 / 0.038) 0.446: 4 ########## (0.011 / 0.049) 0.461: 6 ############## (0.016 / 0.066) 0.476: 4 ########## (0.011 / 0.077) 0.491: 13 ############################### (0.036 / 0.112) 0.506: 5 ############ (0.014 / 0.126) 0.521: 18 ########################################### (0.049 / 0.175) 0.536: 17 ######################################### (0.046 / 0.221) 0.551: 16 ###################################### (0.044 / 0.265) 0.566: 15 #################################### (0.041 / 0.306) 0.581: 18 ########################################### (0.049 / 0.355) 0.596: 15 #################################### (0.041 / 0.396) 0.610: 22 ##################################################### (0.060 / 0.456) 0.625: 25 ############################################################ (0.068 / 0.525) 0.640: 24 ########################################################## (0.066 / 0.590) 0.655: 17 ######################################### (0.046 / 0.637) 0.670: 25 ############################################################ (0.068 / 0.705) 0.685: 25 ############################################################ (0.068 / 0.773) 0.700: 14 ################################## (0.038 / 0.811) 0.715: 15 #################################### (0.041 / 0.852) 0.730: 10 ######################## (0.027 / 0.880) 0.745: 9 ###################### (0.025 / 0.904) 0.760: 9 ###################### (0.025 / 0.929) 0.775: 8 ################### (0.022 / 0.951) 0.789: 11 ########################## (0.030 / 0.981) 0.804: 2 ##### (0.005 / 0.986) 0.819: 2 ##### (0.005 / 0.992) 0.834: 1 ## (0.003 / 0.995) 0.849: 0 (0.000 / 0.995) 0.864: 0 (0.000 / 0.995) 0.879: 0 (0.000 / 0.995) 0.894: 0 (0.000 / 0.995) 0.909: 0 (0.000 / 0.995) 0.924: 2 ##### (0.005 / 1.000) HPD intervals: HPD interval (0.5): 0.57928085447256..0.70810291762841 HPD interval (0.84): 0.49163743364235..0.76199788635690 HPD interval (0.9): 0.48885577055107..0.79692412530860 HPD interval (0.95): 0.43854455251179..0.79943563081298 HPD interval (0.99): 0.36661909211735..0.82897305480622 HPD interval (0.999): 0.34187823456800..0.93869762911401 HPD interval (0.9999): 0.34187823456800..0.93869762911401 HPD interval (0.99999): 0.34187823456800..0.93869762911401 var : rate_diff Probabilities (truncated): 0.614278909105006: 0.0027322404371585 0.551775147464022: 0.0027322404371585 0.53288003594283: 0.0027322404371585 0.529800776638929: 0.0027322404371585 ......... -0.11038087360113: 0.0027322404371585 -0.111878725752178: 0.0027322404371585 -0.118026302716707: 0.0027322404371585 -0.225160235178266: 0.0027322404371585 mean = 0.233294 Histogram: -0.225: 1 ### (0.003 / 0.003) -0.204: 0 (0.000 / 0.003) -0.183: 0 (0.000 / 0.003) -0.162: 0 (0.000 / 0.003) -0.141: 0 (0.000 / 0.003) -0.120: 3 ######## (0.008 / 0.011) -0.099: 3 ######## (0.008 / 0.019) -0.078: 0 (0.000 / 0.019) -0.057: 3 ######## (0.008 / 0.027) -0.036: 8 ##################### (0.022 / 0.049) -0.015: 3 ######## (0.008 / 0.057) 0.006: 9 ####################### (0.025 / 0.082) 0.027: 6 ################ (0.016 / 0.098) 0.048: 8 ##################### (0.022 / 0.120) 0.069: 6 ################ (0.016 / 0.137) 0.090: 17 ############################################ (0.046 / 0.183) 0.111: 18 ############################################### (0.049 / 0.232) 0.132: 13 ################################## (0.036 / 0.268) 0.153: 14 ##################################### (0.038 / 0.306) 0.174: 20 #################################################### (0.055 / 0.361) 0.195: 17 ############################################ (0.046 / 0.407) 0.216: 22 ######################################################### (0.060 / 0.467) 0.237: 14 ##################################### (0.038 / 0.505) 0.258: 18 ############################################### (0.049 / 0.555) 0.279: 22 ######################################################### (0.060 / 0.615) 0.299: 18 ############################################### (0.049 / 0.664) 0.320: 23 ############################################################ (0.063 / 0.727) 0.341: 21 ####################################################### (0.057 / 0.784) 0.362: 21 ####################################################### (0.057 / 0.842) 0.383: 15 ####################################### (0.041 / 0.883) 0.404: 12 ############################### (0.033 / 0.915) 0.425: 5 ############# (0.014 / 0.929) 0.446: 5 ############# (0.014 / 0.943) 0.467: 7 ################## (0.019 / 0.962) 0.488: 7 ################## (0.019 / 0.981) 0.509: 3 ######## (0.008 / 0.989) 0.530: 2 ##### (0.005 / 0.995) 0.551: 1 ### (0.003 / 0.997) 0.572: 0 (0.000 / 0.997) 0.593: 1 ### (0.003 / 1.000) HPD intervals: HPD interval (0.5): 0.16615962570820..0.35740655537574 HPD interval (0.84): 0.01454169071025..0.41182227638169 HPD interval (0.9): 0.00435135048085..0.46469258678364 HPD interval (0.95): -0.04447372259780..0.49398244102205 HPD interval (0.99): -0.11802630271671..0.53288003594283 HPD interval (0.999): -0.22516023517827..0.61427890910501 HPD interval (0.9999): -0.22516023517827..0.61427890910501 HPD interval (0.99999): -0.22516023517827..0.61427890910501 var : diff > 0 Probabilities: 1: 0.9398907103825137 0: 0.0601092896174863 mean = 0.939891 Histogram: 0: 22 #### (0.060 / 0.060) 1: 344 ############################################################ (0.940 / 1.000) HPD intervals: HPD interval (0.5): 1.00000000000000..1.00000000000000 HPD interval (0.84): 1.00000000000000..1.00000000000000 HPD interval (0.9): 1.00000000000000..1.00000000000000 HPD interval (0.95): 0.00000000000000..1.00000000000000 HPD interval (0.99): 0.00000000000000..1.00000000000000 HPD interval (0.999): 0.00000000000000..1.00000000000000 HPD interval (0.9999): 0.00000000000000..1.00000000000000 HPD interval (0.99999): 0.00000000000000..1.00000000000000 */ go ?=> reset_store, run_model(100_000,$model,[show_probs_trunc,mean,show_histogram,show_hpd_intervals, presentation=["rateA","rateB","rate_diff","diff > 0"]]), fail, nl. go => true. model() => NA = 16, % Number of sent mail NB = 16, ObsSA = 6, % Number of signments (observed) ObsSB = 10, RateA = uniform_dist(0,1), % Priors RateB = uniform_dist(0,1), SA = binomial_dist(NA,RateA), % Likelihood SB = binomial_dist(NB,RateB), RateDiff = RateB - RateA, DiffPos = cond(RateB > RateA,1,0), observe(SA == ObsSA, SB == ObsSB), if observed_ok then add("rateA",RateA), add("rateB",RateB), add("rate_diff",RateDiff), add("diff > 0",DiffPos) end.