Explore Confidence Interval, Difference of paired values


The script below provides a way to
  1. Create one population of values with specified means and standard deviations.
  2. Create a parallel population of values with a specified mean and standard deviation to hold the change from the first population to a new population of paired values.
  3. Create that second population of paired values.
  4. Confirm the mean and standard deviation of each population.
  5. Specify the size of a sample to be taken from each population of paired values.
  6. Specify the significance level to use to make confidence intervals.
  7. Specify the number of times to take such samples.
  8. Perform the sampling and, for each sample find the mean difference and the standard deviation of those differences within the pairs of values in the sample.
  9. Based on those mean and standard deviation values, generate confidence intervals.
  10. Keep track of the number of times that the generated confidence interval contains the true mean difference.
  11. Report a dataframe of the very last of the samples taken.
  12. Report the count of the number of generated confidence intervals that contain the true difference in the population.
If we ask for a significant number of samples, say 10,000, we can see that we really do generate about the right number of confidence intervals that contain the true mean.

In the folder containing the function scripts for this course create a new directory, copy the model.R file to that directory, rename the file in the new directory, double click on the file to open Rstudio. Then copy all of the text below the line and paste it into your Rstudio editor pane. Then, you can highlight the entire script and run it to use the default values. After that you can go back and change parameters and run the script again to explore the consequences of those changes. 
  
# We look at a confidence interval for the
# difference of two population means when we 
# are looking at paired values

#  first set up some goal populations
pop_one_mean <- 23.2
pop_one_sd <- 5.3
#  now we want to set up the differences
pop_diff_mean <- 4.3
pop_diff_sd <- 6.2

#  then create the first values as a  population,
#  with 1000 values

#  First generate an approximate standard normal
pop_one <- rnorm( 1000 )

#  Then get its mean and standard deviation
mu_1 <- mean( pop_one)
#  we want the standard deviation of the population
source("../pop_sd.R")
sd_1 <- pop_sd( pop_one )

#  Then create the distribution we want
pop_one <- ( (pop_one-mu_1)/sd_1 )* pop_one_sd+pop_one_mean

#   finally, verify that we have the right population
mean( pop_one )
pop_sd( pop_one )

#  Now do the same thing for pop_diff

#  First generate an approximate standard normal
pop_diff <- rnorm( 1000 )

#  Then get its mean and standard deviation
mu_2 <- mean( pop_diff )

sd_2 <- pop_sd( pop_diff )

#  Then create the distribution we want
pop_diff <- ( (pop_diff-mu_2)/sd_2 ) *
                pop_diff_sd+pop_diff_mean

#   finally, verify that we have the right population
mean( pop_diff )
pop_sd( pop_diff )

#   Now we can create the second values for each
#   each value in the original population

pop_two <- pop_one + pop_diff



#   Now we want to repeat the following process of
#   getting a sample of the pairs of values,  Then,
#   for each pair we can get the difference, in our
#   case we will look at the second value minus the
#   first value (to mimic looking at growth).
#   Then generating the confidence interval
#   for the mean difference of the values.

#   While we are at this, and because we know what that
#   difference is, we can count the number of times 
#   that the true difference is inside our interval.

sig_level <- 0.93
samp_size <- 23

num_reps <- 100
num_success <- 0
num_fail <- 0
true_diff <- pop_diff_mean


#  from the confidence level we want to find the
#  area in the two tails
alpha <- 1 - sig_level
alpha_div_2 <- alpha/2

# We know that we will be using the Student's-t
# distribution.  Therefore, we need to find the degrees
# of freedom.  But this is just one less than the 
# sample size.  Therefore, our t-value will be the
# same for each sample.

t_val <- qt( alpha_div_2, samp_size-1, lower.tail=FALSE)



for( i in (1:num_reps) )
{  # choose samples from pop one get sample mean and sd
   index_1 <- as.integer( runif( samp_size, 1, 1001))
   samp_1 <- pop_one[ index_1 ]
   samp_2 <- pop_two[ index_1 ]
   #  
   diff_values <- samp_2 - samp_1

   xbar_diff <- mean( diff_values )
   sd_diff <- sd( diff_values )
   s_e <- sd_diff/sqrt( samp_size )
   
   # get the confidence interval
   ci_low <- xbar_diff - t_val*s_e
   ci_high <- xbar_diff + t_val*s_e
   in_ci <- (ci_low <= true_diff ) &&
            ( true_diff <= ci_high )
   if( in_ci )
   { num_success <- num_success+1} else
  
   { num_fail <- num_fail + 1}
}

#  look at the last sample as a data frame
data.frame( samp_1, samp_2, diff_values )
#  report the number of successes
num_success