You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
40 lines
1.3 KiB
R
40 lines
1.3 KiB
R
# Goals: Scare the hell out of children with the Cauchy distribution.
|
|
|
|
# A function which simulates N draws from one of two distributions,
|
|
# and returns the mean obtained thusly.
|
|
one.simulation <- function(N=100, distribution="normal") {
|
|
if (distribution == "normal") {
|
|
x <- rnorm(N)
|
|
} else {
|
|
x <- rcauchy(N)
|
|
}
|
|
mean(x)
|
|
}
|
|
|
|
k1 <- density(replicate(1000, one.simulation(20)))
|
|
k2 <- density(replicate(1000, one.simulation(20, distribution="cauchy")))
|
|
|
|
xrange <- range(k1$x, k2$x)
|
|
plot(k1$x, k1$y, xlim=xrange, type="l", xlab="Estimated value", ylab="")
|
|
grid()
|
|
lines(k2$x, k2$y, col="red")
|
|
abline(v=.5)
|
|
legend(x="topleft", bty="n",
|
|
lty=c(1,1),
|
|
col=c("black", "red"),
|
|
legend=c("Mean of Normal", "Mean of Cauchy"))
|
|
# The distribution of the mean of normals collapses into a point;
|
|
# that of the cauchy does not.
|
|
|
|
# Here's more scary stuff --
|
|
for (i in 1:10) {
|
|
cat("Sigma of distribution of 1000 draws from mean of normal - ",
|
|
sd(replicate(1000, one.simulation(20))), "\n")
|
|
}
|
|
for (i in 1:10) {
|
|
cat("Sigma of distribution of 1000 draws from mean of cauchy - ",
|
|
sd(replicate(1000, one.simulation(20, distribution="cauchy"))), "\n")
|
|
}
|
|
|
|
# Exercise for the reader: Compare the distribution of the median of
|
|
# the Normal against the distribution of the median of the Cauchy. |