Gibbs sampler

One sweep replaces each block in turn by a draw from its full conditional, leaving the joint distribution invariant. Each conditional is supplied as a function of the current state and one uniform variate, so the sampler is an inverse-CDF draw; the uniforms come from the van der Corput sequence, which makes the run reproducible and identical across language arms.

Usage

Gibbsm(conditionals, x0, n_iter = 100, burn = 0)

Arguments

conditionals

List of functions f(x, u) returning the new value of block i.

x0

Starting state, one entry per conditional.

n_iter

Number of sweeps.

burn

Sweeps discarded before the summaries are formed.

Value

List with estimate, mean, draws, last, n, method.

Details

Formula: x_i <- conditional_i(x, u) for each block in turn.

References

Geman and Geman (1984), Stochastic relaxation, Gibbs distributions, and the Bayesian restoration of images, IEEE Transactions on Pattern Analysis and Machine Intelligence 6(6):721-741. doi:10.1109/TPAMI.1984.4767596

Examples

conditionals <- list(function(x, u) qnorm(u, 0.5 * x[2], 1),
                     function(x, u) qnorm(u, 0.5 * x[1], 1))
Gibbsm(conditionals, x0 = c(0, 0), n_iter = 50, burn = 10)
#> $estimate
#> [1] -0.5675082
#> 
#> $mean
#> [1] -0.5675082  0.4883692
#> 
#> $draws
#>              [,1]         [,2]
#>  [1,] -0.67448975  0.337244875
#>  [2,] -0.98172694 -0.172224107
#>  [3,] -0.40475142  0.947973672
#>  [4,] -1.06013371 -0.372756170
#>  [5,] -0.67515450  0.549569311
#>  [6,] -0.61236190  0.182595459
#>  [7,] -0.06601295  1.501114067
#>  [8,] -1.11217483 -0.477675004
#>  [9,] -0.81796966  0.367436929
#> [10,] -0.82627170 -0.010885787
#> [11,] -0.24264500  1.196688396
#> [12,] -0.71966670 -0.122631240
#> [13,] -0.46356569  0.778207327
#> [14,] -0.38731810  0.385473113
#> [15,]  0.11432414  1.919893939
#> [16,] -1.19392772 -0.557787777
#> [17,] -0.90499290  0.272017933
#> [18,] -0.94150660 -0.110623408
#> [19,] -0.33300214  1.063357687
#> [20,] -0.88611829 -0.245960063
#> [21,] -0.56807656  0.662743478
#> [22,] -0.49913914  0.283840137
#> [23,]  0.02415019  1.688014820
#> [24,] -0.83193231 -0.298196282
#> [25,] -0.68250785  0.489256955
#> [26,] -0.70215328  0.094019885
#> [27,] -0.15008914  1.342752567
#> [28,] -0.55848248 -0.001550798
#> [29,] -0.36090529  0.897062922
#> [30,] -0.27598292  0.488107551
#> [31,]  0.20487769  2.256313539
#> [32,] -1.28940225 -0.625116838
#> [33,] -0.96266249  0.217952058
#> [34,] -1.00421825 -0.162797518
#> [35,] -0.37950117  0.999413764
#> [36,] -0.97376070 -0.309710527
#> [37,] -0.62168039  0.605716474
#> [38,] -0.55562624  0.233152688
#> [39,] -0.02093706  1.590540136
#> [40,] -0.96640034 -0.385128019
#> [41,] -0.74868960  0.428827764
#> [42,] -0.76348366  0.041834253
#> [43,] -0.19618982  1.268108906
#> [44,] -0.63864419 -0.061928568
#> [45,] -0.41206974  0.837123394
#> [46,] -0.33165368  0.436622614
#> [47,]  0.15952837  2.067192071
#> [48,] -0.95383185 -0.418132989
#> [49,] -0.81151595  0.344457402
#> [50,] -0.87092956 -0.054359326
#> 
#> $last
#> [1] -0.87092956 -0.05435933
#> 
#> $n
#> [1] 50
#> 
#> $method
#> [1] "componentwise full-conditional updates with van der Corput uniforms, Geman & Geman (1984)"
#> 
#> attr(,"class")
#> [1] "morie_rich_result" "list"