Ratio-of-Uniforms Distribution in R -


i have exercise, in have create algorithm follows:

ratio of uniforms based on fact random variable x density f(x) can generate x desired density calculating x = u/v pair (u, v ) uniformly distributed in set

af = {(u,v):0 < v ≤  f(u/v)} 

random points can sampled uniformly in af rejection min- imal bounding rectangle, i.e., smallest possible rectangle contains af . given (u−, u+) × (0, v+) where

v+ = max f(x), x  u− = minx f(x), x  u+ = maxx f(x) 

then ratio-of-uniforms method consists of following simple steps:

  1. generate random number u uniformly in (u−, u+).

  2. generate random number v uniformly in (0, v+).

  3. set x ← u/v .

  4. if v 2 ≤ f(x) accept , return x.

  5. else try again.

my code far:

x <- cnorm(1, mean = 0, sd=1)  myrnorm <- function(pdf){      ## call rou() n times      pdf <- function(x) {exp(-x^2/2)}      } rou <- function(u, v) {    uplus <- 1    vplus <- 1    n <- 100    u <- runif(n, min=0, max=uplus)    v <- runif(n, min=0, max=vplus)    xi <- v/u    while(v < sqrt(xi)) {      if(v^2 <= xi)        return(xi)      } }   myx <- myrnorm(1000) hist(myx) 

but dont know how go on. im ´lost exercise. grateful advise.

following example 1 in page 8 of link , sample code, came solution:

ratiou <- function(nvals) {   h_x = function(x) exp(-x)   # u- b-, u+ b+ , v+ in example:   uminus = 0   uplus = 2/exp(1)   vplus = 1   x.vals <- null   <- 0   repeat {     <- i+1     u <- runif(1,0,vplus)     v <- runif(1,uminus,uplus)     x <- u/v     if(v^2 <= h_x(x)) {       tmp <- x     }     else {       next     }     x.vals <- c(x.vals,tmp)     if(length(x.vals) >= nvals) break   }   answer <- x.vals     answer }  sol = ratiou(1000)  par(mfrow=c(1,2)) hist(sol,breaks=50, main= "using ratiou",freq=f) hist(rexp(1000),breaks = 50, main="using rexp r",freq=f) par(mfrow=c(1,1))  par(mfrow=c(1,2)) plot(density(sol)) plot(density(rexp(1000))) par(mfrow=c(1,1)) 

a lot of code may optimized think enough this purpose. hope helps.


Comments

Popular posts from this blog

php - How to display all orders for a single product showing the most recent first? Woocommerce -

asp.net - How to correctly use QUERY_STRING in ISAPI rewrite? -

angularjs - How restrict admin panel using in backend laravel and admin panel on angular? -