# Problems on EPage 468
#

# Section 2; Question 1:
#
source('chap_9_sect_2_question_1_data.R')

X = DF[ DF$Education==1, ] # funds used for education
n = dim(X)[1]
xbar = mean(X$Percent_Profit)
sx2 = var(X$Percent_Profit)

Y = DF[ DF$Education==0, ]
m = dim(Y)[1]
ybar = mean(Y$Percent_Profit)
sy2 = var(Y$Percent_Profit)

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 1-pt(t, n+m-2) # one-sided test

alpha = 0.01
t_crit = qt( 1-alpha, n+m-2 )

print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 2:
#
n = 77
xbar = -4.7
sx2 = 7.05^2
m = 79
ybar = -1.6
sy2 = 5.36^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = pt(t, n+m-2)

alpha = 0.05
t_crit = qt( alpha, n+m-2 )
print(sprintf('t= %.6f; p-value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 3:
#
n = 476
xbar = 189.0
sx2 = 34.2^2
m = 592
ybar = 177.2
sy2 = 33.3^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 1-pt(t, n+m-2) # one-sided test

alpha = 0.05
t_crit = qt( 1-alpha, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 4:
#
n = 1126
xbar = 491
sx2 = 119^2
m = 5042
ybar = 498
sy2 = 129^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = pt(t, n+m-2) # one-sided test

alpha = 0.05
t_crit = qt( alpha, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 5:
#
n = 93
xbar = 4.17
sx2 = 3.7^2
m = 28
ybar = 4.61
sy2 = 4.28^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 2*(1-pt(abs(t), n+m-2)) # two-sided test

alpha = 0.01
t_crit = qt( 1-alpha/2, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 6:
#
source('chap_9_sect_2_question_6_data.R')

X = DF[ DF$Alcohol==1, ] # funds used for education
n = dim(X)[1]
xbar = mean(X$Age)
sx2 = var(X$Age)

Y = DF[ DF$Alcohol==0, ]
m = dim(Y)[1]
ybar = mean(Y$Age)
sy2 = var(Y$Age)

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = pt(t, n+m-2) # one-sided test

alpha = 0.05

t_crit = qt( alpha, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 7:
#
n = 4
xbar = 1133.0
sx2 = 266.9^2
m = 4
ybar = 1013.5
sy2 = 224.3^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 2*(1-pt(abs(t), n+m-2)) # two-sided test

alpha = 0.05
t_crit = qt( 1-alpha/2, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 8:
#
source('chap_9_sect_2_question_8_data.R')

X = DF[ DF$Sex=='F',]
n = dim(X)[1]
xbar = mean(X$Halflife)
sx2 = var(X$Halflife)

Y = DF[ DF$Sex=='M', ]
m = dim(Y)[1]
ybar = mean(Y$Halflife)
sy2 = var(Y$Halflife)

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 2*(1-pt(abs(t), n+m-2)) # two-sided test

alpha = 0.01
t_crit = qt( 1-alpha/2, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 9:
#
n = 31
xbar = 3.10
sx2 = 1.469^2
m = 57
ybar = 2.43
sy2 = 1.35^2

sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 )

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 2*(1-pt(abs(t), n+m-2)) # two-sided test

alpha = 0.05
t_crit = qt( 1-alpha/2, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 10:
#
n = 10
xbar = 2.1*60 # in minutes
sx2 = 12^2
m = 10
ybar = 1.6*60
sy2 = 16^2
sp2 = ( (n-1)*sx2 + (m-1)*sy2 ) / ( n+m-2 ) # assume sigma_X = sigma_Y

t = ( xbar - ybar )/sqrt( sp2*(1/n + 1/m) )
p_value = 1-pt(t, n+m-2) # one-sided test

alpha = 0.05
t_crit = qt( 1-alpha, n+m-2 )
print(sprintf('t= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 11:
#
n = 6
m = 11
sp2 = 15.3^2

alpha = 0.01
t_crit = qt( 1-alpha/2, n+m-2 )
smallest_abs_diff = sqrt( sp2*(1/n + 1/m) ) * t_crit
print(sprintf('Part (a): |xbar-ybar|_min= %8.3f; t_crit= %8.3f', smallest_abs_diff, t_crit))

n = 13
m = 8
sp2 = 214.9^2

alpha = 0.05
t_crit = qt( 1-alpha, n+m-2 )
smallest_abs_diff = sqrt( sp2*(1/n + 1/m) ) * t_crit
print(sprintf('Part (b): |xbar-ybar|_min= %8.3f; t_crit= %8.3f', smallest_abs_diff, t_crit))


# Section 2; Question 12:
#
n = 10
xbar = 81.6
sx2 = 17.6^2
m = 20
ybar = 79.9
sy2 = 22.9^2

s2 = sx2/n + sy2/m

z = ( xbar - ybar )/sqrt(s2)
p_value = 2*(1-pnorm(abs(z)))

print(sprintf('z= %.6f; p_value= %.6f', z, p_value))


# Section 2; Question 13:
#
xbar = 33 # via interstate
sx2 = 6^2
ybar = 35 # through town
sy2 = 5^2

n = 2
m = ybar*n - xbar*n # the mean of the difference
s = sqrt(n*sx2 + n*sy2) # the standard deviation of the difference
print(sprintf('prob= %.6f', pnorm( -m/s )))


n = 10
m = ybar - xbar # the mean of the difference
s = sqrt(sx2/n + sy2/n) # the standard deviation of the difference
print(sprintf('prob= %.6f', pnorm( -m/s )))


# Section 2; Question 17:
#
severe = c( 640, 80, 1280, 160, 640, 640, 1280, 640, 160, 320, 160 )
asymtomatic = c( 10, 320, 320, 320, 320, 80, 160, 10, 640, 160, 320 )

n = length(severe)
xbar = mean(severe)
sx2 = var(severe)
m = length(asymtomatic)
ybar = mean(asymtomatic)
sy2 = var(asymtomatic)

# the variances are certainly not the same in this case:
#
theta_hat = sx2/sy2
numer = (theta_hat + n/m)^2
denom = theta_hat^2 / (n-1) + (n/m)^2 / (m-1)
nu = round( numer / denom )

w = ( xbar - ybar )/sqrt( sx2/n + sy2/m )
p_value = 1-pt(t, nu) # one-sided test

alpha = 0.05
t_crit = qt( 1-alpha, nu )
print(sprintf('w= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))


# Section 2; Question 18:
#
print(sprintf('nu= %d; n+m-2= %d', nu, n+m-2))


# Section 2; Question 19:
#
source('../Chapter8/chap_8_sect_2_question_2_data.R')

X = DF[ DF$type=='limited', ]
n = dim(X)[1]
xbar = mean(X$deaths)
sx2 = var(X$deaths)

Y = DF[ DF$type=='comprehensive', ]
m = dim(Y)[1]
ybar = mean(Y$deaths)
sy2 = var(Y$deaths)

print(sprintf('n= %d; xbar= %.3f; sx= %.3f; m= %d; ybar= %.3f; sy= %.3f', n, xbar, sqrt(sx2), m, ybar, sqrt(sy2)))


# Section 2; Question 20:
#
democrats = c( 22.4, 24.0, 38.0, 45.7, 21.2, 17.9, 38.2, 33.7, 23.8 )
republican = c( 45.7, 28.6, 14.2, 18.8, 50.3, 40.1, 52.4 )

n = length(democrats)
xbar = mean(democrats)
sx2 = var(democrats)
m = length(republican)
ybar = mean(republican)
sy2 = var(republican)

# the variances are certainly not the same in this case:
#
theta_hat = sx2/sy2
numer = (theta_hat + n/m)^2
denom = theta_hat^2 / (n-1) + (n/m)^2 / (m-1)
nu = round( numer / denom )

w = ( xbar - ybar )/sqrt( sx2/n + sy2/m )
p_value = 1-pt(t, nu) # one-sided test

alpha = 0.01
t_crit = qt( 1-alpha, nu )
print(sprintf('xbar= %.3f; sx= %.3f; ybar= %.3f; sy= %.3f', xbar, sqrt(sx2), ybar, sqrt(sy2)))
print(sprintf('w= %.6f; p_value= %.6f; t_crit= %.6f', t, p_value, t_crit))