The Art of R Programming

> emax

function(nreps) {

x <- rnorm(2*nreps)

maxxy <- pmax(x[1:nreps],x[(nreps+1):(2*nreps)])

return(mean(maxxy))

}

Here, we generated doublenrepsvalues. The firstnrepsvalue simulates

X, and the remainingnrepsvalue represents Y. Thepmax()call then com-

putes the pair-wise maxima that we need. Again, note the contrast here

betweenmax()andpmax(), the latter producing pair-wise maxima.

8.6.2 Obtaining the Same Random Stream in Repeated Runs.............

According to the R documentation, all random-number generators use

32-bit integers for seed values. Thus, other than round-off error, the same

initial seed should generate the same stream of numbers.

By default, R will generate a different random number stream from run

to run of a program. If you want the same stream each time—important in

debugging, for instance—callset.seed(), like this:

> set.seed(8888) # or your favorite number as an argument

8.6.3 Extended Example: A Combinatorial Simulation...................

Consider the following probability problem:

Three committees, of sizes 3, 4 and 5, are chosen from 20 people.

What is the probability that persons A and B are chosen for the

same committee?

This problem is not hard to solve analytically, but we may wish to check

our solution using simulation, and in any case, writing the code will demon-

strate how R’s set operations can come in handy in combinatorial settings.

Here is the code:

1 sim <- function(nreps) {
2 commdata <- list() # will store all our info about the 3 committees
3 commdata$countabsamecomm <- 0
4 for (rep in 1:nreps) {
5 commdata$whosleft <- 1:20 # who's left to choose from
6 commdata$numabchosen <- 0 # number among A, B chosen so far
7 # choose committee 1, and check for A,B serving together
8 commdata <- choosecomm(commdata,5)

Doing Math and Simulations in R 205