Dr. Quan-Hoang Vuong (v19-13) ©2018 SSHPA ∫∫ Vuong & Associates 1 R Codes for Preliminary Checks (Bayesian Method) SSHPA / AISDL 2018 ## column names (variables): scid, bmattitude, bmuse, freqstillgood, capab, bmdiff, ## theodiff, moddiff, resdiff, codediff, plan, trysoft, learntime ## the following codes run well with the hypothetical data set btestdata.csv ## bmdat <read.table("c://dr.vuong/BMSurvey/btestdata.csv", sep=",", header=T) m1 <map( alist( bmuse ~ dbinom(1,p) , logit(p) <a + bfreq*freqstillgood, a ~ dnorm(0,10), bfreq ~ dnorm(0,10) ), data=bmdat) m3 <map( alist( bmuse ~ dbinom(1,p) , logit(p) <a + bfreq*freqstillgood + bcap*capab, a ~ dnorm(0,10), bfreq ~ dnorm(0,10), bcap ~ dnorm(0,10) ), data=bmdat) ## Note for m1: a model with the suspicion that when one considers the old way of ## analysis is still good, one is hesitant to adop the Bayesian method m2 <map( alist( bmattitude ~ dbinom(1,p) , logit(p) <a + bfreq*freqstillgood + bcap*capab, a ~ dnorm(0,10), bfreq ~ dnorm(0,10), bcap ~ dnorm(0,10) ), data=bmdat) ## Note for m2: a model with the suspicion that when one considers the old way ## of analysis is inadequate and one's capability is good enough, ## one will tend to learn and use Bayesian method > compare(m1,m2,m3) WAIC pWAIC dWAIC weight SE dSE m3 28.6 8.6 0.0 0.90 13.79 NA m1 33.1 2.4 4.5 0.09 9.99 10.31 m2 42.3 9.2 13.7 0.00 14.75 8.33 > precis(m3) Dr. Quan-Hoang Vuong (v19-13) ©2018 SSHPA ∫∫ Vuong & Associates 2 Mean StdDev 5.5% 94.5% a -1.71 6.60 -12.26 8.84 bfreq -2.01 6.59 -12.54 8.53 bcap 7.43 3.58 1.71 13.16 > logistic(-1.71) [1] 0.1531637 > logistic(-1.71+7.43) [1] 0.996731 > compare(m3b,m3b.assoc,m3c) WAIC pWAIC dWAIC weight SE dSE m3b 15.3 1.9 0.0 0.61 7.86 NA m3b.assoc 16.2 2.3 0.9 0.39 7.80 0.45 m3c 30.5 9.4 15.1 0.00 14.44 10.58 > plot(compare(m3b,m3b.assoc,m3c)) > precis(m3b.assoc) Mean StdDev 5.5% 94.5% a -3.50 0.85 -4.85 -2.14 bfreq 1.85 3.34 -3.49 7.19 bcap 5.32 2.42 1.46 9.18 # HMC Stan m3b.assoc.stan <map2stan(m3b.assoc, data=bmdat, iter=10000, warmup=1000) # check Gaussian pairs(m3b.assoc.stan) Dr. Quan-Hoang Vuong (v19-13) ©2018 SSHPA ∫∫ Vuong & Associates 3 > precis(m3b.assoc.stan) Mean StdDev lower 0.89 upper 0.89 n_eff Rhat a -3.78 0.94 -5.24 -2.31 4728 1 bfreq 3.02 3.89 -2.91 9.54 4114 1 bcap 7.11 2.98 2.36 11.77 3652 1 > post <extract.samples(m3b.assoc.stan) > str(post) List of 3 $ a : num [1:9000(1d)] -4.42 -4.6 -2.95 -3.94 -2.11 ... $ bfreq: num [1:9000(1d)] 0.298 4.228 4.957 5.496 3.666 ... $ bcap : num [1:9000(1d)] 6.41 5.34 3.91 3.98 3.93 ... > dens(post$a) Dr. Quan-Hoang Vuong (v19-13) ©2018 SSHPA ∫∫ Vuong & Associates 4 > dens(post$bcap) # MCMC Stan for m3b model > m3b.stan <map2stan(m3b, data=bmdat, iter=10000, warmup=1000) > precis(m3b.stan) Mean StdDev lower 0.89 upper 0.89 n_eff Rhat a -1.33 3.21 -6.24 4.01 2489 1 bfreq -2.55 3.20 -7.88 2.29 2480 1 bcap 7.78 2.61 3.55 11.61 2511 1 > pairs(m3b.stan) > post2 <extract.samples(m3b.stan) > dens(post2$a) Dr. Quan-Hoang Vuong (v19-13) ©2018 SSHPA ∫∫ Vuong & Associates 5 > dens(post$bfreq) Stan MCMC performs quite well. References McElreath R. (2018). Statistical Rethinking: A Bayesian Course with Examples in R and Stan. New York: Chapman and Hall/CRC. Vuong QH. (2018). The (ir)rational consideration of the cost of science in transition economies. Nature Human Behaviour, 2(1), 5. Vuong QH, La VP, Vuong TT,... & Ho MT. (2018). Cultural additivity: Behavioural insights from the interaction of Confucianism, Buddhism, and Taoism in folktales. Palgrave Communications, 4(1), 143.