# This script plots data as a histogram
# the data are expected in a file called "indata" . The file need to be in the working directory 
# the first line of the data file should contain characters (no spaces) that give the title of the plot
# the following lines should conatin the observations, one observation per line 

# The script creates two plots 
# one plots the data as given
# plot two plots the normalized data and compares them to data from a normal distribution
# the second panel compares a sample from a normal distribution of the same size as the data

# note: SPACES ARE NOT ALLOWED IN THE NAME
mydata <- read.table("indata", header=TRUE)

######calculat normalized data
mydata_mean <- mean(mydata[[1]])
mydata_sigma <- sqrt(var((mydata[[1]])))
mydata_stand=(mydata-mydata_mean)/mydata_sigma
##################################################

#quartz()
#opens new graphics window on MAC
#x11()
#opens new graphics window on systems with X11 installed 
pdf(file="plot1.pdf")

#### Plot histogram of data #######
layout(matrix(1, 1, 1))
# hist (mydata[[1]],main = paste("Histogram of",colnames(mydata)))  #line n-1
hist (mydata[[1]],main = paste("Histogram of ", colnames(mydata)), breaks=25) #line n
#to change the number of histogram bars change the value for breaks in line n
#to use the program default, uncomment line n-1 and comment line n 
####################################

#quartz()
#opens new graphics window on MAC
#x11()
#opens new graphics window on systems with X11 installed 
pdf(file="plot2.pdf")########If you use the option to create a pdf file, you need to turn the device off
                        ######you can do this by uncommenting the last line

### Plot histogram of normalized data and compare to normal distibutio
### ALso assess the effect of sample size for aparent difference from normal distribution

#puts two plots on top of each other 
#first plot: real data normalized density diagram in red, green histogram of density from 10000 samples from a normal distribution
#second plot compares a sample taken from a normal distribution that has the same size as the original data, to a large sample 
layout(matrix(1:2, 2, 1))
hist (rnorm(100000),freq=F,breaks=100,border="green",,xlim=c(-5,5),ylim=c(-.1,0.6),main = paste("Histogram of", colnames(mydata)))
hist (mydata_stand[[1]],freq=F,breaks=25,add=T,border="red" )
hist (rnorm(100000),freq=F,breaks=100,border="green",,xlim=c(-5,5),ylim=c(-.1,0.6),main = paste("Histogram of rnorm with size of",colnames(mydata)))
hist (rnorm(length(mydata[[1]])),xlim=c(-5,5),ylim=c(-.1,0.6),freq=F,breaks=25,border="red", add=T )

### done
graphics.off()