Re: [R] Question about acception rejection sampling - NOT R question

2007-07-16 Thread Greg Snow
Not a strict proof, but think of it this way:

The liklihood of getting a particular value of x has 2 parts.  1st x has
to be generated from h, the liklihood of this happening is h(x), 2nd the
point has to be accepted with conditional probability f(x)/(c*h(x)).  If
we multiply we get h(x) * f(x)/ ( c* h(x) ) and the 2 h(x)'s cancel out
leaving the liklihood of getting x as f(x)/c.  The /c just indicates
that approximately 1-1/c points will be rejected and thrown out and the
final normalized distribution is f(x), which was the goal.

Hope this helps,

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
[EMAIL PROTECTED]
(801) 408-8111
 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Leeds, 
 Mark (IED)
 Sent: Friday, July 13, 2007 2:45 PM
 To: r-help@stat.math.ethz.ch
 Subject: [R] Question about acception rejection sampling - 
 NOT R question
 
 This is not related to R but I was hoping that someone could 
 help me. I am reading the Understanding the Metropolis 
 Hastings Algorithm
 paper from the American Statistician by Chip and Greenberg, 
 1995, Vol 49, No 4. Right at the beginning they explain the 
 algorithm for basic acceptance  rejection sampling in which 
 you want to simulate a density from f(x) but it's not easy 
 and you are able to generate from another density called 
 h(x). The assumption is that there exists some c such that 
 f(x) = c(h(x) for all x
 
 They clearly explain the steps as follows :
 
 1) generate Z from h(x).
 
 2) generate u from a Uniform(0,1)
 
 3) if u is less than or equal to f(Z)/c(h(Z) then return Z as 
 the generated RV; otherwise reject it and try again.
 
 I think that, since f(Z)/c(h(z) is U(0,1), then u has the 
 distrbution as f(Z)/c(h(Z).
  
 But, I don't understand why the generated and accepted Z's 
 have the same density  as f(x) ?
 
 Does someone know where there is a proof of this or if it's 
 reasonably to explain, please feel free to explain it.
 They authors definitely believe it's too trivial because they 
 don't. The reason I ask is because, if I don't understand 
 this then I definitely  won't understand the rest of the 
 paper because it gets much more complicated.  I willing to 
 track down the proof but I don't know where to look. Thanks.
 
 
 This is not an offer (or solicitation of an offer) to 
 buy/se...{{dropped}}
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


[R] Question about acception rejection sampling - NOT R question

2007-07-13 Thread Leeds, Mark \(IED\)
This is not related to R but I was hoping that someone could help me. I
am reading the Understanding the Metropolis Hastings Algorithm
paper from the American Statistician by Chip and Greenberg, 1995, Vol
49, No 4. Right at the beginning they explain the algorithm for basic
acceptance  rejection sampling in which you want to simulate a density
from f(x) but it's not easy and you are able
to generate from another density called h(x). The assumption is that
there exists some c such that f(x) = c(h(x) for all x

They clearly explain the steps as follows :

1) generate Z from h(x).

2) generate u from a Uniform(0,1)

3) if u is less than or equal to f(Z)/c(h(Z) then return Z as the
generated RV; otherwise reject it and try again.

I think that, since f(Z)/c(h(z) is U(0,1), then u has the distrbution as
f(Z)/c(h(Z).
 
But, I don't understand why the generated and accepted Z's have the same
density  as f(x) ?

Does someone know where there is a proof of this or if it's reasonably
to explain, please feel free to explain it.
They authors definitely believe it's too trivial because they don't. The
reason I ask is because, if I don't understand this then 
I definitely  won't understand the rest of the paper because it gets
much more complicated.  I willing to track down the proof but I don't
know where to look. Thanks.


This is not an offer (or solicitation of an offer) to buy/se...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Question about acception rejection sampling - NOT R question

2007-07-13 Thread Charles C. Berry
On Fri, 13 Jul 2007, Leeds, Mark (IED) wrote:

 This is not related to R but I was hoping that someone could help me. I
 am reading the Understanding the Metropolis Hastings Algorithm
 paper from the American Statistician by Chip and Greenberg, 1995, Vol
 49, No 4. Right at the beginning they explain the algorithm for basic
 acceptance  rejection sampling in which you want to simulate a density
 from f(x) but it's not easy and you are able
 to generate from another density called h(x). The assumption is that
 there exists some c such that f(x) = c(h(x) for all x

 They clearly explain the steps as follows :

 1) generate Z from h(x).

 2) generate u from a Uniform(0,1)

 3) if u is less than or equal to f(Z)/c(h(Z) then return Z as the
 generated RV; otherwise reject it and try again.

 I think that, since f(Z)/c(h(z) is U(0,1), then u has the distrbution as
 f(Z)/c(h(Z).

 But, I don't understand why the generated and accepted Z's have the same
 density  as f(x) ?

 Does someone know where there is a proof of this or if it's reasonably
 to explain, please feel free to explain it.

The original reference to J. von Neumann's work, which no doubt has a 
proof, is in

http://en.wikipedia.org/wiki/Rejection_sampling

along with a suggestive graph.

Here is another graph that may give your intuition a boost:

f - function(x) dnorm(x)/2+dnorm(x,sd=0.1)/2
h - dnorm
Z - rnorm(100)
u - runif(100)
cee - f(0)/h(0) # as is obvious
u.lt.f.over.ch - u  f(Z)/cee/h(Z)
cutpts - seq(-5,5,by=0.1)
midpts - head(cutpts,n=-1)/2 + tail(cutpts,n=-1)/2
tab - table( !u.lt.f.over.ch, cut( Z, cutpts ) )
bp - barplot( tab )
lines( bp, f(midpts)*sum(u.lt.f.over.ch)*0.1,col='red',lwd=2)

Of course, there is some discreteness in this plot.
If you have a 1280x1024 or wider screen or want to zoom in on a pdf(), 
you might try 0.025 or even 0.0125 in place of 0.1.

Turn this graph on its side and think about what u.lt.f.over.ch is doing 
conditionally on Z, i.e. look at one bar and think about why it is split 
where it is.


HTH,

Chuck

 They authors definitely believe it's too trivial because they don't. The
 reason I ask is because, if I don't understand this then
 I definitely  won't understand the rest of the paper because it gets
 much more complicated.  I willing to track down the proof but I don't
 know where to look. Thanks.
 

 This is not an offer (or solicitation of an offer) to buy/se...{{dropped}}

 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


Charles C. Berry(858) 534-2098
 Dept of Family/Preventive Medicine
E mailto:[EMAIL PROTECTED]  UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 92093-0901

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.