Task 1

getwd()
## [1] "C:/Users/Caleb/Documents/Civil Engineering degree coursework/Applied Statistical Methods/Labs/Lab 5"

Task 2

mybin=function(iter=100,n=10, p=0.5){ 
# make a matrix to hold the samples
#initially filled with NA's
sam.mat=matrix(NA,nr=n,nc=iter, byrow=TRUE)
#Make a vector to hold the number of successes in each trial
succ=c()
for( i in 1:iter){
#Fill each column with a new sample
sam.mat[,i]=sample(c(1,0),n,replace=TRUE, prob=c(p,1-p))
#Calculate a statistic from the sample (this case it is the sum)
succ[i]=sum(sam.mat[,i])
}
#Make a table of successes
succ.tab=table(factor(succ,levels=0:n))
#Make a barplot of the proportions
barplot(succ.tab/(iter), col=rainbow(n+1), main="Binomial simulation", xlab="Number of successes")
succ.tab/iter
}
mybin(iter=100,n=10, p=0.7)

## 
##    0    1    2    3    4    5    6    7    8    9   10 
## 0.00 0.00 0.00 0.00 0.04 0.11 0.17 0.27 0.24 0.12 0.05
mybin(iter=200,n=10, p=0.7)

## 
##     0     1     2     3     4     5     6     7     8     9    10 
## 0.000 0.000 0.000 0.010 0.055 0.090 0.160 0.300 0.215 0.150 0.020
mybin(iter=500,n=10, p=0.7)

## 
##     0     1     2     3     4     5     6     7     8     9    10 
## 0.000 0.000 0.000 0.020 0.042 0.112 0.186 0.268 0.236 0.108 0.028
mybin(iter=1000,n=10, p=0.7)

## 
##     0     1     2     3     4     5     6     7     8     9    10 
## 0.000 0.000 0.000 0.007 0.036 0.107 0.189 0.248 0.258 0.125 0.030
mybin(iter=10000,n=10, p=0.7)

## 
##      0      1      2      3      4      5      6      7      8      9     10 
## 0.0000 0.0000 0.0012 0.0088 0.0378 0.1008 0.1965 0.2691 0.2342 0.1229 0.0287
dbinom(0,10,0.7)
## [1] 5.9049e-06
dbinom(1,10,0.7)
## [1] 0.000137781
dbinom(2,10,0.7)
## [1] 0.001446701
dbinom(3,10,0.7)
## [1] 0.009001692
dbinom(4,10,0.7)
## [1] 0.03675691
dbinom(5,10,0.7)
## [1] 0.1029193
dbinom(6,10,0.7)
## [1] 0.2001209
dbinom(7,10,0.7)
## [1] 0.2668279
dbinom(8,10,0.7)
## [1] 0.2334744
dbinom(9,10,0.7)
## [1] 0.1210608
dbinom(10,10,0.7)
## [1] 0.02824752

Yes, each value is approximately correct.

Task 3

bag <- c(1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0)
ind=sample(1:20,5,replace=FALSE)
bag[ind]
## [1] 1 1 1 1 0
ind2=sample(1:20,5,replace=TRUE)
bag[ind2]
## [1] 1 1 0 0 1
myhyper=function(iter=100,N=20,r=12,n=5){
# make a matrix to hold the samples
#initially filled with NA's
sam.mat=matrix(NA,nr=n,nc=iter, byrow=TRUE)
#Make a vector to hold the number of successes over the trials
succ=c()
for( i in 1:iter){
#Fill each column with a new sample
sam.mat[,i]=sample(rep(c(1,0),c(r,N-r)),n,replace=FALSE)
#Calculate a statistic from the sample (this case it is the sum)
succ[i]=sum(sam.mat[,i])
}
#Make a table of successes
succ.tab=table(factor(succ,levels=0:n))
#Make a barplot of the proportions
barplot(succ.tab/(iter), col=rainbow(n+1), main="HYPERGEOMETRIC simulation", xlab="Number of successes")
succ.tab/iter
}
myhyper(iter=100,n=5, N=20,r=12)

## 
##    0    1    2    3    4    5 
## 0.00 0.04 0.17 0.46 0.26 0.07
myhyper(iter=200,n=5, N=20,r=12)

## 
##     0     1     2     3     4     5 
## 0.000 0.060 0.225 0.425 0.220 0.070
myhyper(iter=500,n=5, N=20,r=12)

## 
##     0     1     2     3     4     5 
## 0.004 0.050 0.240 0.382 0.270 0.054
myhyper(iter=1000,n=5, N=20,r=12)

## 
##     0     1     2     3     4     5 
## 0.008 0.059 0.217 0.407 0.251 0.058
myhyper(iter=10000,n=5, N=20,r=12)

## 
##      0      1      2      3      4      5 
## 0.0036 0.0549 0.2334 0.3949 0.2608 0.0524
dhyper(x=0:5, m=12, n=8, k=5)
## [1] 0.003611971 0.054179567 0.238390093 0.397316821 0.255417957 0.051083591

Task 4

I think it forms a random sample between 1:10 an n number of times. Then it runs that process an iter number of times and outputs the relative frequencies of 1-10 onto a barplot.

mysample=function(n, iter=10,time=0.5){
for( i in 1:iter){
#make a sample
s=sample(1:10,n,replace=TRUE)
# turn the sample into a factor
sf=factor(s,levels=1:10)
#make a barplot
barplot(table(sf)/n,beside=TRUE,col=rainbow(10), 
main=paste("Example sample()", " iteration ", i, " n= ", n,sep="") ,
ylim=c(0,0.2)
)

#release the table
Sys.sleep(time)
}
}

mysample(n=1000, iter=30, time=1)

I see a barplot for each iter value that contains the relative frequencies of the random sample from 1-10.

Task 5

choose(8,4)
## [1] 70
ppois(4,2)
## [1] 0.947347
dnbinom(7,3,0.4)
## [1] 0.06449725
pbinom(8,15,0.4)
## [1] 0.9049526

Task 6

MATH4753GRAY::mybin(iter=1000,n=18, p=0.3)

## 
##     0     1     2     3     4     5     6     7     8     9    10    11    12 
## 0.001 0.018 0.052 0.111 0.167 0.181 0.209 0.134 0.078 0.029 0.012 0.007 0.001 
##    13    14    15    16    17    18 
## 0.000 0.000 0.000 0.000 0.000 0.000