9.17.2017 Homework

Lavine. 1.10. Problem 2. Simulating Dice Rolls.

1. simulate 6000 dice rolls. Count the number of 1’s, 2’s, . . . , 6’s.

rolls<-sample(6, 6000, replace=T)
for(i in 1:6) print (sum(rolls==i)) 

[1] 985
[1] 931
[1] 1011
[1] 991
[1] 1062
[1] 1020

simDiceTot

[1] 1012  967  992 1005  996 1028

2. You expect about 1000 of each number. How close was your result to what you expected?

100*(simDiceTot-(rep(1000,6)))/1000

[1]  1.2 -3.3 -0.8  0.5 -0.4  2.8

3. About how often would you expect to get more than 1030 1’s? Run an R simulation to estimate the answer.

rollvec=rep(0,100000)
for(i in 1:100000){ 
  roll=sample(6,6000,replace=T) 
  if(sum(roll==1)>1030) rollvec[i]=1 
}
sum(rollvec)/100000

[1] 0.14533

= 0.14518

simDice<-matrix(NA,nrow=6,ncol=10000)
for(i in 1:10000){
  simDice[,i]<-apply(rmultinom(n=6000, size=1, prob=c(rep(1,6))),1,sum)
} 
length(simDice[1,simDice[1,]>1030])/dim(simDice)[2] 

[1] 0.1475

= 0.141 or 14% of time

Lavine 1.10. Problem 37. Solve analytically or by simulation.

Ecologists are studying salamanders in a forest. There are two types of forest. Type A is conducive to salamanders while type B is not. They are studying one forest but don’t know which type it is. Types A and B are equally likely.

During the study, they randomly sample quadrats. (A quadrat is a square-meter plot.) In each quadrat they count the number of salamanders. Some quadrats have poor salamander habitat. In those quadrats the number of salamanders is 0. Other quadrats have good salamander habitat. In those quadrats the number of salamanders is either 0, 1, 2, or 3, with probabilities 0.1, 0.3, 0.4, and 0.2, respectively. (Yes, there might be no salamanders in a quadrat with good habitat.) In a type A forest, the probability that a quadrat is good is 0.8 and the probability that it is poor is 0.2. In a type B forest the probability that a quadrat is good is 0.3 and the probability that it is poor is 0.7.

1. On average, what is the probability that a quadrat is good?

p(Q=g)=p(g|A)P(A) + p(g|B)p(B) =0.80.5+0.30.5 =0.55

2. On average, what is the probability that a quadrat has 0 salamanders, 1 salamander, 2 salamanders, 3 salamanders?

p(S=x) =p(S=x|g)p(g)+p(S=x|ng)p(ng) p(S=0) = 0.10.55 + 1(0.45) = 0.505 p(S=1) = 0.30.55 + 0(0.45) = 0.165 p(S=2) = 0.40.55 + 0(0.45) = 0.22 p(S=3) = 0.20.55 + 0(0.45) = 0.11

Note: 0.11+0.22+0.165+0.505 = 1

3. The ecologists sample the first quadrat. It has 0 salamanders. What is the probability that the quadrat is good?

P(g|S=0) = p(S=0|g)p(g) / ( p(S=0|g)p(g)+p(S=0|ng)p(ng) ) = (0.10.55) / ( 0.10.55 + 1*0.45 ) = 0.1089109

4. Given that the quadrat had 0 salamanders,what is the probability that the forest is type A?

P(F=A|S=0) = p(S=0|g)p(g|A)P(A)+p(S=0|ng)p(ng|A)P(A) / ( p(S=0|g)p(g|A)P(A)+p(S=0|ng)p(ng|A)P(A) + p(S=0|g)p(g|B)P(B)+p(S=0|ng)p(ng|B)P(B) ) = (0.10.80.5 + 10.20.5) / ( 0.10.80.5 + 10.20.5 + 0.10.30.5 + 10.70.5) = 0.2772277

5. Now the ecologists prepare to sample the second quadrat. Given the results from the first quadrat, what is the probability that the second quadrat is good?

p(Q=g)=p(g|A)P(A)+p(g|B)p(B) =0.80.2772277+0.3(1-0.2772277) =0.4386138

6. Given the results from the first quadrat, what is the probability that they find no salamanders in the second quadrat?

p(S=0) =p(S=0|g)p(g)+p(S=0|ng)p(ng) p(S=0) = 0.10.4386138 + 1(1-0.4386138) = 0.6052476

Lavine 1.10 Problem 44.

This exercise is based on a computer lab that another professor uses to teach the Central Limit Theorem. It was originally written in MATLAB but here it’s translated into R. Enter the following R commands: u <- matrix ( runif(250000), 1000, 250 ) y <- apply ( u, 2, mean )

These create a 1000 x 250 (a thousand rows and two hundred fifty columns) matrix of random draws, called u and a 250-dimensional vector y which contains the means of each column of U.

Now enter the command hist(u[,1]). This command takes the first column of u (a column vector with 1000 entries) and makes a histogram. Print out this histogram and describe what it looks like. What distribution is the runif command drawing from?

Now enter the command hist(y). This command makes a histogram from the vector y. Print out this histogram. Describe what it looks like and how it differs from the one above. Based on the histogram, what distribution do you think y follows?

You generated y and u with the same random draws, so how can they have different distributions? What’s going on here?

runif is sampling from the uniform distribution.

It looks like a normal distribution.

The Central Limit Theorem.

This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Cmd+Shift+Enter.

plot(cars)

Add a new chunk by clicking the Insert Chunk button on the toolbar or by pressing Cmd+Option+I.

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the Preview button or press Cmd+Shift+K to preview the HTML file).