Wednesday 29 July 2009

R stuff for normal distribution fitting

Look at this site for more: http://www.bigre.ulb.ac.be/Users/jvanheld/statistics_bioinformatics/practicals/microarray_fitting_solutions.html


   1:   
   2:  # data loading
   3:  filepath <- system.file("data", "morley.tab" , package="datasets")
   4:  mm <- read.table(filepath)
   5:  m <- mm[,1]
   6:   
   7:  hist(m)
   8:   
   9:  ## install
  10:  library(fBasics)
  11:   
  12:  skewness(m)
  13:   
  14:  kurtosis(m)
  15:   
  16:  plot(density(m))
  17:   
  18:  plot(ecdf(m))
  19:   
  20:  qqnorm(m)
  21:  abline(0,1)
  22:   
  23:  gal <- m
  24:   
  25:  ## Calculate estimators
  26:  m <- mean(gal,na.rm=T)
  27:  s <- sd(gal,na.rm=T)
  28:   
  29:  ## Draw the density histogram of the galactose microarray values
  30:  h <- hist(gal,breaks=100,col='#CCCCFF',border='#CCCCFF',freq=F)
  31:   
  32:  ## On the histogram, draw vertical bars at the following values :
  33:  ## mean, mean + 1*sd, mean -1*sd, mean +2*sd, mean -2*sd
  34:  abline(v=c(m,m-s,m-2*s,m+s,m+2*s),col="#000088",lwd=1)
  35:   
  36:  ## Superimpose the theoretical distribution
  37:  lines(h$mids,dnorm(h$mids,m,s), type="l", lwd=2,col="red")

No comments: