R/rigaussian.R

In mixbox: Observed Fisher Information Matrix for Finite Mixture Model

rigaussian <- function(n, alpha, beta)
{
      Y <- rep( NA, n )
     z2 <- rnorm(n)^2
  ratio <- alpha + alpha^2*z2/(2*beta) - alpha/(2*beta)*sqrt( 4*alpha*beta*z2 + (alpha*z2)^2 )
      y <- alpha^2/ratio
      u <- runif(n)
      j <- which( ratio < alpha*u/(1 - u) )
   Y[j] <- ratio[j]
return( Y )
}