# Section 2; Question 1:
#
p = 1/3250
n = 6000
print( dbinom( 0, n, p ) )
print( sprintf('n= %d; p= %0.6f, n*p= %.3f', n, p, n*p) )
print( dpois( 0, n*p ) )


# Section 2; Question 2:
#
p = 905/289411
n = 10
print( 1 - dbinom( 0, n, p ) )
print( 1 - dpois( 0, n*p ) )


# Section 2; Question 3:
#
p = 1/365
n = 500
print( dpois( 0, n*p ) + dpois( 1, n*p ) )


# Section 2; Question 4:
#
p = 1/10000
n = 20000
print( dpois( 3, n*p ) )
for( k in 1:11 ){
    print( sprintf('P(X >= %2d) = %.6f', k, 1-ppois( k-1, n*p )) )
}


# Section 2; Question 5:
#
p = 0.01
n = 10
print( sprintf('binomia1 prob= %.6f', 1 - dbinom( 0, n, p ) ) )
print( sprintf('poisson prob= %.6f', 1 - dpois( 0, n*p ) ) )



# Section 2; Question 6:
p = 1/150
n = 120
print( sprintf( 'the probability of three or more claims in a given year= %.6f', 1-ppois( 2, n*p ) ) )
#print( sprintf( 'the probability of three or more claims in a given year= %.6f', 1-pbinom( 2, n, p ) ) )


# Section 2; Question 7:
#
p = 1/200
n = 120
print( sprintf( 'the probability loosing two or more bags in a given year= %.6f', 1-ppois( 1, n*p ) ) )
#print( sprintf( 'the probability loosing two or more bags in a given year= %.6f', 1-pbinom( 1, n, p ) ) )


# Section 2; Question 8:
#
p = 1/10^6
n = 9500
print( sprintf( 'the probability of two or more cases (under H0)= %.6f', 1-ppois( 1, n*p ) ) )
#print( sprintf( 'the probability of two or more cases (under H0)= %.6f', 1-pbinom( 1, n, p ) ) )


# Section 2; Question 9:
#
print( log(2)/(100*10^9) )


# Section 2; Question 10:
#
num_of_fatalities = 0:4
freq = c( 109, 65, 22, 3, 1 )
k_bar = sum( num_of_fatalities * freq ) / sum( freq )
empirical = freq / sum( freq )
model = dpois( 0:4, k_bar ) 
print( data.frame( num_of_fatalities=num_of_fatalities, empirical=empirical, model=model ) )


# Section 2; Question 11:
#
num_changes = 0:3
freq = c( 237, 90, 22, 7 )
k_bar = sum( num_changes * freq ) / sum( freq )
empirical = freq / sum( freq )
model = dpois( num_changes, k_bar )
print( data.frame( num_changes=num_changes, empirical=empirical, model=model ) )
print( model * sum( freq ) )


# Section 2; Question 12:
#
source('chap_4_sect_2_question_12_data.R')
T = as.data.frame( table( DF$Bags_Lost ) )
colnames(T) = c('Bags_Lost', 'Frequency')
print(T)
k_bar = mean(DF$Bags_Lost)
empirical = T$Frequency/sum( T$Frequency )
model = dpois( min(DF$Bags_Lost):max(DF$Bags_Lost), k_bar )
print( data.frame( Bags_Lost=T$Bags_Lost, empirical=empirical, model=model ) )


# Section 2; Question 13:
#
number_of_countries = 0:4
freq= c( 82, 25, 4, 0, 2 )
k_bar = sum( number_of_countries * freq ) / sum( freq )
empirical = freq / sum( freq )
model = dpois( number_of_countries, k_bar )
print( data.frame( number_of_countries=number_of_countries, empirical=empirical, model=model ) )
print( model * sum( freq ) ) # the expected frequency


# Section 2; Question 14:
#
number_of_deaths = 0:9
freq = c( 162, 267, 271, 185, 111, 61, 27, 8, 3, 1 )
k_bar = sum( number_of_deaths * freq ) / sum( freq )
empirical = freq / sum( freq )
model = dpois ( number_of_deaths, k_bar )
print( data.frame( number_of_deaths=number_of_deaths, empirical=empirical, model=model ) )
print( model * sum( freq ) ) # the expected frequency


# Section 2; Question 15:
#
num_of_infestations = 0:7
freq = c( 55, 20, 21, 1, 1, 1, 0, 1 )
k_bar = sum( num_of_infestations * freq ) / sum( freq )
empirical = freq / sum( freq )
model = dpois( num_of_infestations, k_bar )
print( data.frame( num_of_infestations= num_of_infestations, empirical=empirical, model= round( model, 2 ) ) )
print( model * sum( freq ) ) # the expected frequency


# Section 2; Question 16:
#
lambda = 1/5 # once in five hours
print( sprintf( 'the probability of less than three breakdowns= %.6f', ppois( 2, lambda*8 ) ) )


# Section 2; Question 17:
#
lambda = 1.5/10
print( sprintf( 'the probability of more than two errors= %.6f', 1-ppois( 2, lambda*30 ) ) )


# Section 2; Question 18:
#
lambda = log(3)
print( sprintf( 'the probability of two or more hits= %.6f', 1-ppois( 1, lambda ) ) )


# Section 2; Question 19:
#
lambda = 1/10
print( sprintf( 'the probability of three or more flaws= %.6f', 1-ppois( 2, lambda*5*8 ) ) )


# Section 2; Question 20:
#
lambda = 482/60 # number of particles per two minute period
print( sprintf( 'Method 1: the probability of three particles in two minutes= %.6f', dpois( 3, lambda*1 ) ) )
lambda = 482/(2*60) # number of particles per one minute period
print( sprintf( 'Method 2: the probability of three particles in two minutes= %.6f', dpois( 3, lambda*2 ) ) )


# Section 2; Question 21:
#
lambda = 0.1 # per day injury rate
print( sprintf( 'the probability of two accidents in five days= %.6f', dpois( 2, lambda*5) ) )
print( sprintf( 'the probability of four accidents in two weeks (10 days)= %.6f', dpois( 4, lambda*10 ) ) )


# Section 2; Question 22:
#
lambda = 2
print( dpois( 4, lambda ) )


# Section 2; Question 26:
#
lambda = 2.5 # crashes per year
print( 1-ppois(3, lambda) )
print( pexp( 3/12, lambda ) )
 

# Section 2; Question 27:
#
lambda = 0.1
print( 1-pexp( 7, rate=lambda ) )


# Section 2; Question 28:
#
lambda = 0.027
1 - exp(-40*lambda) * ( 40*lambda + 1 )


# Section 2; Question 29:
#
lambda = 1.1 # per hundred hours per light
p = pexp( 75/100, lambda ) # bulb fails i.e. P( Y < 75/100 )
n = 50
print( sprintf( 'the expected number to fail (a binomial RV)= %.6f', n*p ) )