population.size = 1000000
population.mean1 = 300
population.std1 = 50
population.mean2 = 500
population.std2 = 50

population.data1 = rnorm(population.size,population.mean1,population.std1)
population.data2 = rnorm(population.size,population.mean2,population.std2)
population.mean1 = mean(population.data1)
population.mean2 = mean(population.data2)

par(mfrow=c(2, 1))

# plot histogram with normal curve
par(mar=c(1,1,1,1))
hist(population.data1, col=rgb(1,0,0,0.5), breaks=20, xlab="Reaction Time",main="",xlim = c(1,800),ylim=c(0,300000)) 
hist(population.data2, col=rgb(0,0,1,0.5), breaks=20, xlab="Reaction Time",main="",xlim = c(1,800),ylim=c(0,300000),add=T)

abline(a = NULL, b = NULL, h = NULL, v = (population.mean1),col=2)
abline(a = NULL, b = NULL, h = NULL, v = (population.mean2),col=4)

my.data = 0
my.data[1] = population.mean1
my.data[2] = population.mean2

barplot(my.data)