R help - Sep 2012 - simulate from conditional distribution

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?
        Thank you very much.
                   hannah

	[[alternative HTML version deleted]]

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.