mu=10
sigma=1
n=300
par(mfrow=c(1,1))
X=rcauchy(n,mu,sigma) # A sample of size n from N(10,1^2).
plot(X)
barX=mean(X)
cat("Based on this sample, sample mean barX is", barX, "\n",
    "where the center can be viewed as mu=", mu)

par(mfrow=c(2,2))
for(n in c(4, 16,25,100))
{
  M=1000
  SampleMean=c()
  for(m in 1:M)
  {
    X<- rcauchy(n,mu,sigma)
    SampleMean=c(SampleMean,mean(X))
  }
  plot(SampleMean)
}

par(mfrow=c(2,2))
par(mar=c(3,5,4,.1))
for(n in c(16,25,36,64))
{
  M=1000
  SampleMedian=c()
  for(m in 1:M)
  {
    X<- rcauchy(n,mu,sigma)
    SampleMean=c(SampleMean,mean(X))
  }
  hist(SampleMean,breaks=c(min(SampleMean)-1,quantile(SampleMean,prob=seq(0.05,0.95,by=0.05)),max(SampleMean)+1),freq=FALSE,xlim=c(mu-3*sigma,mu+3*sigma),ylim=c(0,3),ylab="Sampling distribution of the sample mean", xlab="")
  legend("topright",legend=paste("n=",n))
  lines(seq(mu-3*sigma,mu+3*sigma,length=100),dcauchy(seq(mu-3*sigma,mu+3*sigma,length=100),mu,sigma),col="blue")
  segments(mu,0,mu,3,col="red")
}

#### We should use median

n=300
par(mfrow=c(1,1))
X=rcauchy(n,mu,sigma) # A sample of size n from N(10,1^2).
plot(X)
medX=median(X)
cat("Based on this sample, sample median medX is", medX, "\n",
    "where the center can be viewed as mu=", mu)


par(mfrow=c(2,2))
par(mar=c(3,5,4,.1))
for(n in c(16,36,64,100))
{
  M=1000
  SampleMedian=c()
  for(m in 1:M)
  {
    X<- rcauchy(n,mu,sigma)
    SampleMedian=c(SampleMedian,median(X))
  }
  hist(SampleMedian,freq=FALSE,xlim=c(mu-3*sigma,mu+3*sigma),ylim=c(0,3),ylab="Sampling distribution of the sample median", xlab="")
  legend("topright",legend=paste("n=",n))
  lines(seq(mu-3*sigma,mu+3*sigma,length=100),dcauchy(seq(mu-3*sigma,mu+3*sigma,length=100),mu,sigma),col="blue")
  segments(mu,0,mu,3,col="red")
}