Re: [R] Random number from density()

From: Mark Fisher <mark_at_markfisher.net>
Date: Mon 02 Apr 2007 - 16:56:08 GMT

Arianne ALBERT wrote:
> Hello,
>
> I'm writing some genetic simulations in R where I would like to place
> genes along a chromosome proportional to the density of markers in a
> given region. For example, a chromosome can be presented as a list of
> marker locations:
>
> Chr1<-c(0, 6.5, 17.5, 26.2, 30.5, 36.4, 44.8, 45.7, 47.8, 48.7, 49.2,
> 50.9, 52.9, 54.5, 56.5, 58.9, 61.2, 64.1)
>
> Where the numbers refer to the locations of markers along the
> chromosome. As you can see, there are a lot of markers around 50, but
> they are less dense elsewhere. I would like to take this information to
> randomly select a location on the chromosome to put a gene, where we are
> more likely to pick a location with high marker density (instead of
> using a uniform distribution between 0 and 64.1).
>
> So far I've used density(Chr1) to generate the probability density, but
> can this also be used to generate a random number? All the help I can
> find suggests getting the quantile function (the reciprocal of the
> integral of the pdf), however, it doesn't seem as though density()
> returns a pdf that can be manipulated in this way. Any suggestions?
>
> Thanks in advance,
>
> Arianne
>

I'm new to R and maybe I don't understand what you want, but here's a suggestion:

Chr1 <- c(0, 6.5, 17.5, 26.2, 30.5, 36.4, 44.8, 45.7, 47.8, 48.7, 49.2,

        50.9, 52.9, 54.5, 56.5, 58.9, 61.2, 64.1) n <- length(Chr1)-1
fun <- approxfun(0:n/n, Chr1)
ran <- fun(runif(10000))
hist(ran)

--Mark.



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. Received on Tue Apr 03 03:03:17 2007

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Mon 02 Apr 2007 - 17:30:37 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.