##Simulation: generate 10000 points from N(5,2) and plot the density
data<-rnorm(10000,5,2)
plot(density(data),col="blue", lwd=2,ylim=c(0,0.5),xlim=c(-15,15),main="")

##Simulation: each value in previous data subtracts 5, i.e., y-mu
data1<-data-5
par(new = TRUE)
plot(density(data1),col="black", lwd=2,ylim=c(0,0.5),xlim=c(-15,15),main="")

##Simulation: each value in previous data subtracts 5 and then divided by 2, i.e., (y-mu)/sigma
data2<-(data-5)/2
par(new = TRUE)
plot(density(data2),col="red", lwd=2,ylim=c(0,0.5),xlim=c(-15,15),main="")

##Compare with standard normal curve, i.e., mu=0, std=1 
standard<-function(x) dnorm(x)
curve(standard, col="green", lwd=2,from=-15,to=15,add=T,main="")