stat_data_analysis_rcode

Statistical Data Analysis - Fall 2018
Log | Files | Refs | README

R-IC05.R (3824B)


      1 #*******************************************************#
      2 #*******************************************************#
      3 # 		         STAT320 - IC05 Solutions	      	  	    #
      4 #*******************************************************#
      5 #*******************************************************#
      6 #***************************************************#
      7 #           PVal/Power for Means                    #  
      8 #***************************************************#
      9 
     10 #***************************************************#
     11 #***************************************************#
     12 #***************************************************#
     13 #***************************************************#
     14 #***************************************************#
     15 
     16 #***************************************************#
     17 # Means P-value
     18 #***************************************************#
     19 set.seed(257)
     20 #***************************************************#
     21 
     22 xsamp1=c(96.7, 108.8, 95.1, 110, 104.5, 106.7, 120.5, 113.4, 94.4, 67.6, 81.9, 118.1)
     23 mean(xsamp1)
     24 NSim1=100000
     25 Low1=FALSE
     26 NSize1=12
     27 
     28 mean.pval.f = function(xobs,xsamp,mu0=0,NSize=100,Low = TRUE,NSim=100000)
     29 {
     30   # Rescale sample so as to have mean equal to mu0
     31   rescaledsample = xsamp - mean(xsamp) + mu0
     32   
     33   # Create object where simulated averages will be stored
     34   OXBars = numeric(NSim)
     35   
     36   # Repeat simulations NSim times
     37   for(i in 1:NSim)
     38   {
     39     # Generate sample by resampling original sample WITH replacement 
     40     Sample = sample(x=rescaledsample,size = NSize,replace = TRUE) 
     41     # Save mean of simulated sample
     42     OXBars[i] = mean(Sample)
     43   }
     44   # Compute p-value
     45   if(Low==TRUE)
     46   {
     47     # if alternative Ha is mu < mu0
     48     p = sum(OXBars<=xobs)/NSim
     49   }else
     50   {
     51     # if alternative Ha is mu > mu0
     52     p = sum(OXBars>=xobs)/NSim
     53   }
     54   return(p) 
     55 }
     56 pval = mean.pval.f(xobs = mean(xsamp1), xsamp=xsamp1,mu0=98,NSize=NSize1,Low=Low1,NSim=NSim1)
     57 pval
     58 
     59 mu0=10
     60 
     61 #***************************************************#
     62 # Means Cutoff
     63 #***************************************************#
     64 
     65 mean.pval.cutoff.f = function(a=.05,xsamp,mu0=0,NSize=100,Low = TRUE,NSim=100000)
     66 {
     67   # Rescale sample so as to have mean equal to mu0
     68   rescaledsample = xsamp - mean(xsamp) + mu0
     69   
     70   # Create object where simulated averages will be stored
     71   OXBars = numeric(NSim)
     72   
     73   # Repeat simulations NSim times
     74   for(i in 1:NSim)
     75   {
     76     # Generate sample by resampling original sample WITH replacement 
     77     Sample = sample(x=rescaledsample,size = NSize,replace = TRUE) 
     78     # Save mean of simulated sample
     79     OXBars[i] = mean(Sample)
     80   }
     81   # Compute p-value
     82   if(Low==TRUE)
     83   {
     84     # if alternative Ha is mu < mu0
     85     x = quantile(OXBars,probs = a)
     86   }else
     87   {
     88     # if alternative Ha is mu > mu0
     89     x = quantile(OXBars,probs = 1-a)
     90   }
     91   return(x) 
     92 }
     93 cOff = mean.pval.cutoff.f(a = 0.01,xsamp=xsamp1,mu0=29.5,NSize=18,Low=TRUE,NSim=100000)
     94 cOff
     95 
     96 
     97 #***************************************************#
     98 # Means CI
     99 #***************************************************#
    100 
    101 mean.CI.f = function(xsamp,a=.05,NSize=100,NSim=100000)
    102 {
    103   # Rescale sample so as to have mean equal to mu0
    104   rescaledsample = xsamp - mean(xsamp) + mu0
    105   
    106   # Create object where simulated averages will be stored
    107   OXBars = numeric(NSim)
    108   
    109   # Repeat simulations NSim times
    110   for(i in 1:NSim)
    111   {
    112     # Generate sample by resampling original sample WITH replacement 
    113     Sample = sample(x=rescaledsample,size = NSize,replace = TRUE) 
    114     # Save mean of simulated sample
    115     OXBars[i] = mean(Sample)
    116   }
    117   x = quantile(OXBars,probs = c(a/2,1-a/2))
    118   return(x) 
    119 }
    120 
    121 CI = mean.CI.f(xsamp=xsamp1, a=0.05, NSize=12, NSim=100000)
    122 CI
    123 
    124 pow = mean.pval.f(xobs = mean(xsamp1) ,xsamp=xsamp1,mu0=99,NSize=12,Low=FALSE,NSim=100000)
    125 pow
    126