Greetings, Let me present my modest proposition as function with an input and an output, it is based on the valuable solution of Alberto Monteiro. I have choosen delta superior or equal to 0 wich still provide a REAL solution :
prob <- function(n){ a <- runif(n) b <- runif(n) c <- runif(n) delta <- b^2 - 4 * a * c # a vector of deltas sum(delta >= 0) / length(delta) return(p = sum(delta >= 0) / length(delta)) } # To test > prob(1000) [1] 0.232 > prob(100000) [1] 0.25308 > prob(10000000) [1] 0.2544715 Best regards, Souleymane N'DOYE, M.Sc. Statistician Engineer & Decision Support System Consultant, French. www.labstatconseil.com [EMAIL PROTECTED] P.O. Box 1601 00606 Sarit Center, Nairobi, Kenya Mobile : +254 736 842 478. On Wed, Aug 20, 2008 at 3:56 PM, Alberto Monteiro <[EMAIL PROTECTED]>wrote: > > Jaap Van Wyk wrote: > > > > I would appreciate any help with the following. > > Problem: Suppose A, B and C are independent U(0,1) random variables. > > What is the probability that A(x^2) + Bx + C has real roots? I have > > done the theoretical work and obtained an answer of 1/9 = 0.1111. > > Now I want to show my students to get this in R with simulation. > > Below are two attemps, both giving the answer to be about 0.26. > > > > Could anybody please help me with providing a more elegant way to do > > this? (I am still learning R and trying to get my students to learn > > it as well. I know there must be a better way to get this.) I must > > be doing something wrong ? > > > Always think that R is vector-oriented, so you should think > in terms of vectors. > > A simple solution for your problem would be: > > n <- 10000 > # n <- 10000? Make n <- 1000000 ! > n <- 1000000 > a <- runif(n) # a vector of n unifs > b <- runif(n) > c <- runif(n) > delta <- b^2 - 4 * a * c # a vector of deltas > > # The answer then is how many deltas are non-negative: > > sum(delta > 0) / length(delta) > > # Explanation > b^2 is the vector of the squares of b's > a * c does the pointwise product of vectors > delta > 0 is a logical array with TRUE or FALSE > sum(delta > 0) coerces TRUE to 1 and FALSE to 0 > length(delta) is the length of delta (n) > > Alberto Monteiro > > ______________________________________________ > R-help@r-project.org 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. > [[alternative HTML version deleted]]
______________________________________________ R-help@r-project.org 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.