<h1>Pythonで統計学を学ぶ(3)</h1>
この内容は<A HREF="http://shop.ohmsha.co.jp/shop/shopdetail.html?brandcode=000000001781&search=978-4-274-06710-5&sort=" target="_blank">山田、杉澤、村井(2008)「R」によるやさしい統計学</A>を参考にしています。
<P>
この講義では、「母集団と標本」をとりあげます。
これは、大きな集団から一部を取り出した少数のデータの情報を用いて、もとの集団の性質を推測する推測統計の基本的な理論の学習がねらいです。
特に 標本分布の理解を目的とします。
<p>学習項目です</p>
<ul class="important">
<li class="important"><A HREF="#RS03:popularity">母集団と標本の関係</A></li>
<li class="important"><A HREF="#RS03:estimationt">推測統計</A></li>
<li class="important"><A HREF="#RS03:pointEstimation">点推定</A></li>
<li class="important"><A HREF="#RS03:Basics">推定の基礎</A></li>
<li class="important"><A HREF="#RS03:sampleDistribution">標本分布</A></li>
<LI class="important"><A HREF="#Functions">関数のまとめ</A>
<LI class="important"><A HREF="#RS03:Exercises">演習問題</A>
</ul>

Pythonで統計学を学ぶ(3)

この内容は山田、杉澤、村井(2008)「R」によるやさしい統計学を参考にしています。

この講義では、「母集団と標本」をとりあげます。 これは、大きな集団から一部を取り出した少数のデータの情報を用いて、もとの集団の性質を推測する推測統計の基本的な理論の学習がねらいです。 特に 標本分布の理解を目的とします。

学習項目です

 
<H2><A NAME="RS03:popularity">母集団と標本の関係</A></H2>
<div class="content">
まず重要な用語をまとめておきましょう:
<ul class="important">
  <LI><B>母集団</B>:  対象のデータ全体               
  <LI><B>標本</B>: 母集団から一部のデータを取り出した(抽出した)もの                
  <LI><font class="word"><B>標本抽出</B></font>: 母集団から標本をとりだすこと
  <LI><font class="word"><B>母数</B></font>: 母集団の性質をあらわす統計的指数(比率、平均、分散、相関係数など)
</UL>
例えば、日本の中学生全体のテスト得点データを母集団とすると、名古屋市昭和区の(一部の)中学生のテスト得点データが「標本」になります。そして、平均点や分散などが母数とみなせます。

母集団と標本の関係

まず重要な用語をまとめておきましょう:

  • 母集団: 対象のデータ全体               
  • 標本: 母集団から一部のデータを取り出した(抽出した)もの                
  • 標本抽出: 母集団から標本をとりだすこと
  • 母数: 母集団の性質をあらわす統計的指数(比率、平均、分散、相関係数など)
例えば、日本の中学生全体のテスト得点データを母集団とすると、名古屋市昭和区の(一部の)中学生のテスト得点データが「標本」になります。そして、平均点や分散などが母数とみなせます。

 
<H2><A NAME="RS03:estimation">推測統計</A></H2>
<P>
手元にある標本から母集団の性質(母数)を推測することを推測統計といいます。
推測統計は<font class="word"><B>推定</B></font><font class="word"><B>検定</B></font>に分けられます。
<P>
「推定」は、母数の値に対し具体的な値を求める<font class="word"><B>点推定</B></font>、つまり一つの値で推定の結果を表すものと、<font class="word"><B>区間推定</B></font>、つまりある幅を持った区間で結果を表すものに分けられます。
<P>
また「検定」は母集団について述べた異なる立場の二つの主張(<font class="word">仮説</font>)のうちどれを採用するかを決定するものです。例えば、日本の中学生のテスト得点が5年前と現在とで変化したという仮説と、変化していないという仮説がある場合、どちらを採択するかをきめるのが検定です。
</div>

推測統計

手元にある標本から母集団の性質(母数)を推測することを推測統計といいます。 推測統計は推定検定に分けられます。

「推定」は、母数の値に対し具体的な値を求める点推定、つまり一つの値で推定の結果を表すものと、区間推定、つまりある幅を持った区間で結果を表すものに分けられます。

また「検定」は母集団について述べた異なる立場の二つの主張(仮説)のうちどれを採用するかを決定するものです。例えば、日本の中学生のテスト得点が5年前と現在とで変化したという仮説と、変化していないという仮説がある場合、どちらを採択するかをきめるのが検定です。

 
<H2><A NAME="RS03:pointEstimationt">点推定</A></H2>
<div class="content">
母集団から抽出した標本を用いて、母数の点推定を行う手順についてみていきます。
ここで、標本のデータの個数のことを<font class="word"><B>サンプルサイズ</B></font>
もしくは<font class="word"><B>標本の大きさ</B></font>と言います。
<P>
例えば、10人の17歳男性の身長データがあるとします。これは日本人の17歳男性全体という『母集団』から抽出された『標本』とみなせます。そしてこのとき標本のデータの個数が10ですから、<font class="rule">n = 10</font>と書きます。ここで<font class="word">n</font>はサンプルサイズを表す記号です。

点推定

母集団から抽出した標本を用いて、母数の点推定を行う手順についてみていきます。 ここで、標本のデータの個数のことをサンプルサイズ、 もしくは標本の大きさと言います。

例えば、10人の17歳男性の身長データがあるとします。これは日本人の17歳男性全体という『母集団』から抽出された『標本』とみなせます。そしてこのとき標本のデータの個数が10ですから、n = 10と書きます。ここでnはサンプルサイズを表す記号です。

In [2]:
x
 
import numpy as np
Height=np.array([165.2, 175.9, 161.7, 174.2, 172.1, 163.3, 170.6, 168.4, 171.3])
 
この標本から、「17歳の日本人男性全体の身長」(母集団)の平均(母数の一つ)を点推定してみましょう。
それには、次のようにします。ここで<TT>mean</TT>は、平均値を求めるRの関数でしたね。

この標本から、「17歳の日本人男性全体の身長」(母集団)の平均(母数の一つ)を点推定してみましょう。 それには、次のようにします。ここでmeanは、平均値を求めるRの関数でしたね。

In [3]:
 
print(Height.mean())
169.188888889
<H2><A NAME="RS03:Basics">推定の基礎</A></H2>
<div class="content">
まず用語を確認しましょう:
<UL class="important">
  <LI><font class="word">標本統計量</font>:標本データから計算される平均や分散などの値
  <LI><font class="word">推定値</font>:標本データを用いて計算した母数の推定量
</UL>
母集団、標本、標本統計量、母数の関係を図で表すと以下のようになります: 母集団から抽出された標本データを分析し、標本統計量を得ます、それを使って母数を推定するのです
<CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-1.png" target="_blank">
<IMG class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-1.png" border="0"  width="50%"  hspace="10" /> 
</A>
</CENTER>
ここで、標本統計量からどのように母数を推定するかを表にしておきましょう:
|母数 | 推定量| 関係する関数|オプション|
|:---|:---|:---|:---|
|母平均|標本平均|np.mean| |
|母分散|不偏分散|np.var|ddof=1|
|母標準偏差|不偏分散の正の平方根|np.std|ddof=1|
|母相関係数|標本相関係数|np.corrcoef|対応する要素 |

推定の基礎

まず用語を確認しましょう:

  • 標本統計量:標本データから計算される平均や分散などの値
  • 推定値:標本データを用いて計算した母数の推定量
母集団、標本、標本統計量、母数の関係を図で表すと以下のようになります: 母集団から抽出された標本データを分析し、標本統計量を得ます、それを使って母数を推定するのです

ここで、標本統計量からどのように母数を推定するかを表にしておきましょう:

母数 推定量 関係する関数 オプション
母平均 標本平均 np.mean
母分散 不偏分散 np.var ddof=1
母標準偏差 不偏分散の正の平方根 np.std ddof=1
母相関係数 標本相関係数 np.corrcoef 対応する要素

In [9]:
 
#
 
<h4>課題3-1</h4>
前項で例に出した「身長データ」を用いて、母集団の分散と標準偏差を点推定した値を求めよ。
またそれに用いたコードをあわせて答えよ。

課題3-1

前項で例に出した「身長データ」を用いて、母集団の分散と標準偏差を点推定した値を求めよ。 またそれに用いたコードをあわせて答えよ。

In [ ]:
 
# kadai 3-1
 
<h3>標本抽出に伴う誤差</h3>
<div class="content">
標本から推定値を得る手順よりも大事なこと:
 <UL>
     <LI>実際の母数の値にどの位近い推定値が得られるか?
     <LI>推定の結果はどのくらい信用できるか?
</UL>
<P>
<UL class="important">
    <LI><font class="word"><B>標本誤差</B></font>: 母数の値と推定値とのズレ---
データを沢山集めれば、生じうる誤差は小さくできる
    <LI><font class="word"><B>標本分布</B></font>: どのような推定値が得られる可能性があるか、についての手がかり
</UL>

標本抽出に伴う誤差

標本から推定値を得る手順よりも大事なこと:
  • 実際の母数の値にどの位近い推定値が得られるか?
  • 推定の結果はどのくらい信用できるか?

  • 標本誤差: 母数の値と推定値とのズレ--- データを沢山集めれば、生じうる誤差は小さくできる
  • 標本分布: どのような推定値が得られる可能性があるか、についての手がかり

 
<h3>推定値の信頼性を調べる方法</h3>
<div class="content">
推定値を計算する元になる、標本に含まれる「個々のデータの値」の決定方法について見ていきましょう。
<OL class="important">
    <LI>標本抽出の方法---
単純無作為抽出と呼ばれる方法が前提
  <UL>
    <LI><font class="word">単純無作為抽出</font>:  母集団の中のどのデータにも平等に選ばれる可能性を持っているような標本抽出の方法
    <LI><font class="word">無作為標本</font>: 単純無作為抽出によって得られた標本
   </UL>
    <LI>単純無作為抽出によって得られるデータの性質としての確率変数
    <UL>
        <LI><font class="word">確率変数</font>: 実際に結果が得られるまでどのような値が得られるかが決まっていない変数<BR>
例: サイコロを振った時の出目<BR>
同じ手続でデータを得ても結果に再現性がない⇒
単純無作為抽出によって得られるデータは確率変数といえる
    </UL>
    <LI>確率変数がどのような値をとるのかを示す確率分布
    <UL>
        <LI><font class="word">確率変数の実現値</font>:抽出の結果、確率変数がとる値
        <LI><font class="word">確率分布</font>: 確率変数がどのような値をどのような確率でとるかということを表した分布
            <LI>確率分布を用いた母集団の表現としての母集団分布
    <UL>
        <LI><font class="word">母集団分布</font>:ある変数の母集団における分布<BR>
無作為抽出の場合、標本データの確率分布は母集団分布と同じになります
    </ul>
    <LI>代表的な母集団分布である<font class="word">正規分布</font>
    <LI>Pythonを使って正規分布の母集団から標本を抽出する方法
</OL>

推定値の信頼性を調べる方法

推定値を計算する元になる、標本に含まれる「個々のデータの値」の決定方法について見ていきましょう。
  1. 標本抽出の方法--- 単純無作為抽出と呼ばれる方法が前提
    • 単純無作為抽出: 母集団の中のどのデータにも平等に選ばれる可能性を持っているような標本抽出の方法
    • 無作為標本: 単純無作為抽出によって得られた標本
  2. 単純無作為抽出によって得られるデータの性質としての確率変数
    • 確率変数: 実際に結果が得られるまでどのような値が得られるかが決まっていない変数
      例: サイコロを振った時の出目
      同じ手続でデータを得ても結果に再現性がない⇒ 単純無作為抽出によって得られるデータは確率変数といえる
  3. 確率変数がどのような値をとるのかを示す確率分布
    • 確率変数の実現値:抽出の結果、確率変数がとる値
    • 確率分布: 確率変数がどのような値をどのような確率でとるかということを表した分布
    • 確率分布を用いた母集団の表現としての母集団分布
      • 母集団分布:ある変数の母集団における分布
        無作為抽出の場合、標本データの確率分布は母集団分布と同じになります
    • 代表的な母集団分布である正規分布
    • Pythonを使って正規分布の母集団から標本を抽出する方法

 
[参考: サイコロの出目の確率分布]
理論的には、サイコロの出目の確率分布はつぎのようになります:
<table  border=2 width="700">
 <tr>
    <th width="150">サイコロの出る目</th>
    <th align="center">1</th>
    <th align="center">2</th>
    <th align="center">3</th>
    <th align="center">4</th>
    <th align="center">5</th>
    <th align="center">6</th>
 </tr>
 <tr>
    <th>確率</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
 </tr>
</table>

[参考: サイコロの出目の確率分布] 理論的には、サイコロの出目の確率分布はつぎのようになります:

サイコロの出る目 1 2 3 4 5 6
確率 1/6 1/6 1/6 1/6 1/6 1/6
 
それに対し、シミュレーションにより、サイコロの出る目の度数分布は次のようになりました:
(1) 6回サイコロを振った場合:
<table  border=2 width="700">
 <tr>
    <th width="160">サイコロの目</th>
    <th align="center" width="90">1</th>
    <th align="center" width="90">2</th>
    <th align="center" width="90">3</th>
    <th align="center" width="90">4</th>
    <th align="center" width="90">5</th>
    <th align="center" width="90">6</th>
 </tr>
 <tr>
    <th>度数</th>
    <th align="center">0</th>
    <th align="center">2</th>
    <th align="center">1</th>
    <th align="center">1</th>
    <th align="center">0</th>
    <th align="center">2</th>
 </tr>
 <tr>
    <th>相対度数</th>
    <th align="center">0</th>
    <th align="center">2/6</th>
    <th align="center">1/6</th>
    <th align="center">1/6</th>
    <th align="center">0</th>
    <th align="center">2/6</th>
 </tr>
</table><BR>
(2) 60,000回サイコロを振った場合:
<table  border=2 width="700">
 <tr>
    <th width="160">サイコロの目</th>
    <th align="center" width="90">1</th>
    <th align="center" width="90">2</th>
    <th align="center" width="90">3</th>
    <th align="center" width="90">4</th>
    <th align="center" width="90">5</th>
    <th align="center" width="90">6</th>
 </tr>
 <tr>
    <th>度数</th>
    <th align="center">10,097</th>
    <th align="center">9,918</th>
    <th align="center">9,936</th>
    <th align="center">9,960</th>
    <th align="center">10,050</th>
    <th align="center">10,039</th>
 </tr>
 <tr>
    <th>相対度数</th>
    <th align="center">1.0097</th>
    <th align="center">0.9918</th>
    <th align="center">0.9936</th>
    <th align="center">0.9960</th>
    <th align="center">1.0050</th>
    <th align="center">1.0039</th>
 </tr>
</table>
このように、非常に多くの実現値が得られれば、度数分布は確率分布とほぼ同じになります

それに対し、シミュレーションにより、サイコロの出る目の度数分布は次のようになりました: (1) 6回サイコロを振った場合:

サイコロの目 1 2 3 4 5 6
度数 0 2 1 1 0 2
相対度数 0 2/6 1/6 1/6 0 2/6

(2) 60,000回サイコロを振った場合:

サイコロの目 1 2 3 4 5 6
度数 10,097 9,918 9,936 9,960 10,050 10,039
相対度数 1.0097 0.9918 0.9936 0.9960 1.0050 1.0039
このように、非常に多くの実現値が得られれば、度数分布は確率分布とほぼ同じになります

 
<H3><A NAME="RS03:normalDistribution">正規分布</A></H3>
<div class="content">
正規分布(normal distribution, ガウス分布 Gaussian distributionともいう)は統計学で最もよく用いられる確率分布なので、ここで復習しておこう<P>
正規分布の表現方法: <I>N</I> (μ,σ<sup>2</sup>) --- μは平均値、σ<sup>2</sup>は分散 (σは標準偏差)
<CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-4.png" target="_blank">
<IMG  class="exp" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-4.png" width="40%" hspace="10" />
</A>
</CENTER>
<TABLE >
<TR><TD width="50%">
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-2.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-2.png" border="0" width="60%" hspace="10">
</A></TD>
<TD>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-3.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-3.png" border="0"  width="60%" hspace="10"> 
</A> </TD>
</TR>
<TR><TD  align="center">正規分布 <I>N</I> (0,1)のグラフ</TD>
<TD  align="center">正規分布 <I>N</I> (1,1) のグラフを重ね書き</TD></TR>
</TABLE>
<P>
<ul class="important">
  <LI><font class="word">標準正規分布</font>: 平均0 分散1の正規分布 --- N(0,1)
  <LI>離散変数: 整数など、とびとびの値をとる変数(例:サイコロの出目)
  <LI>連続変数: 実数値をとる変数(例:正規分布に従う確率変数)
  <LI><font class="word">確率密度</font> : 正規分布の図の縦軸は確率密度<BR>
    <font class="rule">確率密度×確率変数の範囲=確率</font>
  <LI><font class="word">確率密度関数</font>:確率密度を確率変数の値の関数として表したもの 
</ul>
正規分布の性質: <UL>
   <LI>μ-σからμ+σ の範囲(偏差値で言えば、40点から60点の範囲)の確率はおよそ 0.683 
   <LI>μ-2σからμ+2σ の範囲(偏差値で言えば、30点から70点の範囲)の確率はおよそ 0.954 
   <LI>μ-3σからμ+3σ の範囲(偏差値で言えば、20点から80点の範囲)の確率はおよそ 0.997 
</UL>
<TABLE >
<TR><TD width="33%">
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-5.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-5.png" border="0" width="80%" hspace="10">
</A></TD>
<TD width="34%">
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-6.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-6.png" border="0"  width="80%" hspace="10">
</A>  </TD>
 <TD width="33%">
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-7.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-7.png" border="0"  width="80%"  hspace="10">
</A> </TD></TR>
<TR><TD width="33%" align="center">μ-σからμ+σ の範囲</TD><TD align="center">μ-2σからμ+2σ の範囲</TD><TD width="33%" align="center">μ-3σからμ+3σの範囲 </TD></TR>
</TABLE>

正規分布

正規分布(normal distribution, ガウス分布 Gaussian distributionともいう)は統計学で最もよく用いられる確率分布なので、ここで復習しておこう

正規分布の表現方法: N (μ,σ2) --- μは平均値、σ2は分散 (σは標準偏差)

正規分布 N (0,1)のグラフ 正規分布 N (1,1) のグラフを重ね書き

  • 標準正規分布: 平均0 分散1の正規分布 --- N(0,1)
  • 離散変数: 整数など、とびとびの値をとる変数(例:サイコロの出目)
  • 連続変数: 実数値をとる変数(例:正規分布に従う確率変数)
  • 確率密度 : 正規分布の図の縦軸は確率密度
       確率密度×確率変数の範囲=確率
  • 確率密度関数:確率密度を確率変数の値の関数として表したもの 
正規分布の性質:
  • μ-σからμ+σ の範囲(偏差値で言えば、40点から60点の範囲)の確率はおよそ 0.683
  • μ-2σからμ+2σ の範囲(偏差値で言えば、30点から70点の範囲)の確率はおよそ 0.954
  • μ-3σからμ+3σ の範囲(偏差値で言えば、20点から80点の範囲)の確率はおよそ 0.997
μ-σからμ+σ の範囲μ-2σからμ+2σ の範囲μ-3σからμ+3σの範囲

In [30]:
 
from scipy import stats as st
import matplotlib.pyplot as plt
x = np.linspace(-4, 4, 10000)
plt.plot(x, st.norm.pdf(x,0,1), label=r'N(%i,%i)$' % (0,1))
plt.axis([-4,4,-0.02,0.5]) 
plt.axhline(0,c='g',ls='-.')
plt.axvline(-2,c='b',ls='--')
plt.axvline(2,c='b',ls='--')
x = np.linspace(-2, 2, 10000)
plt.fill_between(x, stats.norm.pdf(x,0,1))
plt.show()
 
<H4>正規母集団から単純無作為抽出を行う</H4>
<div class="content">
正規母集団から単純無作為抽出を行う場合、特にサンプルサイズが小さいと、
確率変数が正規分布に従っているかどうかは標本データのヒストグラムを見てもわからない(左図)。
しかし大きなサイズになれば、ほぼ正規分布に従うことわかる(右図).

正規母集団から単純無作為抽出を行う

正規母集団から単純無作為抽出を行う場合、特にサンプルサイズが小さいと、 確率変数が正規分布に従っているかどうかは標本データのヒストグラムを見てもわからない(左図)。 しかし大きなサイズになれば、ほぼ正規分布に従うことわかる(右図).

In [ ]:
 
import numpy.random as random
help(random.normal)
In [106]:
 
samplesSmall = random.normal(50, 10.0, 10)  # N(50,10^2)に従う乱数を10個生成
samplesSmall
Out[106]:
array([ 33.09788212,  58.53310324,  46.89549989,  73.32767525,
        59.72873003,  57.4554278 ,  54.04003881,  35.12864043,
        40.43082595,  56.46944824])
In [107]:
 
plt.hist(samplesSmall)
plt.show()
In [108]:
 
samplesLarge = random.normal(50, 10.0, 10000)  # N(50,10^2)に従う乱数を10000個生成
plt.hist(samplesLarge)
plt.show()
 
<H2><A NAME="RS03:sampleDistribution">標本分布</A></H2>
<div class="content">
<font class="important">標本分布</font>とは、標本統計量(標本平均、標本分散など)に関する確率分布のこと
<P>
標本分布は標本における個々のデータの実現値を表した度数分布ではなく、標本統計量の確率分布<BR>
標本分布は母集団分布と標本統計量の種類、そしてサンプルサイズが決まると理論的(数学的)に導かれる。
実際のデータから作成されるわけではない
</div>

標本分布

標本分布とは、標本統計量(標本平均、標本分散など)に関する確率分布のこと

標本分布は標本における個々のデータの実現値を表した度数分布ではなく、標本統計量の確率分布
標本分布は母集団分布と標本統計量の種類、そしてサンプルサイズが決まると理論的(数学的)に導かれる。 実際のデータから作成されるわけではない

 
<H3>標本分布から何がわかるのか</H3>
<div class="content">
標本分布を調べることで、推定値の性質が分かる可能性がある
<BR>
標本分布を調べるときの観点: <OL>
  <LI>標本分布が母数の本当の値を中心として分布しているか⇒平均
  <LI>標本分布が横に大きく広がっていないか⇒標準偏差
</OL>
  (1)が満たされないと推定値と母数との誤差が大きくなる<BR>
  (2) は推定値にどの程度誤差が生じるか調べるもの
</div>

標本分布から何がわかるのか

標本分布を調べることで、推定値の性質が分かる可能性がある
標本分布を調べるときの観点:
  1. 標本分布が母数の本当の値を中心として分布しているか⇒平均
  2. 標本分布が横に大きく広がっていないか⇒標準偏差
  (1)が満たされないと推定値と母数との誤差が大きくなる
  (2) は推定値にどの程度誤差が生じるか調べるもの
 
<H3>標本分布を「経験的」に求める</H3>
<div class="content">
Pythonを用いて「経験的」に標本分布を求める方法---「経験的」とは「現実に得られたデータに基づく」という意味、
理論的な標本分布を「近似的」に再現
<BR>
母集団からサンプルサイズnの標本を何度も繰り返し抽出し、そのたびに標本統計量の実現値を計算して記録する
<BR>
これを実行するときの問題<UL><LI> 母集団分布がどのような分布になるかが不明
 <LI> 母集団分布がわからないと、標本から得られる推定値の信頼性も不明<BR>
⇒ 母集団が正規分布であると仮定して「もし母集団分布がこのような正規分布ならば、このくらいあてになる推定値が得られる」ということを検討
</UL>
</div>

標本分布を「経験的」に求める

Pythonを用いて「経験的」に標本分布を求める方法---「経験的」とは「現実に得られたデータに基づく」という意味、 理論的な標本分布を「近似的」に再現
母集団からサンプルサイズnの標本を何度も繰り返し抽出し、そのたびに標本統計量の実現値を計算して記録する
これを実行するときの問題
  • 母集団分布がどのような分布になるかが不明  
  • 母集団分布がわからないと、標本から得られる推定値の信頼性も不明
    ⇒ 母集団が正規分布であると仮定して「もし母集団分布がこのような正規分布ならば、このくらいあてになる推定値が得られる」ということを検討
 
<H3>正規母集団の母平均の推定</H3>
<div class="content">
仮定:母集団分布は$N(50,10^2)$、 サンプルサイズ n =10

正規母集団の母平均の推定

仮定:母集団分布はN(50,102)、 サンプルサイズ n =10

In [109]:
 
samples = random.normal(50, 10.0, 10)  # N(50,10^2)に従う乱数を10個生成
samples # データの確認
Out[109]:
array([ 48.19415007,  63.26636104,  57.49745353,  39.96272656,
        51.64073697,  60.36259109,  58.71608459,  63.91610309,
        65.66073452,  45.49907224])
In [110]:
 
np.mean(samples)   # 平均
Out[110]:
55.471601370394922
 
この結果、母平均の推定値は 55.5 
[疑問: 55.4716... なのに答は 55.5?]
注意: 55.4716... と表示されたとしても、答として、これをこのまま答えないこと。
<B>有効数字の桁数</B>を考慮して、55.5 と答える(有効数字は2桁、計算途中なので3桁とする)。
母平均は50なので、推定値との差は 5.5 --- この差が標本誤差(誤差)
一回きりの標本抽出で推定値がどれくらいの誤差をもつかは、いろいろ

この結果、母平均の推定値は 55.5

[疑問: 55.4716... なのに答は 55.5?] 注意: 55.4716... と表示されたとしても、答として、これをこのまま答えないこと。

有効数字の桁数を考慮して、55.5 と答える(有効数字は2桁、計算途中なので3桁とする)。

母平均は50なので、推定値との差は 5.5 --- この差が標本誤差(誤差)

一回きりの標本抽出で推定値がどれくらいの誤差をもつかは、いろいろ

 
<H3><A NAME="RS03:makingSample">標本分布を求める</A></H3>
<div class="content">
先の操作を10,000回繰り返してみよう:

標本分布を求める

先の操作を10,000回繰り返してみよう:

In [3]:
 
import numpy as np
import numpy.random as random
meanArray = np.array([ np.mean(random.normal(50, 10.0, 10)) for _ in range(10000)])
 
こうして得られた標本平均のヒストグラムと母集団の正規分布とを、次のコードにより重ねあわせてみると、かなり近い関係にあることが分かります:

こうして得られた標本平均のヒストグラムと母集団の正規分布とを、次のコードにより重ねあわせてみると、かなり近い関係にあることが分かります:

In [4]:
 
from scipy import stats as st
import matplotlib.pyplot as plt
%matplotlib inline
plt.hist(meanArray,normed=True)  # 頻度ではなく割合で表示
x = np.linspace(35, 65, 10000)
plt.plot(x, st.norm.pdf(x,50,10.0**0.5), ls='--',label=r'N(%i,%i)$' % (0,1))
plt.grid()
plt.show()
 
平均$\mu$、分散$\sigma^2$の正規分布 $N(\mu, \sigma^2)$に従う母集団からサンプルサイズ$n$の標本を抽出したとき、<font class="rule">標本平均の標本分布は、理論的に$N(\mu, \sigma^2/n)$ に従います</font>
<P>
誤差の絶対値が5以下、すなわち母平均の推定値が45以上、55以下の範囲で得られた回数を調べてみましょう:

平均μ、分散σ2の正規分布 N(μ,σ2)に従う母集団からサンプルサイズnの標本を抽出したとき、標本平均の標本分布は、理論的にN(μ,σ2/n) に従います

誤差の絶対値が5以下、すなわち母平均の推定値が45以上、55以下の範囲で得られた回数を調べてみましょう:

In [113]:
 
# 45以上55以下のものを1それ以外を0としたリストを作り、そこで0でない(つまり1)の個数を数える
np.count_nonzero([1 if 45 <= x <= 55 else 0 for x in meanArray]) 
Out[113]:
8868
 
全体の89%は本当の母数±5に収まっていることがわかります。理論的にはこのデータは$N(50,10^2)$に従うので、次章の知識を先取りしscipy.statsモジュールの関数cdfを用いた理論的な値と比較してみましょう:

全体の89%は本当の母数±5に収まっていることがわかります。理論的にはこのデータはN(50,102)に従うので、次章の知識を先取りしscipy.statsモジュールの関数cdfを用いた理論的な値と比較してみましょう:

In [63]:
 
norm.cdf(55, loc=50, scale=10.0**0.5) - norm.cdf(45,loc=50,scale=10.0**0.5)
Out[63]:
0.88615370199334187
In [ ]:
 
from scipy.stats import norm
help(norm.cdf)
 
このように、 上位二桁まで一致していることが確かめられました。
<P>
平均$\mu$、分散$\sigma^2$の正規分布 $N(\mu,\sigma^2)$ に従う母集団からサンプルサイズ$n$の標本を抽出したとき、<font class="rule">標本平均の標本分布は、理論的に$N(\mu,\sigma^2/n)$に従う</font>ことは、次のようにして実験的に確かめることができます。

このように、 上位二桁まで一致していることが確かめられました。

平均μ、分散σ2の正規分布 N(μ,σ2) に従う母集団からサンプルサイズnの標本を抽出したとき、標本平均の標本分布は、理論的にN(μ,σ2/n)に従うことは、次のようにして実験的に確かめることができます。

In [114]:
 
import numpy.random as random
for k in [50,60,70,80,90,100,150,200,250,300]:    # サンプルサイズを指定
     meanArray = np.array( [np.mean(random.normal(50, 10.0, k))
                           for _ in range(10000)])  # N(50,10^2)に従うk個の乱数の平均を10,000個
     print("サンプルサイズ=%d, 平均=%f, 標本分散=%f, 標本分散*サンプルサイズ=%f"
               % (k,meanArray.mean(), meanArray.var(),meanArray.var()*k))
サンプルサイズ=50, 平均=50.001157, 標本分散=2.018868, 標本分散*サンプルサイズ=100.943413
サンプルサイズ=60, 平均=49.989521, 標本分散=1.637274, 標本分散*サンプルサイズ=98.236410
サンプルサイズ=70, 平均=50.011332, 標本分散=1.432251, 標本分散*サンプルサイズ=100.257560
サンプルサイズ=80, 平均=50.002751, 標本分散=1.255702, 標本分散*サンプルサイズ=100.456120
サンプルサイズ=90, 平均=50.011754, 標本分散=1.103067, 標本分散*サンプルサイズ=99.276044
サンプルサイズ=100, 平均=49.991631, 標本分散=1.002772, 標本分散*サンプルサイズ=100.277168
サンプルサイズ=150, 平均=49.996742, 標本分散=0.664137, 標本分散*サンプルサイズ=99.620536
サンプルサイズ=200, 平均=50.005825, 標本分散=0.500871, 標本分散*サンプルサイズ=100.174264
サンプルサイズ=250, 平均=50.000568, 標本分散=0.405689, 標本分散*サンプルサイズ=101.422301
サンプルサイズ=300, 平均=49.997763, 標本分散=0.329371, 標本分散*サンプルサイズ=98.811296
 
以上の結果から、サンプルサイズを変えても標本平均の平均が母平均にほぼ一致すること、また標本平均の分散は、母分散/サンプルサイズであることが確かめられました。

以上の結果から、サンプルサイズを変えても標本平均の平均が母平均にほぼ一致すること、また標本平均の分散は、母分散/サンプルサイズであることが確かめられました。

 
<H3>不偏性</H3>
<div class="content">
ある推定量の標本分布の平均が、推定しようとしている母数の値と一致するとき、その推定量は<font class="word">不偏性</font>があるといいます。また、
不偏性がある推定量のことを<font class="word">不偏推定量</font>と言います。
推定量の不偏性は「標本分布が母数の本当の値を中心として分布しているか」に対応した概念です。
</div>

不偏性

ある推定量の標本分布の平均が、推定しようとしている母数の値と一致するとき、その推定量は不偏性があるといいます。また、 不偏性がある推定量のことを不偏推定量と言います。 推定量の不偏性は「標本分布が母数の本当の値を中心として分布しているか」に対応した概念です。
 
<H3><A NAME="RS03:normalError">標準誤差</A></H3>
<div class="content">
「推定量の標本分布の標準偏差(へんさ)」を<font class="word">標準誤差</font>(ごさ)と言います。
これは「標本分布の横への広がり」を表す値です。標準誤差が小さいということは、
標本分布の広がりが小さいことを表し、これは『ばらつきが少ない』ことを意味します。
<P>
例: $N(50,10^2)$の正規母集団から n=10 の標本を抽出したとき、
<font color="green">標本平均の標本分布は$N(50,10^2/n)$に従う</font>ので、標準誤差は$\sqrt{n}$、つまり(n=10なので) およそ3.2となります。
ということは、標本平均から母集団の平均を推定しようとすると、おおよそ68%の確率で46.8から53.2の間の値が得られますが、運が悪いと5以上の誤差が生じることがあることがわかります。
<P>
この例から、(1)母集団の分散(母分散)が大きければ大きいほど標本の平均値は母平均から離れた値を取りやすい(標準誤差が大きくなりやすい)、(2)標本のサイズが小さいほど標準誤差が大きい、ということがわかります。
言い換えれば、標準誤差を小さくするには<font class="rule">サンプルサイズnを大きく</font>するのが有効な方法であると言えます。
<P>
一般に、どのような標本統計量でも、サンプルサイズが大きいほど標準誤差は小さくなります
<P>
<font class="important">中心極限定理</font>: 母集団が正規分布でなくとも、<font class="important">標本平均の標本分布は「サンプルサイズnが大きい時」はほぼ正規分布になる</font>

標準誤差

「推定量の標本分布の標準偏差(へんさ)」を標準誤差(ごさ)と言います。 これは「標本分布の横への広がり」を表す値です。標準誤差が小さいということは、 標本分布の広がりが小さいことを表し、これは『ばらつきが少ない』ことを意味します。

例: N(50,102)の正規母集団から n=10 の標本を抽出したとき、

標本平均の標本分布はN(50,102/n)に従うので、標準誤差はn、つまり(n=10なので) およそ3.2となります。 ということは、標本平均から母集団の平均を推定しようとすると、おおよそ68%の確率で46.8から53.2の間の値が得られますが、運が悪いと5以上の誤差が生じることがあることがわかります。

この例から、(1)母集団の分散(母分散)が大きければ大きいほど標本の平均値は母平均から離れた値を取りやすい(標準誤差が大きくなりやすい)、(2)標本のサイズが小さいほど標準誤差が大きい、ということがわかります。 言い換えれば、標準誤差を小さくするにはサンプルサイズnを大きくするのが有効な方法であると言えます。

一般に、どのような標本統計量でも、サンプルサイズが大きいほど標準誤差は小さくなります

中心極限定理: 母集団が正規分布でなくとも、標本平均の標本分布は「サンプルサイズnが大きい時」はほぼ正規分布になる

 
<div class="exercise">
<h4>課題3-2</h4>
<A HREF="#makingSample">標本分布を求める</A>の項目ではn=10の標本を10,000個作り、その標本平均を求めた。
この課題ではn=100として10,000個の標本平均を求めてみよう。そしてn=10の場合のヒストグラムとn=100の場合のヒストグラムを比較し、標準誤差が小さくなっていることを確かめてみよ。

課題3-2

標本分布を求めるの項目ではn=10の標本を10,000個作り、その標本平均を求めた。 この課題ではn=100として10,000個の標本平均を求めてみよう。そしてn=10の場合のヒストグラムとn=100の場合のヒストグラムを比較し、標準誤差が小さくなっていることを確かめてみよ。

In [115]:
 
%matplotlib inline
import matplotlib.pyplot as plt
import numpy.random as random
SampleAverage10 = np.array([np.mean(random.normal(50, 10.0, 10)) for _ in range(10000)])
SampleAverage100 = np.array([np.mean(random.normal(50, 10.0, 100)) for _ in range(10000)])
plt.hist(SampleAverage10, color="yellow")     # n=10の方を黄色で表示
plt.hist(SampleAverage100, color="green")        # n=100の方を緑色で表示
plt.title('Histogram of sample average')
plt.xlabel('x')
plt.ylabel('Frequency')
plt.show()
print("std of sample10=%f, sample100=%f"%(SampleAverage10.std(), SampleAverage100.std()))
std of sample10=3.137151, sample100=0.990142
 
<H3>標本平均以外の標本分布</H3>
<P>
「平均」以外の母数のうち、分散と中央値の標本分布について述べます。
<P>
<h4>標本分散と不偏分散の標本分布</H4>
<P>
母分散の推定量としての<font class="word">標本分散と不偏分散の違い</font>を実験で示しましょう。

標本平均以外の標本分布

「平均」以外の母数のうち、分散と中央値の標本分布について述べます。

標本分散と不偏分散の標本分布

母分散の推定量としての標本分散と不偏分散の違いを実験で示しましょう。

In [117]:
 
import numpy as np
import numpy.random as random
SampleVars = []
UnbiasedVars = []
for i in range(10000):
    Sample = random.normal(50, 10.0,10)  # N(50, 10^2)に従う乱数10個
    SampleVars = SampleVars+[Sample.var()]  # 標本分散
    UnbiasedVars = UnbiasedVars + [Sample.var(ddof=1.0)]
SampleVars = np.array(SampleVars)
UnbiasedVars = np.array(UnbiasedVars)
print(SampleVars.mean())
print(UnbiasedVars.mean())
90.2677729196
100.297525466
 
平均的に母分散の推定値として不偏分散のほうが近い値であることがわかります。次にその標準偏差を求めてみましょう。

平均的に母分散の推定値として不偏分散のほうが近い値であることがわかります。次にその標準偏差を求めてみましょう。

In [118]:
 
print(SampleVars.std())
print(UnbiasedVars.std())
42.4670221473
47.1855801637
 
不偏分散のほうが推定値のばらつきが大きいことがわかります。

不偏分散のほうが推定値のばらつきが大きいことがわかります。

In [119]:
 
plt.subplot(2,2,1)
plt.hist(SampleVars)
plt.title(u'Sample Variance')
plt.subplot(2,2,2)
plt.hist(UnbiasedVars)
plt.title(u'Unbiased Variance')
plt.show()
In [120]:
 
plt.subplot(2,2,1)
plt.hist(SampleVars,bins=50)
plt.title(u'Sample Variance')
plt.subplot(2,2,2)
plt.hist(UnbiasedVars,bins=50)
plt.title(u'Unbiased Variance')
plt.show()
 
<UL>
  <LI>不偏分散は不偏性をもつという点では望ましい性質をもっているが、
<FONT class="rule">標本分散に比べて過大な推定値が得られる可能性が高い</font>
   <LI>不偏分散の平方根は母標準偏差の不偏推定量ではない
</UL>
  • 不偏分散は不偏性をもつという点では望ましい性質をもっているが、 標本分散に比べて過大な推定値が得られる可能性が高い
  • 不偏分散の平方根は母標準偏差の不偏推定量ではない
 
<h4>中央値の標本分布</H4>
<P>
実験: 標本中央値の標本分布を求め、標本平均と比較する(どちらが母平均の推定量としてよいか?)

中央値の標本分布

実験: 標本中央値の標本分布を求め、標本平均と比較する(どちらが母平均の推定量としてよいか?)

In [121]:
 
import numpy.random as random
SampleAverage = []
SampleMedian = []
for i in range(10000):
    Sample = random.normal(50, 10.0, 10)  # N(50, 10^2)に従う乱数10個
    SampleAverage = SampleAverage+[np.mean(Sample)]  # 標本平均
    SampleMedian = SampleMedian + [np.median(Sample)]   # 中央値
SampleAverage = np.array(SampleAverage)
SampleMedian = np.array(SampleMedian)
print(SampleAverage.mean())
print(SampleMedian.mean())
50.0321506814
50.0576939024
 
標本中央値も標本平均も母平均とほぼ一致---母平均の推定量としての候補として優劣が付けられない

標本中央値も標本平均も母平均とほぼ一致---母平均の推定量としての候補として優劣が付けられない

In [122]:
 
print(np.std(SampleAverage))
print(np.std(SampleMedian))
3.15787180632
3.6923757235
 
標本平均のほうが標準誤差が小さい---標本平均は母平均の不偏推定量の中で</font><font class="word">標準誤差が最小</font>

標本平均のほうが標準誤差が小さい---標本平均は母平均の不偏推定量の中で標準誤差が最小

In [123]:
 
plt.subplot(2,2,1)
plt.hist(SampleAverage,bins=50)
plt.title(u'Sample Mean')
plt.subplot(2,2,2)
plt.hist(SampleMedian,bins=50)
plt.title(u'Sample Median')
plt.show()
 
---
<H2><A NAME="Functions">関数のまとめ</A></H2>
注:  numpyをnp, np.randomを random、 matplotlib.pyplotをplt、pandasをpd、scipy.statsをstと略記する
|目的 | 関数名とモジュール | 使い方|
|:---|:---|:---|
|一様分布に従う乱数(復習)| random.random(個数) | random.random(10) # [0,1)の乱数を10個|
|正規分布に従う乱数 | random.normal(平均,標準偏差,個数) | random.normal(50,10,10)  # N(50, 10^2)に従う乱数10個|
|指定の範囲の数の生成 | np.linspace(最小値,最大値,個数) | np.linspace(35, 65, 10000) # [35,65]の範囲の数を10000個|
|正規分布の密度関数| st.norm.pdf(xデータ,平均,標準偏差) | plt.plot(x, st.norm.pdf(x,50,10.0**0.5)) # $N(50,\sqrt{10})$の正規分布の関数描画|

関数のまとめ

注: numpyをnp, np.randomを random、 matplotlib.pyplotをplt、pandasをpd、scipy.statsをstと略記する

目的 関数名とモジュール 使い方
一様分布に従う乱数(復習) random.random(個数) random.random(10) # [0,1)の乱数を10個
正規分布に従う乱数 random.normal(平均,標準偏差,個数) random.normal(50,10,10) # N(50, 10^2)に従う乱数10個
指定の範囲の数の生成 np.linspace(最小値,最大値,個数) np.linspace(35, 65, 10000) # [35,65]の範囲の数を10000個
正規分布の密度関数 st.norm.pdf(xデータ,平均,標準偏差) plt.plot(x, st.norm.pdf(x,50,10.0**0.5)) # N(50,10)の正規分布の関数描画
 
<hr noshade="noshade" size="5">
<H2><A NAME="RS03:Exercises">演習問題3</A></H2>
<div class="exercise">
<h4>演習問題3-1</h4>
<P>
 $N(50,10^2)$の正規母集団から$n = 20$の標本抽出を5,000回繰り返すことにより、標本平均の経験的な標本分布を求めよ。またこの経験的な標本分布から得られるグラフ(ヒストグラム)に、理論的に導かれる標本分布を重ねあわせて表示せよ。<BR>
この結果下のようなグラフが得られればよい:
<CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-A.png" target="_blank">
<IMG class="fig"  SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-A.png" border="0"  width="50%"  hspace="10"> 
</A>
</CENTER>


演習問題3

演習問題3-1

N(50,102)の正規母集団からn=20の標本抽出を5,000回繰り返すことにより、標本平均の経験的な標本分布を求めよ。またこの経験的な標本分布から得られるグラフ(ヒストグラム)に、理論的に導かれる標本分布を重ねあわせて表示せよ。

この結果下のようなグラフが得られればよい:

In [ ]:
 
 
<h4>演習問題3-2</h4>
<P>
理論的な標本分布について、サンプルサイズを$n=1, 4, 9, 16, 25$と変化させた時に、標本分布の形状がどのように代わるかを調べてみよ。ただし、母集団分布は標準正規分布$N(0, 1^2)$とする。
<BR>
これによって次のような図が得られれば良い。
<CENTER>
<A HREF=http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-B.png" target="_blank">
<IMG class="fig"  SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-B.png" border="0"  width="50%"  hspace="10"> 
</A>
</CENTER>

演習問題3-2

理論的な標本分布について、サンプルサイズをn=1,4,9,16,25と変化させた時に、標本分布の形状がどのように代わるかを調べてみよ。ただし、母集団分布は標準正規分布N(0,12)とする。


これによって次のような図が得られれば良い。

<A HREF=http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-B.png" target="_blank">

In [ ]:
 
 
<h4>演習問題3-3</h4>
<P>
最小値0、最大値100の一様分布から$n=20$の標本抽出を5,000回繰り返すことにより、標本平均の経験的な標本分布を求めよ。中心極限定理からこれは正規分布に近づくはずである。そこでこのヒストグラムに、標本平均データの平均と標準偏差による正規分布のグラフを重ねせて表示せよ。<BR>
この結果、次のようなグラフが得られればよい:
<CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-C.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-3-C.png" border="0"  width="50%"  hspace="10"> 
</A>
</CENTER>

演習問題3-3

最小値0、最大値100の一様分布からn=20の標本抽出を5,000回繰り返すことにより、標本平均の経験的な標本分布を求めよ。中心極限定理からこれは正規分布に近づくはずである。そこでこのヒストグラムに、標本平均データの平均と標準偏差による正規分布のグラフを重ねせて表示せよ。
この結果、次のようなグラフが得られればよい:

In [ ]: