先日、"Using R for Introductory Econometrics"を勉強するゼミをやっていて、
教科書中にRスクリプトが公開されていない図（Figure 1.19とFigure 1.21）があったので、作ってみました。
たまに見返したくなるので、ここにあげておきます。
##スクリプト（Figure 1.19）
set.seed(123456)
ybar = numeric(10000)
for (i in c(10,50,100,1000)){
#"i" means the different simple sizes:
for (j in 1:10000){
sample = rnorm(i,10,2)
ybar[j] = mean(sample)
}
#change the name of each vector:
if (i == 10){
ybar10 = ybar
}
if (i == 50){
ybar50 = ybar
}
if (i == 100){
ybar100 = ybar
}
else {
ybar1000 = ybar
}
}
plot(density(ybar10), xlim=c(8.5,11.5), ylim=c(0,2), main="density of mean Y with different sample sizes", lty=1)
par(new=T)
plot(density(ybar50), xlim=c(8.5,11.5), ylim=c(0,2), main="", xlab="", lwd=1.5, lty=2)
par(new=T)
plot(density(ybar100), xlim=c(8.5,11.5), ylim=c(0,2), main="", xlab="", lwd=2, lty=3)
par(new=T)
plot(density(ybar1000), xlim=c(8.5,11.5), ylim=c(0,2), main="", xlab="", lwd=2.5, lty=4)
legend("topright", c("n=10", "n=50", "n=100", "n=1000"), lwd=c(1, 1.5, 2, 2.5), lty=1:4)
##結果（Figure 1.19）
###fig1.19.Rのスクリプトの説明
ここでは、正規分布（平均10、標準偏差2）に基づいてサンプルサイズの異なる（n=10~1,000）標本をそれぞれ10,000個生み出して平均を計算し、その密度分布を図にプロットしました。
なお、シード値を123456に固定しましたが、異なる値に変えると異なる結果が得られます。まあ、結局だいたい同じ図が描けると思いますが。
###Figure 1.19の説明
標本のサンプルサイズが大きくなればなるほど、標本平均の分散が0に近づくという"law of large numbers"をモンテカルロシミュレーションで示すために使われた図です（p61, p62）。
注意点としては、これは「無数にある標本の平均」の分布であり、標本内の数値の分布ではないということが挙げられます。
##スクリプト（Figure 1.21）
#definition of function, named "comp":
comp = function(x, i){
set.seed(i)
ybar = numeric(10000)
for (j in 1:10000){
sample = rnorm(x, 1, sqrt(2))
ybar[j] = mean(sample)
}
ybarN = ybar
for (j in 1:10000){
sample = rchisq(x, 1)
ybar[j] = mean(sample)
}
ybarC = ybar
#definition of the y limitation:
if (x<=10){
lims = 1
}
if (10<x&x<=100){
lims = 3
}
if (100<x&x<=1000){
lims = 10
}
if (1000<x&x<=10000){
lims = 30
}
#how to make the plot:
plot(density(ybarN), xlim=c(mean(ybar)-3*sd(ybar),mean(ybar)+3*sd(ybar)), ylim=c(0,lims), main="density of mean Y, norm and chi distributions", lwd=1, lty=1)
par(new=T)
plot(density(ybarC), xlim=c(mean(ybar)-3*sd(ybar),mean(ybar)+3*sd(ybar)), ylim=c(0,lims), main="", xlab="", lwd=2, lty=2)
legend("topright",c("normal", "chi2"), lwd=1:2, lty=1:2)
}
print("please type → comp(size, seed). if you do not have any ideas, please type → comp(100,12345). Then, a figure should appear.")
comp(100, 2134)
##結果（Figure 1.21）
###fig1.21のスクリプトの説明
fig1.19の目的は、同一の正規分布に従うサンプルサイズだけ異なる標本平均の分布を比べることだったので、1つの図の中に4つの分布を重ねる形で表示しました。
しかし、ここでは分布が異なる2つの標本の平均の分布を比べることが目的なので、1つの図の中にまとめるのはやめました。
その代わり、comp(size, random seed)という関数を新しく作って、好きなサンプルサイズとシード値で2つの分布を比べられるようにしました。
なお、2つの分布は、正規分布（平均10、標準偏差sqrt(2)）、カイ2乗分布（自由度1）です。
正規分布の標準偏差を上記のように設定した理由は、少し数式で説明する必要があります。以前、勉強したんですが、今は忘れちゃって説明できないので、気が向いたときに書き足しておきます。
###Figure 1.21の説明
ここでは、サンプルがどのような分布に従っていたとしても、サンプルサイズが大きくなればなるほど、サンプル平均の分布は正規分布に従うという"central limit theorem"をモンテカルロシミュレーションで示すために使われた図です（p61, p63）。
上に出力した図では、サンプルサイズを100に設定しましたが、確かに同じような分布になってます。
このように中心極限定理がどのような場合でも成り立つなら、「平均値が正規分布に従う」ことを仮定した定理などの信頼性が自分の中で増します。たしか、OLSにおける推定量の標準誤差を推定する際にもどこかで正規分布が仮定されていたような気がするんですが、、、
まあ、これもまた勉強した時に書き足しておきます。
以上です。ありがとうございました。
##参考文献
[1]Florian Heiss(2016) "Using R for Introductory Econometrics", Amazon Digital Services