Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

STAT 2011 Probability and Estimation Theory – Semester 1, 2022

Computer Practical Sheet Week 2

Producing a report

Please carefully read the document on preparing reports using RStudio on the Canvas page before attempting this computer practical.

Computer Problems for Week 2

1. We are going to investigate the following question:

What is more likely, getting  at least once in five rolls of one die, or getting  at least once in 25 rolls of two dice?

(a)  The following R command creates a pseudo-random sample of size 5000  with re-

placement from the population {1, 2, 3, 4, 5, 6}:

set.seed(2018)

rolls1 =  sample(x=1:6,  size=5000, replace=TRUE)

To convince yourself the correct thing has been created, also execute table(rolls1) to see a frequency table of the 5000 values just generated and look at the histogram through hist(rolls1,breaks=0:6)

(b) Form the vector rolls1 into a 1,000-by-5 matrix called five.rolls using the matrix()

function.  For example, if y=1:10 then M=matrix(y,nrow=2,ncol=5) forms it into the matrix

1 3  5 7  9

2 4  6 8  10

while M=matrix(y,nrow=2,ncol=5,byrow=TRUE) gives

1  2 3 4  5

6 7 8  9  10

although either way is fine in this case (note for future reference though!).

(c) Each row of your matrix represents five rolls of a die. We wish to count how many rows have at least one 1 in them (this corresponds to a ).  Equivalently, we can compute the minimum value of each row, and count how many of these are 1.  Use the apply() function to create a vector called min.roll consisting of the 1,000 row minima of the matrix four.rolls.  Hint:  if M is a matrix then apply(M,1,sum) gives the row sums of M, apply(M,2,min) the column minima of M, etc..

(d)  Counting how many and/or what proportion of elements of min.roll are equal to

1 can be achieved using a logical comparison.   The vector min.roll==1 consists of 1,000 TRUE’s and FALSE’s, with TRUE corresponding to a minimum of 1.  These TRUE’s and FALSE’s are subsequently interpreted as 1’s and 0’s respectively, if need be. So sum(min.roll==1) gives us the count.

(e) We now do something similar for the second event.

i. create a vector called rolls2 consisting of a random sample of size 50,000 with replacement from {1, 2, 3, 4, 5, 6};

ii. form the vector rolls2 into a 25,000-by-2 matrix called two.rolls;

iii.  occurs if and only if the sum is 2; thus form a vector called  sum.rolls consisting of 25,000 row sums of the matrix two.rolls;

iv. form this vector of sums into a 25-by-1,000 matrix called two.dozen;

v. each column now corresponds to 25 rolls of two dice; form a vector min.pair of column minima, and

vi. finally count how many of these are equal to 2.

(f)  Convert your counts obtained here and in part (d) into estimates of the corresponding

probabilities we are interested in; call them p1.est and p2.est.  Based on these, what can you say in response to the question posed at the beginning of the exercise?

(g) How accurate are your estimates? We can get an idea of the accuracy of the estimates

using a little simulation. Replicate the whole procedure 25 times using a for-loop:

results1  <-  0

results2  <-  0

for  (i  in  1:25){

...                                      # your  commands from q1a to q1c go here

results1[i]  <-  sum(min.roll==1)

...                                      # your  commands from q1e parts  (i)-(v) go here

results2[i]  <-  sum(min.pair==2)

}

(h)  Define

prob.ests1  <- results1/1000

and similarly prob.ests2; also

se1  <-  sd(prob.ests1)

and similarly se2. These give an idea of how accurate the procedures are in general, and also how reliable our single estimates p1.est and p2.est are in particular. We call these the (estimated) standard errors of our estimates.

(i)  Compute the actual probabilities of these two events, assuming each possible se- quence of outcomes is equally likely in each of the two cases.  (Hint:  consider the complement). Call them p1 and p2.

(j)  Compute

abs(p1.est  - p1)/se1

and do the same for case 2.  In both cases your estimate should be within 2 or 3 standard errors of the true value. Did this occur?