###################################################
# Computer Lab 2 – F71SM
###################################################
## In the tasks below, some (but not all) R code is given
## for guidance. In some places the code is also incomplete,
## so you will need to figure out how to fill in the details.
## Make sure you use the help menus in R (e.g.”?mean” will
## open up a help window for “mean()”).
###################################################
## Task 1
# Tutorial 1, Question 2:
# Repeat the calculations
# The data can be found in file T1Q2.RData
load(“T1Q2.RData”)
objects()
sum(x)
…
q = quantile(x,probs=c(0.25,0.75))
…
# Note that R by default uses a different definition of quantiles
# The definition in lectures notes can be obtained using
# type=6 in quantile()
###################################################
## Task 2
# Tutorial 1, Question 7:
# Confirm the answers by computing the mean, sd and skewness
# of the two samples
a = rep(c(1:4),c(80,60,40,20))
b = …
…
# For skewness, use function from lab 1
…
###################################################
## Task 3
# Worked example 2.4:
# Confirm the answers for all values of k
p = numeric(5)
# Note that p[0] is not defined for a vector in R
for(k in 0:4){
p[k+1] = choose(4,k)*choose(7,(6-k))/choose(11,6)
}
p
###################################################
## Task 4
# Simulate 10,000 values throws of a fair six-sided die
# Let X be the r.v. giving the score of the die
# Compute (approximately) the mean and variance of X
# (this is also worked example 4.1).
# A single value of X can be simulated as (check):
trunc(6*runif(1, 0, 1)) + 1
# So 10^4 values simulated as
x = …
mean(x); var(x)
# Compute (approximately) P(X<=3)
length(x[x<=3]) / length(x)
# Repeat these calculations to see how the anwsers
# change slightly
###################################################
## Task 5
# Worked example 3.3:
# Simulate 10,000 values of r.v. Y, where Y is
# number of throws until we get a "6"
# Confirm the mean and variance of Y
# A sigle value of y can be simulated as follows:
# Generate throws of the die as in Task 3
# and stop when a 6 is thrown
y = 0
x = 0
while(x!=6){
x = trunc(6*runif(1, 0, 1)) + 1
y = y+1
}
y
# Use for loop to simulate 10^4 values of Y
y = numeric(10^4)
for (i in 1:10^4){
...
}
mean(y)
var(y)
# Again, repeat these calculations to see how the anwsers
# change slightly
###################################################
## Task 6
# Worked examples 3.5, 3.6:
# Plot the pdf and cdf of r.v. X
x1 = seq(1.001,4,by=0.01)
f.x1 = ...
F.x1 = ...
# Change the plotting window layout
par(mfrow=c(1,2))
plot(x1,f.x1,type="l",xlab="x", ylab="f(x)", main="PDF")
plot(x1,F.x1,type="l",xlab="x", ylab="F(x)", main="CDF")
# Can also change the limits of the x axis:
plot(x1,f.x1,type="l",xlim=c(0,3.3), xlab="x", ylab="f(x)", main="PDF")
plot(x1,F.x1,type="l",xlim=c(-1,3.3),xlab="x", ylab="F(x)", main="CDF")
###################################################
## Task 7
# Worked example 4.2
# X ~ Binomial(10,0.2)
# Calculate pmf f(3), cdf F(3)
...
#Simulate 10^3 values of X and confirm E(X), Var(X)
x = rbinom(...)
mean(x); var(x)
###################################################
## Task 8
# Worked example 4.3
# X ~ Poisson(7)
# Calculate pmf f(9), P(X>=8)
dpois(…)
1-ppois(…)
# Simulate N=10^4 values of X and confirm E(X), Var(X)
# Use rpois()
…
# Repeat above with N=10^2. Comment.
###################################################
## Task 9
# Worked example 4.5
# X ~ Binomial(200,0.02)
# Calcualate P(X = 5 or 6) using the binomial above and a Poisson approximation
pbinom(…) – pbinom(…)
ppois(…) – ppois(…)
###################################################
## Task 10
# Worked example 4.5
# X ~ Geo(0.3)
# Calcualate P(X = 1, 2 or 3)
…
###################################################
## Task 11
# Worked example 4.6
# X ~ Exp(0.2)
# Calcualate P(6 < X < 10), P(6 <= X < 10)
# Use pexp()
...
# Simulate N=10^4 values of X and confirm E(X), Var(X), P(6 < X < 10)
# Use rexp() ...
...
# Plot the pdf and cdf of r.v. X
x1 = seq(...,by=0.01)
f.x1 = ...
F.x1 = ...
# Change the plotting window layout
par(mfrow=c(1,2))
plot(...)
plot(...)
###################################################
## Task 12
# Worked example 4.7
# X ~ N(500,10)
# Calculate P(493 < X < 505), 99th quantile
# Use pnorm(), qnorm()
...
###################################################
## Task 13
# X ~ Gamma(4,3)
# Calculate 0.95 quantile
# Use qgamma()
...
# Y = 6*X ~ Chi-square(8)
# Calculate 0.95 quantile of X and confirm that it is the same as above
# Use qchisq()