Skip to content

Rejection sampling to draw from distributions

2 messages · Anup Nandialath, Charles C. Berry

#
On Fri, 14 Mar 2008, Anup Nandialath wrote:

            
Usually.

It really helps to "provide commented, minimal, self-contained, 
reproducible code" as you are asked, but you have not shown what dal() is.

If dal() is **vectorized**, then a vectorized approach like this will 
often work well:

rout <- matrix(0,n)
rejected <- seq(n)
k <- n
while( k > 0 ){
   rout[ rejected ] <- rnorm( k )
   ratio <- dal(val=rout[ rejected ], p=p)/(s*dnorm(rout[ rejected ]))
   alpha <- runif( k , min=0, max=1 )
   rejected <- rejected[ alpha >= ratio ]
   k <- length( rejected )
}

Obviously, I cannot test this without dal().

HTH,

Chuck

p.s. Putting a limit on iterations is usually a good idea, too.
Charles C. Berry                            (858) 534-2098
                                             Dept of Family/Preventive Medicine
E mailto:cberry at tajo.ucsd.edu	            UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 92093-0901