An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20120913/4e6f2c23/attachment.pl>
simulate from conditional distribution
3 messages · Li Li, Thomas Lumley, Petr Savicky
On Fri, Sep 14, 2012 at 9:02 AM, li li <hannah.hlx at gmail.com> wrote:
Dear all,
Y, X are bivariate normal with all the parameters known.
I would like to generate numbers from the distribution Y | X > c
where c is a constant.
Does there exist an R function generating
random numbers from such a distribution?
Not directly, as far as I know, but you can easily simulate X|X>c by
transforming uniform random numbers using the inverse CDF, and Y|X=x
is univariate Normal with mean linear in x and variance independent of
x.
-thomas
Thomas Lumley Professor of Biostatistics University of Auckland
On Thu, Sep 13, 2012 at 05:02:54PM -0400, li li wrote:
Dear all,
Y, X are bivariate normal with all the parameters known.
I would like to generate numbers from the distribution Y | X > c
where c is a constant.
Does there exist an R function generating
random numbers from such a distribution?
Hi.
One of the options, how to generate such numbers, is the rejection
method. This works well, if Pr(X > c) is not too small, but may be
very bad otherwise.
For generating a single number, try the following
library(mvtnorm)
sigma <- rbind(c(3, 1), c(1, 2))
cc <- 2
while (1) {
v <- rmvnorm(1, sigma=sigma)
if (v[2] > cc) break
}
x <- v[1]
For a larger number, try something like the following
n <- 500
for (k in 2^(1:10)) {
v <- rmvnorm(k*n, sigma=sigma)
x <- v[v[, 2] > cc, 1]
if (length(x) >= n) break
}
stopifnot(length(x) >= n)
x <- x[1:n]
Hope this helps.
Petr Savicky.