Grapham: Graphical Models with Adaptive Random Walk Metropolis Algorithms

GRAPHA M: GRAPHIC AL MODELS WITH AD APTIVE RANDOM W ALK METROPOLIS ALGORITH MS MA TTI VIHOLA Abstract. Recen tly developed adaptive Mark ov chain Mon te Car lo (MCMC) metho ds hav e been applied success fully to man y problems in Bay esian statistics. Grapham is a new open source implemen tation covering several such methods, with emphasis o n g raphical mo dels f or directed acyclic gr aphs. The implemen ted algo- rithms include the seminal Adaptive Metr o p olis algo r ithm a djusting the prop osa l cov ar iance according to the histor y of the chain and a Metrop o lis alg orithm adjust- ing the prop o s al scale bas ed on the obser ved acceptance probability . Different v ar i- ants of the algor ithms allow one, for example, to use these tw o algor ithms together, employ delay ed rejection and adjust several parameters of the algorithms. T he im- plement ed Metrop olis-within-Gibbs update allows arbitrary sampling blo cks. The softw are is written in C and uses a simple extensio n langua ge Lua in configur a tion. 1. Intr oduction Mark o v chain Monte Carlo (MCMC) is a g eneral framew ork for computing ex- p ectations o v er complicated distributions in general state s paces. The metho ds are based on constructing a Mark ov chain ( X n ) n ≥ 1 so that the ergo dic av erages I N = N − 1 P N k =1 f ( X k ) con v erge to I = R f ( x ) π ( x )d x as N → ∞ , where π is the target distribution of interes t. Suc h a c hain is often easy to construct using the Metrop olis-Hastings a lgorithm; see, for example, Rob ert a nd Casella (199 9). De- p ending on π , ho w ev er, it ma y b e difficult to design a practical algorithm so that I N w ould appro ximate I w ell with a mo derate n um b er of samples N . Recen tly prop osed adaptiv e MCMC algorit hms adjust the pa r ameters of the algo- rithm (the propo sal distribution) on-the-fly , aiming to allo w efficien t sim ulation. They ha v e attracted increasing attention in the last few y ears, af t er Haario et al. (2001) presen ted the seminal Adaptiv e Metrop olis (AM) algorit hm, and Andrieu and Rob ert (2001) related ada ptive MCMC to the g eneral con text of the Robbins-Monro sto chas - tic appro ximation. After that, sev eral authors hav e pro p osed new algorithms and v a riations, and pro vided theoretical v alidation of the metho ds (Haario et al., 2005, 2006; A tc had ´ e and Rosen thal, 2005; Andrieu and Moulines, 2006; Roberts a nd Ro senthal, 2009, 2007; Sa ksman and Vihola, 2008; Atc had ´ e and F ort, 2008; Bai et al., 2008; Viho la , 2009); see also the recen t review b y Andrieu and Thoms (2008) and references therein. Date : September 2, 2009. Key wor ds and phr ases. Adaptive Marko v c hain Mon te Carlo , Bay esian statistics, computational metho d, graphical mode l, s ta tistical soft ware. 1 2 MA TTI V IHOLA Grapham is a n op en source implemen tation of sev eral adaptive MCMC algorithms based on the random walk Metrop olis sampler. The purp ose o f Gra pha m is to pro vide an exp erimen tal to ol for ev aluating the p erfo rmance of suc h a lgorithms with practical problems, esp ecially in Bay esian statistics. The source co de of the soft w are and addi- tional documen tatio n a r e a v ailable for downloading in http://iki.fi/mvi hola/grapham/ . Rosen thal (2007) describ es another adaptiv e MCMC implemen tatio n: AMCMC. It is an R interface to one adaptiv e MCMC algorithm (referred to as ‘ASCM’ in Section 2 b elo w). Grapham differs fr o m AMCMC in that it relies on a hierarc hical mo del sp ecification and provides more alternativ e algorithms. Unlik e AMCMC, G rapham also prov ides a set of ready-made standard distribution functions the user can emplo y as a part of their mo del sp ecification. This is in tended to allow faster dev elopmen t while p ermitting the user to define arbitrary distributions easily . The mo dels are sp ecified in Gra pham by defining a set of v ariables with their conditional distributions. Suc h mo dels are often referred to as ‘graphical mo dels’; see, for example, Lauritzen (1996) and references therein. This underlying philosophy of Grapham reminds that of BUGS (Spiegelhalter et al., 1996–200 8); see also the review Murph y (200 7) of other softw are for graphical mo dels. The adv an tage of Grapham o ve r BUGS is that the adaptiv e MCMC algorithms can b e muc h more efficien t than the non-adaptive (Metrop olis-within-)Gibbs algorithms of BUGS. One should, ho w ev er, notice that Grapham is an exp erimen tal pro ject not offering the v ersatilit y and maturity of BUGS. It is also lik ely that BUGS p erforms b etter than Grapham with man y simpler mo dels. 2. Algorithms The general form of the algorithms implemen ted in Gr apham can b e described as follo ws. Let X 0 ≡ x 0 ∈ R d b e a given starting p oint f or the state c hain, and θ 0 and L 0 stand for the initial scaling parameter and the (low er-diagonal with non-zero diagonal) shap e matrix, respective ly . F or n ≥ 1, the recursion fo llo ws: (S1) form a prop o sal Y n = X n − 1 + θ n − 1 L n − 1 W n , where W n is an indep enden t sample from a symmetric prop osal distribution q 0 , (S2) with probability α n = min { 1 , π ( Y n ) /π ( X n − 1 ) } , the prop osal is accepted and X n = Y n ; otherwise, the prop osal is rejected and X n = X n − 1 , and (S3) up date the scaling parameter θ n − 1 → θ n > 0 and the shap e L n − 1 → L n ∈ R d × d according to the selected adaptiv e algorithm. The steps ( S1 ) and (S2) implemen t an iteration o f the random-w alk Metrop olis al- gorithm with the prop osal distribution q 0 scaled b y the f a ctor θ n − 1 L n − 1 . Step (S3) implemen ts the adaptation, c hanging the scaling parameters θ n and L n based on t he history of the c hain. Examples of suc h up dates are give n b elow. Instead of applying the iteration (S1)– ( S3) at once to all the elemen ts of the v ec- tor X n , o ne ma y use Metrop o lis-within-Gibbs a nd apply t he iteratio n sequen tially to subsets of the elemen ts of X n , as in the single comp onen t AM algorithm suggested b y Haario et al. (2005) . These sampling blo c ks can b e selected freely in Gra pham. GRAPHAM: GRAPHICAL MODELS WITH ADAPTIVE METR OPOLIS 3 The pro p osal distribution q 0 in (S1) can a lso b e ch osen. Grapham curren tly imple- men ts (m ultiv ariate) Gaussian, studen t, uniform (in a cub e) and (a d - fold pro duct of ) Lapla ce prop osal distributions. The a daptation of (S3) depends o n the selec ted algorithm. The Adaptiv e Metrop- olis (AM) algorithm of Haario et a l. (2001 ) implies constan t scaling θ n = θ 0 for all n ≥ 1. The shap e matrix L n is the Cholesky factor of a cov ariance estimate of the c hain. In part icular, L n L T n = C n with a p ositiv e definite C 0 ∈ R d × d and defined through M n = (1 − η n ) M n − 1 + η n X n and (1) C n = (1 − η n ) C n − 1 + η n ( X n − M n − 1 )( X n − M n − 1 ) T , (2) with M 0 ≡ x 0 . The we ight sequence η n ∈ (0 , 1) can b e selected arbitrarily , but it is recommended to c ho ose η n deca ying to zero. F or example, setting η n = η 0 > 0 for all n ≥ 1 results in an algo rithm similar to the Adaptive Prop osal (AP) a lg orithm (Haario et al., 1999). This alg orithm do es no t , in general, prov ide v alid sim ulatio n; see the example in Haario et al. (2001). The original AM alg orithm emplo ys the default v alue η n = ( n + 1) − 1 , in whic h case M n and C n coincide with the av erage and (a symptotically) the sample co v a riance of X 0 , . . . , X n , resp ective ly . The up dated Cholesky factor L n +1 of C n +1 is computed efficien tly f rom L n b y a rank one up date requiring O ( d 2 ) o p erations (Donga rra et al., 1979). Observ e that the same order of op erations is needed when forming the prop o sal Y n in (S1). The a daptiv e scaling Metrop olis (ASCM) algorithm as prop osed b y A tc had ´ e and Rosenthal (2005) and Rob erts and Rosen thal (2007, 2 009) lea v es t he shap e matrix constan t L n = L 0 for all n ≥ 1. The scaling parameter θ n is up da ted according to the ob- serv ed acceptance probability . The default up date in Grapham is (3) θ n = θ n − 1 h 1 + η n  α n α ∗ − 1 i , where α ∗ is the desired acceptance probabilit y . The default v alues for α ∗ are 0 . 44 in dimension one a nd 0 . 234 otherwise fo llo wing Ro b erts and Rosen thal (20 09). The user can also supply an alternativ e, arbitrar y up date function easily , a s exemplified in Section 4. These t w o algorithms, AM and ASCM, can b e used sim ultaneously , as suggested in A tc had´ e and F ort (200 8) and Andrieu and Thoms (2008). Additional fla v ours to the algorithms include a Rao -Blac kw ellised v ersion of AM (Andrieu and Thoms, 2008) mo difying the up date form ulae (1) and (2) to M n = (1 − η n ) M n − 1 + η n [ α n Y n + (1 − α n ) X n − 1 ] and (4) C n = (1 − η n ) C n − 1 + η n  α n ( Y n − M n − 1 )( Y n − M n − 1 ) T (5) +(1 − α n )( X n − 1 − M n − 1 )( X n − 1 − M n − 1 ) T  . There is a p ossibilit y to use (t w o-stage) dela y ed rejection (DR) with AM (Haario et al., 2006). DR can also b e applied when using ASCM, so that only the first-stage a ccep- tance probability α n is emplo ye d in (3). 4 MA TTI V IHOLA Grapham implemen ts three different burn-in strategies f or adaptat ion. The default ‘greedy’ strategy p erforms con tinuous adaptatio n during the whole MCMC run. The ‘traditional’ strategy as prop osed in Haario et a l. (2001) uses a fixed prop osal fo r the burn-in and then p erforms contin uous adaptation during the rest of the sim ulation. One may also apply a ‘freeze’ strategy adapting only during the burn- in a nd k eeping the obtained parameters fixed during the estimation run. It is po ssible to employ a mixture of tw o prop osal density comp onents , a fixed and an adaptiv e one (Rob erts and Rosenthal, 2 009). This is implemen ted in G rapham so that, with probability p mix , the initial pa rameters L 0 and θ 0 are used in (S1) instead of the adapted v a lues θ n − 1 and L n − 1 . The user may define also a non-constant mixing probabilit y p ( n ) mix ∈ [0 , 1]. This feature can b e used, for example, to in tro duce a ‘gradual burn-in,’ by defining a deca ying sequence p ( n ) mix → 0 . 3. Implement a tion Grapham does not ha v e an in teractiv e ‘use r in terface.’ It is simply executed from the command prompt (shell) with input file names a s parameters. The in- put files con tain the mo del sp ecification a nd the sim ulation parameters. It is also p ossible to define the functional of in terest in the input files. F o r more compli- cated functionals, ho w ev er, it may b e con v enien t to store (a subset o f ) the samples sim ulated b y Gr a pham and pro cess them in ano t her env ironment. The samples can b e sa ve d into a file in the CSV (comma separated v alues) f ormat or in a sim- ple binary format. The former allow s the samples to b e easily imp o rted to man y other en vironmen ts. There are r eady-made functions for loading the binary data files in to R (R Deve lopmen t Core T eam, 2 0 09), Matlab R  (The MathW orks, Natick , Massac h usetts) and Octav e (Eato n, 2002) env ironments . The core of Grapha m is implemen ted in C. It includes some n umerical F ortran subroutines from the Netlib rep o sitory (Browne et al., 1995) a nd can optionally b e compiled with the dSFMT ra ndo m n um b er generator of Saito and Matsumoto (200 8) instead of using the ra ndo m num ber generators provide d b y the C standard libra ries. The configuration of Grapham is done using the small and public ly av ailable extension language Lua (Ierusalimsc h y et al., 199 6). While minimalistic, Lua is in fact a full- featured programming language offering a great flexibilit y . F or example, the user can supply functions as configurat io n parameters and apply data from external files in the mo del. In fa ct, Grapham includes some functions written in Lua, for example for reading data files in the CSV format. The Numeric Lua pac k age (Carv alho, 2005) can also b e compiled with Grapham to allow easy w orking with ve ctor- v alued v ariables. There ar e n umerous ready-made distribution functions a v ailable fo r defining the conditional densities asso ciated with the v ariables. The densities can also b e defined arbitrarily as L ua functions. Like wise, the functional of interes t may b e written in Lua. Ho we v er, to allow optimal p erformance, G rapham allo ws the user to supply densities and functionals in a separate C library with ease. GRAPHAM: GRAPHICAL MODELS WITH ADAPTIVE METR OPOLIS 5 µ a t 1 t 2 t 18 · · · y 1 y 2 y 18 Figure 1. The graphical represen tation of the baseball mo del. The no des with observ ed v alues (‘data’) are sho wn in grey . 4. An Example Se ssion Consider the baseball mo del of Rosenthal (1996) used as an example also with AM- CMC (Rosen thal , 2007). It consists of 38 real-v alued v ar ia bles, defined hierar chically as depicted in Fig. 1 . The file sp ecifying this mo del in Gr a pham is shown in Fig. 2 . The mo del is defined in the Lua table model , defined in lines 4–20. Each v ariable is defined by a n en try containing a loga rithmic densit y , conditional on the paren t v a riables. The v ariables µ , t and y in the example ha v e standard distributions: µ has (an improp er) uniform distribution ov er R , while t and y a re conditionally Gaussian with means µ and t a nd v ariances a and v , resp ectiv ely . The recipro cal of the v ariable a is exp onen tially distributed; this is defined through a Lua function defined in lines 16–18, calling dexp , the exp onen tial distribution f unction. The mo del is, in fact, then mo dified by t he function repeat_bloc k . The blo ck o f v a riables ( y , t ) in the mo del is replicated 18 times to obtain blo c ks ( y 1 , t 1 ) , . . . , ( y 18 , t 18 ). A t the same time, t he function repeat_bl ock sets the v alues of y i to t he 18 v alues read fro m the CSV file baseball.da ta using the function read_csv . The following sho ws an example run o f G rapham with the mo del sp ecification o f Fig. 2. $ ./grapham models/baseba ll.lua Functional average = [ 0.392507 0.267393 0.318917 ] Acceptance rates: ( a ): 44.03% ( t7 ): 43.97% ( t9 ): 44.00% The pa r t of the output sho wn ab o ve con tains the computed estimate of the exp ected v a lue of the f unctiona l sp ecified in lines 23–25 of Fig . 2, that is, the mean o f the v ector [ t 1 , µ, a ], giving a similar estimate fo r t 1 as obta ined b y AMCMC (Rosen thal , 2007). Moreov er, the a v erage acceptance rate of the no des w as approximately 44%, whic h is the default v alue of the desired acceptance probabilit y α ∗ . The run consisted of 40000 (of whic h 10 0 00 burn-in) iterations using the ASCM algo rithm for each real- v a lued v ariable a t a time. This algor ithm is v ery similar t o the one implemen ted in AMCMC. The running time of Grapham w as a ppro ximately 1.0 seconds with In tel P entium 4 at 2.8 0GHz. As a comparison, the same run with AMCMC (with b oth the densit y and the functional sp ecified in C for o ptimal p erformance) to ok appro ximately 3.8 seconds . The faster sim ulation sp eed of Grapham is explained b y 6 MA TTI V IHOLA 1 const = { 2 v = 0.00434 3 } 4 model = { 5 mu = { 6 density = "duni form" 7 }, 8 t = { 9 parents = {"mu" ,"a"}, den sity = "dnorm" 10 }, 11 y = { 12 parents = {"t", "v"}, density = "dnorm " 13 }, 14 a = { 15 init_val = 1, 16 density = functi on(a_) 17 return dexp(1/a _, 1/2) 18 end 19 }, 20 } 21 _, y = read_csv( "models/ baseball.da ta") 22 repeat_bl ock({"y" ,"t"}, y [1]) 23 function functio nal() 24 return {t1, mu, a} 25 end 26 para = { 27 niter = 30000, nburn = 10000, algorith m = "as cm", 28 } Figure 2. The Lua co de in the file models/basebal l.lua sp ecifying the mo del of Fig. 1 in Grapham. the hierarc hical mo del setup, whic h Grapham can take a dv an tage of. That is, only part of the conditional densitie s in the target distribution need to b e ev alua t ed when eac h v ariable is up dated. Let us mo dify t he a b ov e example, b y adding the lines show n in Fig. 3 to the mo del sp ecification of F ig. 2. The supplied function para.scaling_adapt replaces the default up date in ( 3), and in fact implemen ts exactly the scaling adaptation algorithm of AMCMC ( Rosen thal, 20 0 7; Rob erts and Ro sen thal, 2009). The v a lue set to the parameter para.dr means that dela y ed rejection is used, with a 0.1 times do wn-scaled prop osal in the second stage. Moreo v er, instead of the default Gaussian distribution, the prop o sal samples are draw n from a Studen t’s t -distribution. $ ./grapham models/baseba ll.lua models/amcmc _dr.lua Functional average = [ 0.392465 0.266204 0.321466 ] Acceptance rates: ( a ): 70.26% (47.49%/22.77%) ( t7 ): 70.66% GRAPHAM: GRAPHICAL MODELS WITH ADAPTIVE METR OPOLIS 7 1 para.sca ling_ada pt = funct ion(sc, alpha, dim, k) 2 if alpha >0.44 then 3 delt a = 1 4 else 5 delt a = -1 6 end 7 retu rn sc*exp(delta *min(0.0 1, 1/sqrt(k+ 1))) 8 end 9 para.dr = 0.1; para.pro posal = "s tudent" Figure 3. The Lua co de in the file models/amcmc dr.lua . (47.55%/23. 10%) ( t9 ): 71.08% (47.62%/23.46%) In this case, the total acceptance rate of eac h blo c k is around 70%, of whic h roughly t w o thirds are accepted in the first stage and one third in the second, delay ed rejection stage. The estimates obtained for t 1 , µ and a a pp ear similar to the first run. Finally , to exemplify how the data sim ulated by Grapham can b e used in ot her en vironmen ts, let us run G rapham with the command line $ ./grapham models/baseba ll.lua -e "para.outfile= ’bb.bin’" This command includes the c huc k of Lua co de para.outfile=’bb.bi n’ after r eading the file baseball.lua . As a consequence , the sim ulated samples are written in the file bb.bin . In R , o ne could, fo r example, write > source("tools /grapham_read.r") > data <- grapham_read("b b.bin", nthin=10) > plot(data$a, data$t1) whic h w ould plot ev ery tenth o f the 30000 sim ulated samples of ( a, t 1 ) in t he same figure. 5. Conclusions Grapham provides a flexible op en-source t est b ed for ev aluat ing the p erformance of differen t adaptive random w alk Metrop olis alg o rithms, especially with hierarc hi- cal mo dels o ften encoun tered in Ba ye sian statistics. It provide s a fa irly simple and general wa y of determining mo dels and functionals and for incorp orating data in to the mo del. The simulation sp eed of Grapham is go o d, even in a relatively high- dimensional setting, as the implemen ted algorithms inv olv e at most a quadrat ic n um b er of op erations with resp ect to the dimension. The user has extensiv e con- trol o v er the v arious parameters of the algorithms, enabling a thorough testing of differen t a daptation strategies. Moreo ve r, ne w adapt ive algorithms of the similar random walk t yp e can b e easily added to Gr a pham. 8 MA TTI V IHOLA A cknowl e dgements The author thanks Professor Antti Pen ttinen and the referees f or helpful commen ts on the man uscript. This w ork was supp ort ed by the Academ y o f F inla nd, pro jects no. 1 1 0599 and 201392 , by the Finnish Academ y of Science and Letters, the Vilho, Y rj¨ o and Kalle V¨ ais¨ al¨ a F oundatio n, b y the F innish Cen tre of Excellence in Analysis and Dynamics Researc h and b y the F innish G raduate School in Sto chas tics and Statistics. Reference s C. Andrieu and ´ E. Moulines. On the ergo dicit y prop erties of some adapt ive MCMC algorithms. A nn. Appl. Pr ob ab. , 16(3):1 4 62–1505 , 2 006. C. Andrieu and C. P . Rob ert. Con trolled MCMC for optimal sampling. T ec hnical Rep ort Ceremade 012 5, Univ ersit ´ e P aris Dauphine, 20 0 1. C. Andrieu and J. Thoms. A tutorial on adaptiv e MCMC. Statist. Co mput. , 18(4): 343–373, D ec. 2008 . Y. A tc had´ e and G . F ort . Limit theorems for some adaptiv e MCMC algor it hms with subgeometric k ernels. URL http://arxiv .org/abs/0807.2952 . Preprin t, July 2008. Y. F. A tc had´ e and J. S. Rosen thal. On a daptiv e Mark ov c hain Monte Carlo alg o- rithms. Bernoul li , 11(5 ):815–828 , 2005. Y. Bai, G. O. Rob erts, and J. S. Rosen thal. On the con tainmen t condition f o r adaptiv e Marko v chain Monte Carlo a lgorithms. URL http://prob ability.ca/jeff/research.html . Preprin t, July 2008. S. Brow ne, J. Do ng arra, E. G rosse, and T. Row a n. The Netlib mathematical softw are rep ository . T echnic al rep o rt, Corp oration for National Researc h Initiativ es, 19 9 5. URL http://netlib.org/ . L. Carv a lho . Numeric Lua pro ject page, 2005. URL http://luaf orge.net/projects/numlua/ . Visited 1 1th June 2009. J. J. Dongar r a , J. R . Bunch, C. B. Moler, and G. W. Stew art. LINP ACK Users’ Guide . So ciet y for Industrial and Applied Mathematics, 1979. ISBN 0-89871- 1 72- X. J. W. Eaton. GNU Octave Manual . Net w ork Theory Limited, 2002. ISBN 0 -95416 1 7- 2-6. URL http://www.gnu .org/software/octave/ . H. Haario, E. Saksman, a nd J. T amminen. Adaptiv e prop osal distribution for r a ndom w alk Metrop olis a lgorithm. Co mput. Statist. , 1 4(3):375– 3 95, 1999. H. Haario, E. Saksm an, and J. T amminen. An adaptiv e M etrop olis algorithm. Bernoul li , 7 (2):223–2 4 2, 2001. H. Haario , E. Saksman, and J. T a mminen. Comp onen twis e adaptation for high dimensional MCMC. Comput. Statist. , 20(2 ):265–274 , 2 005. H. Haario, M. Laine, A. Mira, and E. Saksman. DRAM: Efficien t adaptiv e MCMC. Statist. Comp ut. , 16(4):339 –354, Dec. 200 6. GRAPHAM: GRAPHICAL MODELS WITH ADAPTIVE METR OPOLIS 9 R. Ierusalimsc h y , L. H. de Figueiredo, and W. C. F ilho. L ua — an extensible extension language. Softwar e: Pr actic e and Exp erienc e , 26(6):63 5–652, June 1 996. URL http://www. lua.org/ . S. L. Lauritzen. Gr aphic al Mo dels . Oxford Univ ersit y Press, 1996. ISBN 0 -19-85 2219- 3. K. Murph y . Soft w are for graphical mo dels: A review. ISBA Bul letin , 14(4), Dec. 2007. R D ev elopmen t Core T eam. R: A L anguage and Envir onment for Statistic al Com- puting . R F oundation for Statistical Computing, Vienna, Austria, 2009. URL http://www. R- project.org . C. P . Rob ert and G . Casella. Monte Carlo Statistic a l Metho ds . Springer-V erlag, New Y ork, 1999. ISBN 0-387-98 707-X. G. O. Rob erts and J. S. R osen thal. Coupling and ergo dicit y of adaptiv e Marko v c hain Mon te Carlo algorithms. J. Appl. Pr ob a b. , 44(2) :458–475 , 20 07. G. O. Rob erts and J. S. Ro sen thal. Examples o f adaptiv e MCMC. J. Comput. Gr aph. Statist. , 18(2 ):349–367 , 2 009. J. S. Rosen thal. AMCMC: An R interface for adaptiv e MCMC. Comput. Statist. Data A nal. , 51(12) , Aug. 2007. URL http://www.proba bility.ca/amcmc/ . J. S. Rosen thal. Analysis of the G ibbs sampler for a mo del related to James-Stein estimators. Statist. C omput. , 6(3):26 9–275, 1996 . M. Saito and M. Matsumoto. SIMD-orien ted fa st Mersenne t wister: a 128-bit pseu- dorandom num ber generator. In Mon te Carlo and Quasi-Monte Carlo Metho ds 2006 , pages 607–622. Springer, 2008. E. Saksman and M. Vihola. On the ergo dicit y of the adaptiv e Metrop olis algorithm on unbounded do mains. URL http://arxiv. org/abs/0806.2933 . Preprint, June 2008. D. J. Spiegelhalter, A. Thomas, N. G. Best, W. R. Gilks, and D. Lunn. BUGS: Bay esian inference using Gibbs sampling, 1996–200 8. URL http://www. mrc- bsu.cam.ac . uk/bugs/ . M. Vihola. On the stabilit y and ergo dicity of an adaptiv e scaling Metrop olis algo- rithm. URL http://arx iv.org/abs/0903.4061 . Preprint, Mar. 2009. Ma tti Vihola, Dep ar tment of Ma thema tics and St a tistics, P.O. Box 35 (MaD), FI-40014 University of Jyv ¨ askyl ¨ a, Finland E-mail addr ess : matt i.vih ola@i ki.fi URL : http:// iki.f i/mvihola/