exp.gibbs: Gibbs Sampling for Two Truncated Exponential Variables,...

In SMPracticals: Practicals for Use with Davison (2003) Statistical Models

Description

Performs Gibbs sampling for problem with two truncated exponential variables. See Practical 11.3 of Davison (2003) for details.

Usage

exp.gibbs(u1 = NULL, u2 = NULL, B, I = 100, S = 100)

Arguments

u1	Initial values for variable 1
u2	Initial values for variable 2
B	Value at which exponential distribution is truncated
I	Number of iterations of sampler
S	Number of replicates of sampler

Details

This is provided simply so that readers spend less time typing. It is not intended to be robust and general code.

Value

A 2 x S x I array containing the values of the variables for the successive iterations

Author(s)

Anthony Davison ([e-mail])

References

Davison, A. C. (2003) Statistical Models. Cambridge University Press.Practical 11.3.

Examples

add.exp.lines <- function( exp.out, i, B=10)
{
  dexp.trunc <- function( u, lambda, B ) 
     dexp(u, rate=lambda)/(1-exp(-lambda*B))
  S <- dim(exp.out)[2]
  I <- dim(exp.out)[3]
  u <- seq(0.0001,B,length=1000)
  fu <- rep(0,1000)
  for (s in 1:S) fu <- fu + dexp.trunc(u,exp.out[3-i,s,I],B)/S
  lines(u,fu,col="red")
  invisible()
}
par(mfrow=c(3,2))
B <-10; I <- 15; S <- 500
exp.out <- exp.gibbs(B=B,I=I,S=S)
hist(exp.out[1,,I],prob=TRUE,nclass=15,xlab="u1",ylab="PDF",xlim=c(0,B),ylim=c(0,1))
add.exp.lines(exp.out,1)
hist(exp.out[2,,I],prob=TRUE,nclass=15,xlab="u2",ylab="PDF",xlim=c(0,B),ylim=c(0,1))
add.exp.lines(exp.out,2)

Example output

Loading required package: ellipse

Attaching package: 'ellipse'

The following object is masked from 'package:graphics':

    pairs

SMPracticals documentation built on May 2, 2019, 11:12 a.m.