#Copyright Marius Hofert

#remove previously used objects and load libs
remove(list=objects()) 

#initialise MT rng
set.seed(1)

#====user specifications====

#set path to working directory
mypath="/Users/mhofert/mhofert/latex/projects/job/2009\ cop/lecture/r"

#====exercise 1 (c)====

#create obscure.dat
myobscuredata=rchisq(n=1000,df=3)
obscurepath=paste(mypath,"/obscure.dat",sep="")
write.table(myobscuredata,obscurepath,row.names=F,col.names=F) #or use write.matrix() from the library MASS

#read data obscure.dat
obscure=read.table(obscurepath)

#apply the standard exponential distribution
mydata=dexp(obscure[,1]) #by exercise 1 (a) mydata should follow a U[0,1] distribution
plot(mydata) #does not look uniform! (can be tested with Anderson-Darling, Kolmogorov-Smirnov etc.)
plot(runif(1000)) #this are 1000 random variates from a U[0,1] distribution

#====exercise 1 (d)====

#(i)

#define the density of a double exponential distribution
f=function(x,mu,lambda){
	return(exp(-abs(x-mu)/lambda)/(2*lambda))
}

#define x values
x=seq(-10,10,length.out=5000)

#plot to pdf device
pdf(paste(mypath,"/doubleexpdensities.pdf",sep="")) #open pdf device for plotting
plot(x,f(x,0,1),type="l",col="1",ylab="f(x)",main="Double exponential densities") #first plot
lines(x,f(x,0,2),type="l",col="2") #add another plot
lines(x,f(x,1,1),type="l",col="3") #add another plot
legend(x="topleft",legend=c(expression(paste(mu,"=0, ",lambda,"=1",sep="")),expression(paste(mu,"=0, ",lambda,"=2",sep="")),expression(paste(mu,"=1, ",lambda,"=1",sep=""))),col=c(1,2,3),lty=c(1,1,1)) #add a legend
dev.off() #close plotting device

#(ii)

#define F^- (this actually coincides with F^-1 in this case)
generalizedF=function(y,mu,lambda){
	return(mu-lambda*sign(y-0.5)*log(1-2*abs(y-0.5)))
}

#generate U[0,1] random variates
u=runif(1000)

#apply F^-
x=generalizedF(u,1,2)

#plot data (just for fun)
plot(x)

#compute sample mean and sample variance
mean(x) #should be close to mu=1
var(x) #should be close to 2*lambda^2=8