[R] R project

From: Nicola Smith <nicola_1988_at_hotmail.com>
Date: Sun, 13 Apr 2008 16:30:13 +0000

Hi,I am currently doing a project in which we are to investigate the size and power of three different one sample tests over three different distributions using a number of different sample sizes and values for mu1. I have written a function and was trying to get my answer for each test into the right position in an array so the output is the power of each combination of test, distribution, sample size and mu1, however my code is not working and i keep getting an error message returned. My code is: proj1.sim<-function(){
  tests<-c("t","Wilcoxon","bootstrap")
  distns<-c("Normal","Uniform","Beta")
  sample.sizes<-c("5","10","25","50","100")   mu1s<-c("0","0.5","1","1.5","2")
  res<-array(0,c(length(tests),length(distns),length(sample.sizes),length(mu1s)))   dimnames(res)<-list(tests,distns,sample.sizes,mu1s)   for(test in tests){
    for(distn in distns){

      for(sample.size in sample.sizes){
        for(mu1 in mu1s){
          result<-hypoth.test(sample.size,distn,test,mu1)
          res[test,distn,sample.size,mu1]<-result
  }}}}
 output(res)
}  

hypoth.test<-function(sample.size,distn,test,mu1){
#Purpose:
#Samples n values from the chosen distribution
#Inputs:
#n - sample size of data to generate
#distn - the distribution to generate data from
#test - which test to use
#H0 - the null hypothesis
#mu.true - the true value of the mean under the chosen distribution
#sigma.true - the true value of the variance under the chosen distribution
#alpha - trimming level
#K - number of simulations to perform
#Outputs:
#reject.null - the number of times we reject the null hypothesis
 

  H0<-0
  alpha<-0.05
  K<-1100
  sigma<-1
  #create vectors for the results
  reject<-numeric(K)
  mu1<-numeric(5)
  sample.size<-numeric(5)
  #loop through the observations the required number of K simulations   for(i in 1:K){
    #determine whether the data comes from a Normal or Uniform distribution     if(distn=="Normal"){

      #Normal distribution
      dat<-rnorm(n=sample.size,mean=mu1,sd=sigma)
    } else if(distn=="Uniform"){
      dat<-runif(n=sample.size,min=-sqrt(3)+mu1, max=sqrt(3)+mu1)
    } else {
      z<-rbeta(n=sample.size,shape1=0.8,shape2=7.2)
      dat<-((z-0.1)*(sigma/0.1))+H0

    }
    #determine which test is required to sample the observations     if(test=="t"){
      #t-test
      res<-t.test(dat,alternative="two.sided",mu=H0)
    } else if(test=="Wilcoxon"){
      res<-wilcox.test(dat,alternative="two.sided",mu=H0)
    } else {
      nboot<-(K-1)
      n<-length(dat)
      bootmean<-NULL
      for(i in 1:nboot){
        rdata<-sample(dat,n,replace=T)
        bootmean[i]<-mean(rdata)
      }
      nlo<-round((nboot+1)*(alpha/2))
      nhi<-round((nboot+1)*((1-alpha)/2))
      bootmean<-sort(bootmean)
      low<-bootmean[nlo]
      high<-bootmean[nhi]
      res<-c(low,high)

    }
    #record whether we reject the null or not i.e. whether the p-value is less     #than or equal to the given alpha level(reject), or not (do not reject)     if(test=="bootstrap"){
      reject[i]<-(H0<low|H0>high)
    } else {
      reject[i]<-(res$p.value<=alpha)
    }
    #sum together the total number of null hypothesis rejected     reject.null<-sum(reject)
  }
  power<-(reject.null)/K
  return(power)
} Would you be able to tell me where I am going wrong? Thanks,Nicola Smith

[[elided Hotmail spam]]

        [[alternative HTML version deleted]]



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Sun 13 Apr 2008 - 17:36:25 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sun 13 Apr 2008 - 18:30:29 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive