xxxxxxxxxx
 
<h1>Pythonで統計学を学ぶ(4)</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>
この講義では、「統計的仮説検定」をとりあげます。
これは、統計的仮説検定の手順の理解と用語の習熟がねらいです。
また、代表的な統計的仮説検定、つまり標準正規分布を用いた検定、t分布を用いた検定、無相関検定、カイ二乗検定について学びます。
<p>学習項目です:
<ul class="important">
<li class="important"><A HREF="#RS04:necessity">統計的仮説検定の必要性</A></li>
<li class="important"><A HREF="#RS04:testProcedure">統計的仮説検定の手順と用語</A></li>
<li class="important"><A HREF="#RS04:averageTest">標準正規分布を用いた検定(1つの平均値の検定: 母分散が既知)</A></li>
<li class="important"><A HREF="#RS04:tTest">t分布を用いた検定(1つの平均値の検定: 母分散が未知)</A></li>
<li class="important"><A HREF="#RS04:correlationTest">相関係数の検定(無相関検定)</A></li>
<li class="important"><A HREF="#RS04:chisquareTest">独立性の検定(カイ二乗検定)</A></li>
<li class="important"><A HREF="#Functions">関数のまとめ</A></li>
<li class="important"><A HREF="#Exercises">演習問題</A></li>
</ul>

Pythonで統計学を学ぶ(4)

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

この講義では、「統計的仮説検定」をとりあげます。 これは、統計的仮説検定の手順の理解と用語の習熟がねらいです。 また、代表的な統計的仮説検定、つまり標準正規分布を用いた検定、t分布を用いた検定、無相関検定、カイ二乗検定について学びます。

学習項目です:

 
<h2><A NAME="RS04:necessity">統計的仮説検定の必要性</A></h2>
<div class="content">
下の散布図を見てください(青色の円は、点の分布の状態を表すために描いたものです):
<CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-1.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-1.png" width="60%" hspace="10" />
</A>
</CENTER>
これをみると、この2つの変数 a と bの間には相関関係がないようにみえます。実際
$corrcoef(a,b)= -0.034$なので相関なしといえます。
ところが、a, bそれぞれから30点ずつ無作為抽出したデータ xa, xb (下にその散布図を示す)は、
ときに 0.378 という『弱い相関』を示すことがあります。
<BR>
参考: 無相関の母集団から相関するデータを作る<BR>
次が無相関の母集団から相関するデータを作った方法です:

統計的仮説検定の必要性

下の散布図を見てください(青色の円は、点の分布の状態を表すために描いたものです):

これをみると、この2つの変数 a と bの間には相関関係がないようにみえます。実際 corrcoef(a,b)=0.034なので相関なしといえます。 ところが、a, bそれぞれから30点ずつ無作為抽出したデータ xa, xb (下にその散布図を示す)は、 ときに 0.378 という『弱い相関』を示すことがあります。


参考: 無相関の母集団から相関するデータを作る
次が無相関の母集団から相関するデータを作った方法です:

In [185]:
 
import matplotlib.pyplot as plt
fig = plt.figure()
a = random.normal(50,10,500)
b = random.normal(50,10,500)
plt.scatter(a, b,c='w')
plt.xlabel('a')
plt.ylabel('b')
ax = fig.add_subplot(111)
# 中心(50,50)で半径25の円を描画
circle = plt.Circle((50,50),25, fill=False, color='b')
ax.add_patch(circle)
aspect = (ax.get_xlim()[1] - ax.get_xlim()[0]) / (ax.get_ylim()[1] - ax.get_ylim()[0])                     
ax.set_aspect(aspect)  # 縦横の縮尺を調整
plt.show()
In [34]:
 
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import numpy.random as random
a = random.normal(50,10,500)
b = random.normal(50,10,500)
for _ in range(100):
    xa = random.choice(a,30)
    xb = random.choice(b,30)
    if (abs(np.corrcoef(xa,xb)[0,1]) > 0.30):
        break
plt.scatter(xa, xb)
plt.xlabel('xa')
plt.ylabel('xb')
x = np.linspace(10, 90, 10000)
lm = np.polyfit(xa, xb, 1)
plt.plot(x, lm[0]*x+lm[1],"g")
print(np.corrcoef(xa,xb)[0,1])
0.306203680441
 
これは、実際には作為的に作られたデータです。しかし、あなたが発明した機器の有効性を示す論文や、あなたが作った薬品の効果を示す論文において、都合の良いデータだけを集めたのではないか、もしくは作為がないとしても、このデータは本当に偶然の結果であって多数のデータを取ればこのようなグラフにはならないのではないか、という疑いがかけられることがあります。
<P>
そのような疑いや批判には、(前の章で学んだように)『標本抽出が無作為抽出であること』(都合の良いデータを集めたわけではない)、そして(本章で学ぶように)『母集団に全く相関がないとしたら、抽出した標本からこのような結果が得られる可能性が非常に小さいということ』(多数のデータを集めても同じような結果が得られる確率が高い)を示さなければなりません。そして、
<font class="word"><B>統計的仮説検定</B></font>は確率に基づき、後者の主張を行うための方法です。(前者の「無作為抽出」は、統計による分析の大前提です)

これは、実際には作為的に作られたデータです。しかし、あなたが発明した機器の有効性を示す論文や、あなたが作った薬品の効果を示す論文において、都合の良いデータだけを集めたのではないか、もしくは作為がないとしても、このデータは本当に偶然の結果であって多数のデータを取ればこのようなグラフにはならないのではないか、という疑いがかけられることがあります。

そのような疑いや批判には、(前の章で学んだように)『標本抽出が無作為抽出であること』(都合の良いデータを集めたわけではない)、そして(本章で学ぶように)『母集団に全く相関がないとしたら、抽出した標本からこのような結果が得られる可能性が非常に小さいということ』(多数のデータを集めても同じような結果が得られる確率が高い)を示さなければなりません。そして、

統計的仮説検定は確率に基づき、後者の主張を行うための方法です。(前者の「無作為抽出」は、統計による分析の大前提です)

<h2><A NAME="RS04:testProcedure">統計的仮説検定の手順と用語</A></h2>
<div class="content">
統計的仮説検定の一般的な手順を次の表に示します:
|手順|やること|
|:---|:---|
|1 | 母集団に関する<font class="word">帰無仮説</font><font class="word">対立仮説</font>を設定|
|2 |<font class="word">検定統計量</font>を選択|
|3 | <font class="word">有意水準 α</font>の値を決定|
|4 | データを収集した後、データから<font class="word">検定統計量の実現値</font>を求める|
|5| 結論: <font class="word">検定統計量の実現値</font><font class="word">棄却域</font>に入れば<font class="word">帰無仮説</font>を棄却し、<font class="word">対立仮説</font>を採択する。そうでなければ、<font class="word">帰無仮説</font>を採択する|

統計的仮説検定の手順と用語

統計的仮説検定の一般的な手順を次の表に示します:
手順 やること
1 母集団に関する帰無仮説対立仮説を設定
2 検定統計量を選択
3 有意水準 αの値を決定
4 データを収集した後、データから検定統計量の実現値を求める
5 結論: 検定統計量の実現値棄却域に入れば帰無仮説を棄却し、対立仮説を採択する。そうでなければ、帰無仮説を採択する

<H3>1. 帰無仮説と対立仮説</H3>
<ul class="important">
   <LI><font class="important">帰無仮説</font>:提案する手法が従来の手法と「差がない」、
提案する手法は「効果がない」という仮説---<font color="green">本来主張したいこととは逆の仮説</font>
<font class="important">この仮説が棄却されることを目標として仮説検定を行う</font></LI>
具体的には、母平均$\mu = 0$ (母平均は0である), 母相関係数$\rho =0$ (相関がない), 母平均の差$\mu_1 - \mu_2 = 0$ (差がない)というような仮説
    <LI> <font class="important">対立仮説</font>: 帰無仮説が棄却されたときに採択される仮説---
帰無仮説とは逆の仮説であり、実験などで<font class="important">示したい・主張したいことを表したもの</font></LI>
具体的には、母平均$\mu \neq 0$(母平均は<B>0でない</B>), 母相関係数 $\rho \neq 0$ (<B>相関がある</B>), 母平均の差$\mu_1 - \mu_2 \neq 0$ (<B>差がある</B>)というような仮説
</ul>
対立仮説の設定により、検定は次のどちらかで行う(両側検定の方がより厳しい条件であり、<font class="important">普通は両側検定</font>で行う):
<ul class="important">
   <LI><font class="important">片側検定</font>:対立仮説が、母平均$\mu > 0$(もしくは$\mu < 0$ )、
母相関係数 $\rho > 0$(もしくは$\rho < 0$)、、母平均の差$\mu_1 > \mu_2$ (もしくは$\mu_1 < \mu_2$ )、の場合</LI>
   <LI><font class="important">両側検定</font>:対立仮説が、母平均μ≠0、母相関係数ρ≠0 、母平均の差 μ<sub>1</sub> - μ<sub>2</sub>≠ 0の場合</LI>
 要するに、両側検定では、例えば母平均$\mu \neq 0$を調べるには、母平均$\mu > 0$$\mu < 0$ の両方を調べなければならない
</ul>
<font class="important">帰無仮説が正しいものとして分析を行う</font>
実際に得られたデータから計算された検定統計量の値によって採択を判断する。
<P>
<font class="rule">帰無仮説が正しいとしたとき</font>、検定統計量が、ほぼ起こり得ない値(それほど極端な値)であれば、帰無仮説を棄却する(つまり、<font color="green">本来の主張を表す対立仮説が採択される</font>)。そうでなければ(確率的に十分起こりうるような値であれば、帰無仮説を採択する(この場合は、<font color="green">本来主張したかった対立仮説が棄却されてしまう</font>)。

1. 帰無仮説と対立仮説

  • 帰無仮説:提案する手法が従来の手法と「差がない」、 提案する手法は「効果がない」という仮説---本来主張したいこととは逆の仮説

    この仮説が棄却されることを目標として仮説検定を行う

  • 具体的には、母平均μ=0 (母平均は0である), 母相関係数ρ=0 (相関がない), 母平均の差μ1μ2=0 (差がない)というような仮説

  • 対立仮説: 帰無仮説が棄却されたときに採択される仮説--- 帰無仮説とは逆の仮説であり、実験などで示したい・主張したいことを表したもの
  • 具体的には、母平均μ0(母平均は0でない), 母相関係数 ρ0 (相関がある), 母平均の差μ1μ20 (差がある)というような仮説

対立仮説の設定により、検定は次のどちらかで行う(両側検定の方がより厳しい条件であり、普通は両側検定で行う):

  • 片側検定:対立仮説が、母平均μ>0(もしくはμ<0 )、 母相関係数 ρ>0(もしくはρ<0)、、母平均の差μ1>μ2 (もしくはμ1<μ2 )、の場合
  • 両側検定:対立仮説が、母平均μ≠0、母相関係数ρ≠0 、母平均の差 μ1 - μ2≠ 0の場合
  • 要するに、両側検定では、例えば母平均μ0を調べるには、母平均μ>0μ<0 の両方を調べなければならない

帰無仮説が正しいものとして分析を行う。 実際に得られたデータから計算された検定統計量の値によって採択を判断する。

帰無仮説が正しいとしたとき、検定統計量が、ほぼ起こり得ない値(それほど極端な値)であれば、帰無仮説を棄却する(つまり、本来の主張を表す対立仮説が採択される)。そうでなければ(確率的に十分起こりうるような値であれば、帰無仮説を採択する(この場合は、本来主張したかった対立仮説が棄却されてしまう)。

xxxxxxxxxx
 
<H3>2. 検定統計量</H3>
<ul class="important">
    <LI><font class="important">検定統計量</font>:統計的仮説検定のために用いられる標本統計量のこと。代表的な検定統計量の例: $t, \chi^2、F$</LI>
    <LI><font class="important">検定統計量の実現値</font>:実際のデータ(手に入った標本)を基に計算してえられる具体的な値のこと</LI>
</ul>
検定統計量の実現値は、<font class="rule">対立仮説に合うほど 0から離れた値を示す</font>

2. 検定統計量

  • 検定統計量:統計的仮説検定のために用いられる標本統計量のこと。代表的な検定統計量の例: t,χ2F
  • 検定統計量の実現値:実際のデータ(手に入った標本)を基に計算してえられる具体的な値のこと
検定統計量の実現値は、対立仮説に合うほど 0から離れた値を示す

xxxxxxxxxx
 
<H3>3. 有意水準と棄却域</H3>
対立仮説を採択するか決定するときに基準になるのが<font class="word">有意水準</font>(αで表されます)
<BR>
有意水準は<font class="rule">5%または1%</font>(α=0.05またはα=0.01)に設定することが多い(つまり、標本が100回に5回(5%の場合)以下にしか現れないデータであった---こんなことは偶然では起こりえない---、だから帰無仮説が成り立たないと考えて良いのではないか、という判断基準)
<BR>
帰無仮説が正しいものとして考えた時の標本分布を<font class="word">帰無分布</font>という---帰無分布に基づいて確率計算される
<BR>
帰無仮説のもとで、非常に生じにくい検定統計量の値の範囲を<font class="word">棄却域</font>という---帰無仮説が棄却される領域(だから、この範囲に入るのが『望ましい』)
<UL>
    <LI><font class="word">採択域</font>: 棄却域以外の部分---「帰無仮説が採択される領域」</LI>
    <LI><font class="word">臨界値</font>: 棄却域と採択域の境目の値</LI>
</UL>
棄却域に検定統計量の実現値が入ったら、帰無仮説を棄却する---本来主張したかったことが採択される!
<CENTER>
<TABLE>
<TR><TD align="center">
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-3.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-3.png" width="70%" hspace="10" />
</A>
</TD></TR>
 <TR><TD align="center">正規分布を帰無分布とした時の棄却域
 </TD></TR></TABLE>
</CENTER>

3. 有意水準と棄却域

対立仮説を採択するか決定するときに基準になるのが有意水準(αで表されます)


有意水準は5%または1%(α=0.05またはα=0.01)に設定することが多い(つまり、標本が100回に5回(5%の場合)以下にしか現れないデータであった---こんなことは偶然では起こりえない---、だから帰無仮説が成り立たないと考えて良いのではないか、という判断基準)


帰無仮説が正しいものとして考えた時の標本分布を帰無分布という---帰無分布に基づいて確率計算される


帰無仮説のもとで、非常に生じにくい検定統計量の値の範囲を棄却域という---帰無仮説が棄却される領域(だから、この範囲に入るのが『望ましい』)

  • 採択域: 棄却域以外の部分---「帰無仮説が採択される領域」
  • 臨界値: 棄却域と採択域の境目の値
棄却域に検定統計量の実現値が入ったら、帰無仮説を棄却する---本来主張したかったことが採択される!

正規分布を帰無分布とした時の棄却域
 
<H3>4 &amp; 5. 統計的仮説検定の結果の報告</H3>
<div class="content">
検定統計量の実現値が棄却域に入った場合、「差がない」という帰無仮説を棄却し、
「差がある」という対立仮説を採択する。
<ul class="important">
検定結果は5% (または1%)水準で有意である<BR>
または
<I>p</I> &lt; .05 (または <I>p</I> &lt; .01 )で<font class="important">有意差</font>が見られた 」
と記述する。
</UL>
帰無仮説が棄却できない場合は、「検定の結果、差が有意でなかった」または「有意差が認められなかった」と書く。
</div>

4 & 5. 統計的仮説検定の結果の報告

検定統計量の実現値が棄却域に入った場合、「差がない」という帰無仮説を棄却し、 「差がある」という対立仮説を採択する。
    検定結果は5% (または1%)水準で有意である
    または 「p < .05 (または p < .01 )で有意差が見られた 」 と記述する。
帰無仮説が棄却できない場合は、「検定の結果、差が有意でなかった」または「有意差が認められなかった」と書く。
 
<h4>課題4-1</h4>
あなたはランダムに配置された対象物(例えば地雷や石油や埋蔵金など)を衛星からのセンサーデータを元に限定された時間(例えば1時間)内に検出する機器を作成した。100個のデータに対し検出率は0.70であった。そして、その性能が従来の製品(検出率は0.60と宣伝されている)よりも優れていることを統計的仮説検定の手法により示したい。<BR>
どのような帰無仮説と対立仮説をたてればよいか、また検定方法は片側か両側か、有意水準はどのくらいに設定したらよいか、考えを述べよ。

課題4-1

あなたはランダムに配置された対象物(例えば地雷や石油や埋蔵金など)を衛星からのセンサーデータを元に限定された時間(例えば1時間)内に検出する機器を作成した。100個のデータに対し検出率は0.70であった。そして、その性能が従来の製品(検出率は0.60と宣伝されている)よりも優れていることを統計的仮説検定の手法により示したい。
どのような帰無仮説と対立仮説をたてればよいか、また検定方法は片側か両側か、有意水準はどのくらいに設定したらよいか、考えを述べよ。

In [ ]:
 
<H3>p値</H3>
<div class="content">
<font class="important">p値</font><B>帰無仮説が正しいという仮定のもとで</B>、標本から計算した検定統計量の実現値以上の値が得られる確率
p値が有意水準より小さい時に帰無仮説を棄却する
[参考: p値が小さいことの意味]
p値の大きさが対立仮説を採択する(帰無仮説を棄却する)決め手となります。p値が小さいということは、
『帰無仮説が正しいとすると』確率的にほとんど起こりえないことが起きた(有意水準が5%なら100回中5回以下、
1%なら100回中1回以下)ということを意味します。逆にp値が大きいということは、確率的にはよくあることが起きた(だから、この結果では差があるとはいえない)、
ということになります。

p値

p値帰無仮説が正しいという仮定のもとで、標本から計算した検定統計量の実現値以上の値が得られる確率

p値が有意水準より小さい時に帰無仮説を棄却する

[参考: p値が小さいことの意味] p値の大きさが対立仮説を採択する(帰無仮説を棄却する)決め手となります。p値が小さいということは、 『帰無仮説が正しいとすると』確率的にほとんど起こりえないことが起きた(有意水準が5%なら100回中5回以下、 1%なら100回中1回以下)ということを意味します。逆にp値が大きいということは、確率的にはよくあることが起きた(だから、この結果では差があるとはいえない)、 ということになります。

xxxxxxxxxx
<H3>第1種の誤りと第2種の誤り</H3>
<UL>
 <LI><font class="word">第1種の誤り α</font>: 「<font class="rule">帰無仮説が真のとき</font>、これを<font class="rule">棄却</font>してしまう」誤りのこと</LI>
   この種の誤りを犯す確率が「有意水準」または「危険率」
 <LI><font class="word">第2種の誤り β</font>:「<font class="rule">帰無仮説が偽のとき</font>、これを<font class="rule">採択</font>する(棄却できない)」誤りのこと</LI>
本当は差があるのに「差がない」と判断してしまう誤り
</UL>

第1種の誤りと第2種の誤り

  • 第1種の誤り α: 「帰無仮説が真のとき、これを棄却してしまう」誤りのこと
  • この種の誤りを犯す確率が「有意水準」または「危険率」
  • 第2種の誤り β:「帰無仮説が偽のとき、これを採択する(棄却できない)」誤りのこと
  • 本当は差があるのに「差がない」と判断してしまう誤り
xxxxxxxxxx
 
<H3>検定力</H3>
<font class="important">検定力</font><font class="rule">帰無仮説が偽の場合</font>、全体の確率1から第2種の誤りの確率(1 - β)を引いた確率 <BR>「第2種の誤りを犯さない確率」とも、つまり間違っている帰無仮説を正しく棄却できる確率のこと 

検定力

検定力帰無仮説が偽の場合、全体の確率1から第2種の誤りの確率(1 - β)を引いた確率
「第2種の誤りを犯さない確率」とも、つまり間違っている帰無仮説を正しく棄却できる確率のこと 

<hr noshade size=3>
<h2><A NAME="RS04:averageTest">標準正規分布を用いた検定(1つの平均値の検定: 母分散が既知)</A></h2>
<div class="content">
正規母集団 $N(\mu, \sigma^2)$ から無作為に標本を抽出する(サンプルサイズを$n$ とする)と<BR>
  標本平均の分布も正規分布<BR>
  標本平均の平均は [ア] 、分散は [イ ] (<font color="red">問題</font>:ア、イに当てはまる記号を書け---課題4-2)<BR>
これを標準化したものを<font class="word">検定統計量</font>とする($\bar{X}$は標本データの平均):
$$Z = \frac{\bar{X}-\mu}{\sigma/\sqrt{n}}$$


標準正規分布を用いた検定(1つの平均値の検定: 母分散が既知)

正規母集団 N(μ,σ2) から無作為に標本を抽出する(サンプルサイズをn とする)と
  標本平均の分布も正規分布
  標本平均の平均は [ア] 、分散は [イ ] (問題:ア、イに当てはまる記号を書け---課題4-2)

これを標準化したものを検定統計量とする(X¯は標本データの平均):

Z=X¯μσ/n

<h4>課題4-2</h4>
  正規母集団 $N(\mu, \sigma^2)$から無作為に標本を抽出したとき、
理論的に標本平均の平均と、分散がそれぞれどのように表されるか、書きなさい(つまり、上の[ア], [イ]の箇所を補うこと)。
またこれを標準化して得られる検定統計量がZで表されている理由を答えなさい。
[ヒント]
<A HREF="Rstatistics-03.html#makingSample">標本分布を求める</A>の項を読みなおしてください。また、標準化については<A HREF="Rstatistics-01.html#RS01:normalization">標準化</A>の項を見てください。

課題4-2

正規母集団 N(μ,σ2)から無作為に標本を抽出したとき、 理論的に標本平均の平均と、分散がそれぞれどのように表されるか、書きなさい(つまり、上の[ア], [イ]の箇所を補うこと)。

またこれを標準化して得られる検定統計量がZで表されている理由を答えなさい。

[ヒント]

標本分布を求めるの項を読みなおしてください。また、標準化については標準化の項を見てください。

xxxxxxxxxx
 
<!--- 解答用 --->
 
---
<H4>Pythonを使った実習</H4>
<B>例題:</B>「心理学テスト」が$N(12, 10)$の正規分布に従うものとする。
次のデータ(<A NAME="SampleData">「指導法データ」</A>と呼ぶ)はこの母集団から無作為抽出した標本と考えてよいかどうかを判定せよ

Pythonを使った実習

例題:「心理学テスト」がN(12,10)の正規分布に従うものとする。 次のデータ(「指導法データ」と呼ぶ)はこの母集団から無作為抽出した標本と考えてよいかどうかを判定せよ

In [36]:
 
from __future__ import division
import numpy as np
SampleData = np.array([13,14,7,12,10,6,8,15,4,14,9,6,10,12,5,12,8,8,12,15])
次のステップで行う:<OL class="other">
 <LI>帰無仮説と対立仮説をたてる: 帰無仮説は「無作為抽出した標本と考えて良い」、つまり$\mu = 12$</LI>
     対立仮説は「無作為抽出した標本ではない」、つまり$\mu \neq 12$
 <LI>検定統計量の選択: 標本データを標準化した値(Zで表す)</LI>
 <LI>有意水準の決定: 両側検定で、有意水準 5%、つまり$\alpha=0.05$</LI>
 <LI>検定統計量の実現値の計算:</LI></OL>

次のステップで行う:

  1. 帰無仮説と対立仮説をたてる: 帰無仮説は「無作為抽出した標本と考えて良い」、つまりμ=12
  2. 対立仮説は「無作為抽出した標本ではない」、つまりμ12
  3. 検定統計量の選択: 標本データを標準化した値(Zで表す)
  4. 有意水準の決定: 両側検定で、有意水準 5%、つまりα=0.05
  5. 検定統計量の実現値の計算:

In [38]:
 
z =  (np.mean(SampleData) - 12) / (10.0/len(SampleData))**0.5             # 標準化
z
Out[38]:
-2.8284271247461898
xxxxxxxxxx
 
! 5. 帰無仮説の棄却か採択かの決定: <B>帰無仮説によればこの標本は正規分布に従う</B>。そこで<font class="word">scipy.stats</font>モジュールのppf関数で棄却の臨界値を求める、
もしくは<font class="word">cdf</font>関数でp値を求める
   <UL class="other">
    <LI><font class="word">下側確率</font>:標準正規分布に従う確率変数$Z$を例にとると、$Z$がある値$\alpha$以下となる確率 $Prob(Z \leq \alpha)$</LI>
    <LI><font class="word"> 上側確率</font>:標準正規分布に従う確率変数$Z$を例にとると、$Z$がある値$\alpha$より大きくなる確率  $Prob(Z > \alpha)$</LI>
   </UL>

! 5. 帰無仮説の棄却か採択かの決定: 帰無仮説によればこの標本は正規分布に従う。そこでscipy.statsモジュールのppf関数で棄却の臨界値を求める、 もしくはcdf関数でp値を求める

  • 下側確率:標準正規分布に従う確率変数Zを例にとると、Zがある値α以下となる確率 Prob(Zα)
  • 上側確率:標準正規分布に従う確率変数Zを例にとると、Zがある値αより大きくなる確率 Prob(Z>α)

In [45]:
 
import scipy.stats as st
st.norm.ppf(0.025)    # 下側確率0.05/2 = 0.025となるzの値を求める
# 下側確率であるから、この値よりもZ値が小さければ棄却される
Out[45]:
-1.9599639845400545
In [46]:
 
# 上側確率 1 - 0.05/2 = 0.975となるzの値を求める
st.norm.ppf(0.975)  
# 上側確率であるから、この値よりもZ値が大きければ棄却される
Out[46]:
1.959963984540054
In [4]:
 
help(st.norm.ppf)
Help on method ppf in module scipy.stats._distn_infrastructure:

ppf(self, q, *args, **kwds) method of scipy.stats._continuous_distns.norm_gen instance
    Percent point function (inverse of `cdf`) at q of the given RV.
    
    Parameters
    ----------
    q : array_like
        lower tail probability
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    x : array_like
        quantile corresponding to the lower tail probability q.

 
この結果、棄却域は $Z < -1.959964$ または $Z >  1.959964$となるので、$Z$の値は棄却域に入る。
よって、結論「有意水準5%において、指導法データは心理学テスト(という母集団)から無作為抽出した標本とはいえない」。
<BR>
なお、関数cdfを用いて、直接p値を求めることもできる:

この結果、棄却域は Z<1.959964 または Z>1.959964となるので、Zの値は棄却域に入る。 よって、結論「有意水準5%において、指導法データは心理学テスト(という母集団)から無作為抽出した標本とはいえない」。


なお、関数cdfを用いて、直接p値を求めることもできる:

In [47]:
 
st.norm.cdf(-2.828427)   # 下側確率
# 下側確率とすれば、p値は0.0023という小さな値(<  0.05)
Out[47]:
0.0023388684020295768
In [48]:
 
st.norm.cdf(2.828427)   # 上側確率
Out[48]:
0.99766113159797043
In [49]:
 
# 両側検定なので2倍する
2*st.norm.cdf(-2.828427)
# 両側検定であるから2倍したp値は0.0047という小さな値(< 0.05)
Out[49]:
0.0046777368040591535
In [3]:
 
help(st.norm.cdf)
Help on method cdf in module scipy.stats._distn_infrastructure:

cdf(self, x, *args, **kwds) method of scipy.stats._continuous_distns.norm_gen instance
    Cumulative distribution function of the given RV.
    
    Parameters
    ----------
    x : array_like
        quantiles
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    cdf : ndarray
        Cumulative distribution function evaluated at `x`

 
<h4>課題4-3</h4>
  標準正規分布のグラフを書き、有意水準5%の棄却域を黄色で表し、例題のZ値がどこに位置するかを重ね書きした図を作成せよ。
[ヒント]
前節<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Python-statistics3.html#RS03:normalDistribution">正規分布</A>の「正規分布グラフに領域を表示する関数」で紹介した関数を拡張修正して用いる。
Z値以下の領域をオレンジ色で表すと次のような図が得られる:<BR>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-E1.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-E1.png" width="30%" hspace="10" align="center" /></A>

課題4-3

標準正規分布のグラフを書き、有意水準5%の棄却域を黄色で表し、例題のZ値がどこに位置するかを重ね書きした図を作成せよ。

[ヒント] 前節正規分布の「正規分布グラフに領域を表示する関数」で紹介した関数を拡張修正して用いる。 Z値以下の領域をオレンジ色で表すと次のような図が得られる:

 
<h2><A NAME="RS04:tTest">t分布を用いた検定(1つの平均値の検定: 母分散が未知)</A></h2>
<div class="content">
正規母集団からの無作為標本であっても、母集団の分散σ<sup>2</sup>がわからない場合、
先の方法が使えません---先の検定で用いた検定統計量が計算できないからです
<BR>
そこで母分散の平方根σ の代わりに、標本から求められる不偏分散の平方根$\hat{\sigma}$を使い、
$$ t = \frac{\bar{X}-\mu}{\hat{\sigma}/\sqrt{n}} $$
を検定統計量とする。
これは<font class="word">自由度(df)</font> $n-1$<font class="word">t分布</font>に従う
<UL class="important">
<LI><font class="important">t分布</font>:統計学でよく利用される、正規分布の形に似た左右対称・山形の確率分布
<LI><font class="important">自由度(df)</font>:t分布の形状を決める
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-5.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-5.png" width="70%" hspace="0" align="right" />
</A>

t分布を用いた検定(1つの平均値の検定: 母分散が未知)

正規母集団からの無作為標本であっても、母集団の分散σ2がわからない場合、 先の方法が使えません---先の検定で用いた検定統計量が計算できないからです


そこで母分散の平方根σ の代わりに、標本から求められる不偏分散の平方根σ̂ を使い、

t=X¯μσ̂ /n
を検定統計量とする。 これは自由度(df) n1t分布に従う

  • t分布:統計学でよく利用される、正規分布の形に似た左右対称・山形の確率分布

  • 自由度(df):t分布の形状を決める

 
<H4>Pythonを使った実習</H4>
 <div class="content">
<B>例題:</B>「心理学テスト」が平均12の正規分布に従うものとする(<FONT class="other">分散は未知!</font>)。
<A HREF="#SampleData">前項にあげた「指導法データ」</A> (SampleData) が、この母集団から無作為抽出した標本と考えてよいかどうかを判定せよ
次のステップで行う:<OL class="other">
 <LI>帰無仮説と対立仮説をたてる: 帰無仮説は「無作為抽出した標本と考えて良い」、つまり$\mu=12$<BR>
     対立仮説は「無作為抽出した標本ではない」、つまり$\mu \neq 12$
 <LI>検定統計量の選択: <font class="rule">標本の不偏分散の平方根 $\hat{\sigma}$ を用い</font>
    $$ t = \frac{\bar{X}-\mu}{\hat{\sigma}/\sqrt{n}} $$を検定統計量とする
 <LI>有意水準の決定: 両側検定で、有意水準 5%、つまり$\alpha=0.05$
 <LI>検定統計量の実現値の計算:
 ```
 t =  (np.mean(SampleData) - 12) / (np.var(SampleData, ddof=1)/len(SampleData))**0.5             # 検定量
 ```
 <LI> 帰無仮説の棄却か採択かの決定: <B>帰無仮説によればこの検定統計量は自由度$df=n-1=19$のt分布に従う</B>
 ```
st.t.ppf(0.025,19)     # df=19、下側確率0.05/2 = 0.025となるtの値を求める(scipy.stats.tモジュールのppf関数
⇒ -2.0930240544082634
# 下側確率であるから、この値よりもt値が小さければ棄却される)
st.t.ppf(0.975,19)     # df=19、上側確率1-0.05/2 = 0.975となるtの値を求める(scipy.stats.tモジュールのppf関数
# 上側確率であるから、この値よりもt値が大きければ棄却される
⇒ 2.093024054408263
```
この結果、棄却域は $t < -2.093024$ または $t > 2.093024$となるので、tの値は棄却域に入る。<BR>
関数 cdf を用いて、直接p値を求めることもできる:
```
st.t.cdf(-2.616648,19)   # 下側確率 (scipy.stats.tモジュールのcdf関数)
# 下側確率とすれば、p値は0.0085 という小さな値(< 0.05)
 0.0084854604500293768
print(st.t.cdf(2.616648,19))  # 上側確率
 0.99151453955
print(2*st.t.cdf(-2.616648,19))    # 両側検定なので2倍する
 0.0169709209001
# 両側検定より2倍したp値は0.017という小さな値(< 0.05)
```
<LI>よって、結論「有意水準5%において、指導法データは心理学テスト(という母集団)から無作為抽出した標本とはいえない」。

Pythonを使った実習

例題:「心理学テスト」が平均12の正規分布に従うものとする(分散は未知!)。

前項にあげた「指導法データ」 (SampleData) が、この母集団から無作為抽出した標本と考えてよいかどうかを判定せよ

次のステップで行う:

  1. 帰無仮説と対立仮説をたてる: 帰無仮説は「無作為抽出した標本と考えて良い」、つまりμ=12
    対立仮説は「無作為抽出した標本ではない」、つまりμ12
  2. 検定統計量の選択: 標本の不偏分散の平方根 σ̂  を用い
    t=X¯μσ̂ /n
    を検定統計量とする
  3. 有意水準の決定: 両側検定で、有意水準 5%、つまりα=0.05
  4. 検定統計量の実現値の計算:
     t =  (np.mean(SampleData) - 12) / (np.var(SampleData, ddof=1)/len(SampleData))**0.5             # 検定量
    

  5. 帰無仮説の棄却か採択かの決定: 帰無仮説によればこの検定統計量は自由度df=n1=19のt分布に従う

    st.t.ppf(0.025,19)     # df=19、下側確率0.05/2 = 0.025となるtの値を求める(scipy.stats.tモジュールのppf関数
    ⇒ -2.0930240544082634
    # 下側確率であるから、この値よりもt値が小さければ棄却される)
    st.t.ppf(0.975,19)     # df=19、上側確率1-0.05/2 = 0.975となるtの値を求める(scipy.stats.tモジュールのppf関数
    # 上側確率であるから、この値よりもt値が大きければ棄却される
    ⇒ 2.093024054408263
    

    この結果、棄却域は t<2.093024 または t>2.093024となるので、tの値は棄却域に入る。
    関数 cdf を用いて、直接p値を求めることもできる:

    st.t.cdf(-2.616648,19)   # 下側確率 (scipy.stats.tモジュールのcdf関数)
    # 下側確率とすれば、p値は0.0085 という小さな値(< 0.05)
    ⇒ 0.0084854604500293768
    print(st.t.cdf(2.616648,19))  # 上側確率
    ⇒ 0.99151453955
    print(2*st.t.cdf(-2.616648,19))    # 両側検定なので2倍する
    ⇒ 0.0169709209001
    # 両側検定より2倍したp値は0.017という小さな値(< 0.05)
    

  6. よって、結論「有意水準5%において、指導法データは心理学テスト(という母集団)から無作為抽出した標本とはいえない」。

In [186]:
 
# 以上の実行
t =  (np.mean(SampleData) - 12) / (np.var(SampleData, ddof=1)/len(SampleData))**0.5             # 検定量
print("t = %f" % t)
t = -2.616648
In [187]:
 
import scipy.stats as st
st.t.ppf(0.025,19)     # df=19、下側確率0.05/2 = 0.025となるtの値を求める(scipy.stats.tモジュールのppf関数
# 下側確率であるから、この値よりもt値が小さければ棄却される)
Out[187]:
-2.0930240544082634
In [188]:
 
st.t.ppf(0.975,19)     # df=19、上側確率1-0.05/2 = 0.975となるtの値を求める(scipy.stats.tモジュールのppf関数
# 上側確率であるから、この値よりもt値が大きければ棄却される
Out[188]:
2.093024054408263
In [189]:
 
st.t.cdf(-2.616648,19)   # 下側確率 (scipy.stats.tモジュールのcdf関数)
# 下側確率とすれば、p値は0.0085 という小さな値(< 0.05)
Out[189]:
0.0084854604500293768
In [190]:
 
print(st.t.cdf(2.616648,19))  # 上側確率
print(2*st.t.cdf(-2.616648,19))    # 両側検定なので2倍する
# 両側検定より2倍したp値は0.017という小さな値(< 0.05)
0.99151453955
0.0169709209001
In [5]:
 
help(st.t.ppf)
Help on method ppf in module scipy.stats._distn_infrastructure:

ppf(self, q, *args, **kwds) method of scipy.stats._continuous_distns.t_gen instance
    Percent point function (inverse of `cdf`) at q of the given RV.
    
    Parameters
    ----------
    q : array_like
        lower tail probability
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    x : array_like
        quantile corresponding to the lower tail probability q.

In [6]:
 
help(st.t.cdf)
Help on method cdf in module scipy.stats._distn_infrastructure:

cdf(self, x, *args, **kwds) method of scipy.stats._continuous_distns.t_gen instance
    Cumulative distribution function of the given RV.
    
    Parameters
    ----------
    x : array_like
        quantiles
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    cdf : ndarray
        Cumulative distribution function evaluated at `x`

 
----
Pythonでt検定するための関数: 
以上のことをすべてやってくれる関数がscipy.statsモジュールのttest_1samp関数である。

Pythonでt検定するための関数:

以上のことをすべてやってくれる関数がscipy.statsモジュールのttest_1samp関数である。

In [70]:
 
import scipy.stats as st
help(st.ttest_1samp)
Help on function ttest_1samp in module scipy.stats.stats:

ttest_1samp(a, popmean, axis=0, nan_policy='propagate')
    Calculates the T-test for the mean of ONE group of scores.
    
    This is a two-sided test for the null hypothesis that the expected value
    (mean) of a sample of independent observations `a` is equal to the given
    population mean, `popmean`.
    
    Parameters
    ----------
    a : array_like
        sample observation
    popmean : float or array_like
        expected value in null hypothesis, if array_like than it must have the
        same shape as `a` excluding the axis dimension
    axis : int or None, optional
        Axis along which to compute test. If None, compute over the whole
        array `a`.
    nan_policy : {'propagate', 'raise', 'omit'}, optional
        Defines how to handle when input contains nan. 'propagate' returns nan,
        'raise' throws an error, 'omit' performs the calculations ignoring nan
        values. Default is 'propagate'.
    
    Returns
    -------
    statistic : float or array
        t-statistic
    pvalue : float or array
        two-tailed p-value
    
    Examples
    --------
    >>> from scipy import stats
    
    >>> np.random.seed(7654567)  # fix seed to get the same result
    >>> rvs = stats.norm.rvs(loc=5, scale=10, size=(50,2))
    
    Test if mean of random sample is equal to true mean, and different mean.
    We reject the null hypothesis in the second case and don't reject it in
    the first case.
    
    >>> stats.ttest_1samp(rvs,5.0)
    (array([-0.68014479, -0.04323899]), array([ 0.49961383,  0.96568674]))
    >>> stats.ttest_1samp(rvs,0.0)
    (array([ 2.77025808,  4.11038784]), array([ 0.00789095,  0.00014999]))
    
    Examples using axis and non-scalar dimension for population mean.
    
    >>> stats.ttest_1samp(rvs,[5.0,0.0])
    (array([-0.68014479,  4.11038784]), array([  4.99613833e-01,   1.49986458e-04]))
    >>> stats.ttest_1samp(rvs.T,[5.0,0.0],axis=1)
    (array([-0.68014479,  4.11038784]), array([  4.99613833e-01,   1.49986458e-04]))
    >>> stats.ttest_1samp(rvs,[[5.0],[0.0]])
    (array([[-0.68014479, -0.04323899],
           [ 2.77025808,  4.11038784]]), array([[  4.99613833e-01,   9.65686743e-01],
           [  7.89094663e-03,   1.49986458e-04]]))

 
「指導法データ」 (SampleData) を用いてその使い方を示す: <T>st.ttest_1samp(データ,$\mu$)</T>

「指導法データ」 (SampleData) を用いてその使い方を示す: st.ttest_1samp(データ,μ)

In [72]:
 
st.ttest_1samp(SampleData,12.0)
Out[72]:
Ttest_1sampResult(statistic=-2.616648017377738, pvalue=0.016970920269563441)
 
この表示から、t値が-2.617、p値が0.017 (両側検定)であることが得られる。

この表示から、t値が-2.617、p値が0.017 (両側検定)であることが得られる。

 
---
<h2><A NAME="RS04:correlationTest">相関係数の検定(無相関検定)</A></h2>
<div class="content">
<font class="important">無相関検定</font>:「母集団において相関が0である」と設定して行う検定
<BR>
母集団相関係数(母相関)に関する検定を行うときは、標本相関係数rから次を求めて検定統計量とする:
$$ t = \frac{r\sqrt{n-2}}{\sqrt{1-r^2}}$$

相関係数の検定(無相関検定)

無相関検定:「母集団において相関が0である」と設定して行う検定


母集団相関係数(母相関)に関する検定を行うときは、標本相関係数rから次を求めて検定統計量とする:

t=rn21r2

 
<H4>Pythonを使った実習</H4>
<div class="content">
<B>例題:</B>以下で与えられる「統計学テスト1」(StatTest1)と「統計学テスト2」(StatTest2)の得点の相関係数の検定を行え。有意水準は5%とする。

Pythonを使った実習

例題:以下で与えられる「統計学テスト1」(StatTest1)と「統計学テスト2」(StatTest2)の得点の相関係数の検定を行え。有意水準は5%とする。

In [73]:
 
import numpy as np
StatTest1 = np.array([6,10,6,10,5,3,5,9,3,3,11,6,11,9,7,5,8,7,7,9])
StatTest2 = np.array([10,13,8,15,8,6,9,10,7,3,18,14,18,11,12,5,7,12,7,7])
 
次のステップで行う:<OL>
<LI>帰無仮説と対立仮説をたてる: 帰無仮説は$\rho = 0$、つまり母相関$=0$<BR>
     対立仮説は「$\rho \neq 0$」、つまり母相関$\neq 0$
 <LI>検定統計量の選択: $t = \frac{r\sqrt{n-2}}{\sqrt{1-r^2}}$
 <LI>有意水準の決定: 両側検定で、有意水準 5%、つまり$\alpha = 0.05$
 <LI>検定統計量の実現値の計算:
```
SampleCorr = np.corrcoef(StatTest1, StatTest2)[0,1]
print("Sample Correlation = %f" % SampleCorr)  # 標本相関
⇒ Sample Correlation = 0.749659
SampleSize = len(StatTest1)
tDividend =  SampleCorr * (SampleSize - 2.0)**0.5
tDivider = (1.0 - SampleCorr**2)**0.5
t = tDividend/tDivider
⇒ t statistics = 4.805707
```
<LI>帰無仮説の棄却か採択かの決定: <B>帰無仮説によればこの検定統計量は自由度$df=n-2=18$のt分布に従う</B>
```
import scipy.stats as st
print(st.t.ppf(0.025,18))    # df=18、下側確率0.05/2 = 0.025となるtの値を求める
⇒ -2.10092204024
# 下側確率であるから、この値よりもt値が小さければ棄却される
print(st.t.ppf(0.975,18))    # df=18、上側確率 1 - 0.05/2 = 0.975となるtの値を求める
⇒ 2.10092204024
# 上側確率であるから、この値よりもt値が大きければ棄却される
```
<LI>この結果、棄却域は $t < -2.101$ または $t > 2.101$となるので、tの値は棄却域に入る。
よって結論「統計学テスト1(StatTest1)と統計学テスト2(StatTest2)は有意水準5%において強い相関(相関係数0.75)がある」。
</OL>

次のステップで行う:

  1. 帰無仮説と対立仮説をたてる: 帰無仮説はρ=0、つまり母相関=0
    対立仮説は「ρ0」、つまり母相関0
  2. 検定統計量の選択: t=rn21r2
  3. 有意水準の決定: 両側検定で、有意水準 5%、つまりα=0.05
  4. 検定統計量の実現値の計算:

    SampleCorr = np.corrcoef(StatTest1, StatTest2)[0,1]
    print("Sample Correlation = %f" % SampleCorr)  # 標本相関
    ⇒ Sample Correlation = 0.749659
    SampleSize = len(StatTest1)
    tDividend =  SampleCorr * (SampleSize - 2.0)**0.5
    tDivider = (1.0 - SampleCorr**2)**0.5
    t = tDividend/tDivider
    ⇒ t statistics = 4.805707
    

  5. 帰無仮説の棄却か採択かの決定: 帰無仮説によればこの検定統計量は自由度df=n2=18のt分布に従う

    import scipy.stats as st
    print(st.t.ppf(0.025,18))    # df=18、下側確率0.05/2 = 0.025となるtの値を求める
    ⇒ -2.10092204024
    # 下側確率であるから、この値よりもt値が小さければ棄却される
    print(st.t.ppf(0.975,18))    # df=18、上側確率 1 - 0.05/2 = 0.975となるtの値を求める
    ⇒ 2.10092204024
    # 上側確率であるから、この値よりもt値が大きければ棄却される
    

  6. この結果、棄却域は t<2.101 または t>2.101となるので、tの値は棄却域に入る。 よって結論「統計学テスト1(StatTest1)と統計学テスト2(StatTest2)は有意水準5%において強い相関(相関係数0.75)がある」。

In [83]:
 
from __future__ import division
import numpy as np
SampleCorr = np.corrcoef(StatTest1, StatTest2)[0,1]
print("Sample Correlation = %f" % SampleCorr)  # 標本相関
SampleSize = len(StatTest1)
tDividend =  SampleCorr * (SampleSize - 2.0)**0.5
tDivider = (1.0 - SampleCorr**2)**0.5
t = tDividend/tDivider
print("t statistics = %f" % t)
Sample Correlation = 0.749659
t statistics = 4.805707
In [86]:
 
import scipy.stats as st
print(st.t.ppf(0.025,18))    # df=18、下側確率0.05/2 = 0.025となるtの値を求める
# 下側確率であるから、この値よりもt値が小さければ棄却される
print(st.t.ppf(0.975,18))    # df=18、上側確率 1 - 0.05/2 = 0.975となるtの値を求める
# 上側確率であるから、この値よりもt値が大きければ棄却される
-2.10092204024
2.10092204024
 
なおscipy.stats.tモジュールのcdf関数を用いて、直接p値を求めることもできる:

なおscipy.stats.tモジュールのcdf関数を用いて、直接p値を求めることもできる:

In [93]:
 
print(st.t.cdf(t,18))  # 上側確率
print( (1.0  - st.t.cdf(t,18))*2.0 )  # 両側検定なので2倍する
# 両側検定により2倍したp値は0.00014という小さな値(< 0.05)
0.999929188559
0.000141622881555
 
---
Pythonで無相関検定するための関数: scipy.stat.pearsonr

Pythonで無相関検定するための関数: scipy.stat.pearsonr

In [92]:
 
help(st.pearsonr)
Help on function pearsonr in module scipy.stats.stats:

pearsonr(x, y)
    Calculates a Pearson correlation coefficient and the p-value for testing
    non-correlation.
    
    The Pearson correlation coefficient measures the linear relationship
    between two datasets. Strictly speaking, Pearson's correlation requires
    that each dataset be normally distributed, and not necessarily zero-mean.
    Like other correlation coefficients, this one varies between -1 and +1
    with 0 implying no correlation. Correlations of -1 or +1 imply an exact
    linear relationship. Positive correlations imply that as x increases, so
    does y. Negative correlations imply that as x increases, y decreases.
    
    The p-value roughly indicates the probability of an uncorrelated system
    producing datasets that have a Pearson correlation at least as extreme
    as the one computed from these datasets. The p-values are not entirely
    reliable but are probably reasonable for datasets larger than 500 or so.
    
    Parameters
    ----------
    x : (N,) array_like
        Input
    y : (N,) array_like
        Input
    
    Returns
    -------
    r : float
        Pearson's correlation coefficient
    p-value : float
        2-tailed p-value
    
    References
    ----------
    http://www.statsoft.com/textbook/glosp.html#Pearson%20Correlation

In [95]:
 
import scipy.stats as st
SampleCorr = st.pearsonr(StatTest1, StatTest2)
print(SampleCorr)
(0.74965896482424488, 0.00014162288155448246)
 
pearsonr関数の出力の第一要素は標本相関係数(0.75)、第二要素は両側検定によるp値である。

pearsonr関数の出力の第一要素は標本相関係数(0.75)、第二要素は両側検定によるp値である。

xxxxxxxxxx
---
<H2><A NAME="RS04:chisquareTest">独立性の検定(カイ自乗検定)</A></H2>
<div class="content">
2つの質的変数が独立かどうかを確かめる---
<font class="word">独立</font>とは「2つの質的変数に連関がない」こと
<ul class="important">
 <LI><font class="word">独立性の検定</font>:2つの質的変数間の連関の有意性を調べる検定<BR>
  <LI><font class="word">期待度数</font>:2つの変数の間に連関がない(独立である)という帰無仮説のもとで、帰無仮説が正しければ(連関がなければ)これくらいの度数をとるだろうと期待される度数<BR>
    クロス集計表におけるセルの期待度数 = (セルが属する行の周辺度数 ×セルが属する列の周辺度数)÷総度数
<LI><font class="word">$\chi^2$(カイ2乗)</font>という確率分布を利用するため、カイ自乗(2乗)検定ともいう。
   <LI><font class="word">独立性の検定における検定統計量の式</font> 
$$ \chi^2 = \frac{(O_1-E_1)^2}{E_1}+\frac{(O_2-E_2)^2}{E_2}+\cdots+\frac{(O_k-E_k)^2}{E_k}$$
<BR>
   $O_1 \sim O_k$は観測度数、$E_1 \sim E_k$は期待度数
  <LI>カイ二乗分布: <CENTER>
<A HREF="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-8.png" target="_blank">
<IMG  class="fig" SRC="http://lang.sist.chukyo-u.ac.jp/Classes/PythonProbStat/Figs/RS-4-8.png" align="top" width="40%" hspace="5" align="right"  /></A>
</CENTER>
</ul>

独立性の検定(カイ自乗検定)

2つの質的変数が独立かどうかを確かめる---

独立とは「2つの質的変数に連関がない」こと

  • 独立性の検定:2つの質的変数間の連関の有意性を調べる検定
  • 期待度数:2つの変数の間に連関がない(独立である)という帰無仮説のもとで、帰無仮説が正しければ(連関がなければ)これくらいの度数をとるだろうと期待される度数
    クロス集計表におけるセルの期待度数 = (セルが属する行の周辺度数 ×セルが属する列の周辺度数)÷総度数
  • χ2(カイ2乗)という確率分布を利用するため、カイ自乗(2乗)検定ともいう。
  • 独立性の検定における検定統計量の式
    χ2=(O1E1)2E1+(O2E2)2E2++(OkEk)2Ek

    O1Okは観測度数、E1Ekは期待度数
  • カイ二乗分布:

 
<H4>Pythonを使った実習</H4>
<div class="content">
<B>例題:</B>20名の学生に対し数学(Math)と統計学(Stat)の好き嫌いをアンケート調査した結果が以下。このことから、一般に数学と統計学の好き嫌いの間に有意な連関があるといえるかどうか、有意水準5%で検定せよ。
```
Math = np.array(["嫌い","嫌い","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","好き","好き",
    "嫌い","好き","嫌い","嫌い","好き","嫌い","嫌い","嫌い","嫌い"])
Stat = np.array(["好き","好き","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","嫌い","好き",
    "好き","好き","嫌い","好き","嫌い","嫌い","嫌い","嫌い","嫌い"])
```
このクロス集計表は以下:
<CENTER>
<table border=4 width=600 align=center>
  <tr bgcolor="#cccccc">
    <th width=25%>  </th>
    <th width=25%>統計学好き</th>
    <th width=25%>統計学嫌い</th>
    <th width=25%>合計</th>
  </tr>
  <tr>
   <TD align="center">数学好き</TD>   <TD>Expectation_11</TD>     <TD>Expectation_12</TD>    <TD align="center">6</TD>
  </tr>
    <tr>
   <TD align="center">数学嫌い</TD>   <TD>Expectation_21</TD>     <TD>Expectation_22</TD>    <TD align="center">14</TD>
  </tr>
  <tr>
   <TD  align="center"></TD>   <TD align="center">8</TD>     <TD align="center">12</TD>    <TD align="center">20</TD>
  </tr>
 </TABLE>
</CENTER>
マスのことを<font class="word">セル</font>、セルに書かれた数値を<font class="word">観測度数</font>
観測度数を各々、列・行で合計したものを<font class="word">周辺度数</font>、周辺度数の合計を<font class="word">総度数</font>と呼ぶ
<BR>
<font class="rule">自由度 df = (行の数-1)×(列の数-1)</font>

Pythonを使った実習

例題:20名の学生に対し数学(Math)と統計学(Stat)の好き嫌いをアンケート調査した結果が以下。このことから、一般に数学と統計学の好き嫌いの間に有意な連関があるといえるかどうか、有意水準5%で検定せよ。

Math = np.array(["嫌い","嫌い","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","好き","好き",
    "嫌い","好き","嫌い","嫌い","好き","嫌い","嫌い","嫌い","嫌い"])
Stat = np.array(["好き","好き","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","嫌い","好き",
    "好き","好き","嫌い","好き","嫌い","嫌い","嫌い","嫌い","嫌い"])

このクロス集計表は以下:

統計学好き 統計学嫌い 合計
数学好き Expectation_11 Expectation_12 6
数学嫌い Expectation_21 Expectation_22 14
8 12 20
マスのことをセル、セルに書かれた数値を観測度数、 観測度数を各々、列・行で合計したものを周辺度数、周辺度数の合計を総度数と呼ぶ


自由度 df = (行の数-1)×(列の数-1)

In [144]:
 
Math = np.array(["嫌い","嫌い","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","好き","好き",
    "嫌い","好き","嫌い","嫌い","好き","嫌い","嫌い","嫌い","嫌い"])
Stat = np.array(["好き","好き","好き","好き","嫌い","嫌い","嫌い","嫌い","嫌い","嫌い","好き",
    "好き","好き","嫌い","好き","嫌い","嫌い","嫌い","嫌い","嫌い"])
import pandas as pd
data = pd.DataFrame({ 'Stat':Stat, 'Math':Math})
table = pd.crosstab(data.Math,data.Stat,margins=True)  # クロス集計表を作る
table
Out[144]:
Stat 好き 嫌い All
Math
好き 4 2 6
嫌い 4 10 14
All 8 12 20
 
次のステップで行う:<OL class="other">
<LI>帰無仮説と対立仮説をたてる: 帰無仮説H<sub>0</sub>は、「数学と統計学の2つの変数は独立(連関なし)」<BR>
     対立仮説H<sub>1</sub>は「「数学と統計学の2つの変数は独立(連関なし)」
 <LI>検定統計量の選択: $$ \chi^2 = \frac{(O_1-E_1)^2}{E_1}+\frac{(O_2-E_2)^2}{E_2}+\cdots+\frac{(O_k-E_k)^2}{E_k}$$
 <LI>有意水準の決定: 5%とする(片側検定---カイ二乗検定は棄却域が一つしかない)
 <LI>検定統計量の実現値の計算:

次のステップで行う:

  1. 帰無仮説と対立仮説をたてる: 帰無仮説H0は、「数学と統計学の2つの変数は独立(連関なし)」
    対立仮説H1は「「数学と統計学の2つの変数は独立(連関なし)」
  2. 検定統計量の選択:
    χ2=(O1E1)2E1+(O2E2)2E2++(OkEk)2Ek
  3. 有意水準の決定: 5%とする(片側検定---カイ二乗検定は棄却域が一つしかない)
  4. 検定統計量の実現値の計算:

In [145]:
 
Expectaion_11 = 8*6/20.0
Expectaion_12 = 12*6/20.0
Expectaion_21 = 8*14/20.0
Expectaion_22 = 12*14/20.0
ExpectedFrequency = np.array([Expectaion_11, Expectaion_21,Expectaion_12,Expectaion_22])
ObservedFrequency = np.array([4,4,2,10])
ChiSqElements = (ObservedFrequency - ExpectedFrequency)**2 / ExpectedFrequency
ChiSq = np.sum(ChiSqElements)   #  検定統計量
print(ChiSq)
2.53968253968
 
<BR>5. 帰無仮説の棄却か採択かの決定: <B>帰無仮説によればこの検定統計量は自由度$df=1$$\chi^2$分布に従う</B>


5. 帰無仮説の棄却か採択かの決定: 帰無仮説によればこの検定統計量は自由度df=1χ2分布に従う

In [104]:
 
import scipy.stats as st
help(st.distributions.chi2.ppf)
Help on method ppf in module scipy.stats._distn_infrastructure:

ppf(self, q, *args, **kwds) method of scipy.stats._continuous_distns.chi2_gen instance
    Percent point function (inverse of `cdf`) at q of the given RV.
    
    Parameters
    ----------
    q : array_like
        lower tail probability
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    x : array_like
        quantile corresponding to the lower tail probability q.

In [8]:
 
st.distributions.chi2.ppf(0.95,1)  # 自由度1のカイ二乗分布で確率0.95となるχ2の値を求める
# これが棄却域を定める---この値よりもχ2値が大きければ棄却される
# カイ二乗分布は『上側』のみ
Out[8]:
3.8414588206941236
 
この結果、棄却域は  $\chi^2 >  3.84$ となり、この例題における$\chi^2$の値(=2.54)は棄却域に入っていない、つまり帰無仮説は棄却されず、採択される。
よって結論「有意水準5%において、数学と統計学の2つの変数は<B>独立ではない(連関がない)</B>」。
<BR>
なお、cdf関数を用いて、直接p値を求めることもできる:

この結果、棄却域は χ2>3.84 となり、この例題におけるχ2の値(=2.54)は棄却域に入っていない、つまり帰無仮説は棄却されず、採択される。 よって結論「有意水準5%において、数学と統計学の2つの変数は独立ではない(連関がない)」。


なお、cdf関数を用いて、直接p値を求めることもできる:

In [107]:
 
help(st.distributions.chi2.cdf)
Help on method cdf in module scipy.stats._distn_infrastructure:

cdf(self, x, *args, **kwds) method of scipy.stats._continuous_distns.chi2_gen instance
    Cumulative distribution function of the given RV.
    
    Parameters
    ----------
    x : array_like
        quantiles
    arg1, arg2, arg3,... : array_like
        The shape parameter(s) for the distribution (see docstring of the
        instance object for more information)
    loc : array_like, optional
        location parameter (default=0)
    scale : array_like, optional
        scale parameter (default=1)
    
    Returns
    -------
    cdf : ndarray
        Cumulative distribution function evaluated at `x`

In [108]:
 
1.0 - st.distributions.chi2.cdf(ChiSq,1)  # 上側確率
# p値が有意水準0.05よりも大きいので、帰無仮説は棄却されない
Out[108]:
0.11101710551218613
 
---
Pythonでカイ二乗検定するための関数<font class="rule">chisquare</font>(scipy.statsモジュール):  

Pythonでカイ二乗検定するための関数chisquare(scipy.statsモジュール):

In [109]:
 
help(st.chisquare)
Help on function chisquare in module scipy.stats.stats:

chisquare(f_obs, f_exp=None, ddof=0, axis=0)
    Calculates a one-way chi square test.
    
    The chi square test tests the null hypothesis that the categorical data
    has the given frequencies.
    
    Parameters
    ----------
    f_obs : array_like
        Observed frequencies in each category.
    f_exp : array_like, optional
        Expected frequencies in each category.  By default the categories are
        assumed to be equally likely.
    ddof : int, optional
        "Delta degrees of freedom": adjustment to the degrees of freedom
        for the p-value.  The p-value is computed using a chi-squared
        distribution with ``k - 1 - ddof`` degrees of freedom, where `k`
        is the number of observed frequencies.  The default value of `ddof`
        is 0.
    axis : int or None, optional
        The axis of the broadcast result of `f_obs` and `f_exp` along which to
        apply the test.  If axis is None, all values in `f_obs` are treated
        as a single data set.  Default is 0.
    
    Returns
    -------
    chisq : float or ndarray
        The chi-squared test statistic.  The value is a float if `axis` is
        None or `f_obs` and `f_exp` are 1-D.
    p : float or ndarray
        The p-value of the test.  The value is a float if `ddof` and the
        return value `chisq` are scalars.
    
    See Also
    --------
    power_divergence
    mstats.chisquare
    
    Notes
    -----
    This test is invalid when the observed or expected frequencies in each
    category are too small.  A typical rule is that all of the observed
    and expected frequencies should be at least 5.
    
    The default degrees of freedom, k-1, are for the case when no parameters
    of the distribution are estimated. If p parameters are estimated by
    efficient maximum likelihood then the correct degrees of freedom are
    k-1-p. If the parameters are estimated in a different way, then the
    dof can be between k-1-p and k-1. However, it is also possible that
    the asymptotic distribution is not a chisquare, in which case this
    test is not appropriate.
    
    References
    ----------
    .. [1] Lowry, Richard.  "Concepts and Applications of Inferential
           Statistics". Chapter 8. http://faculty.vassar.edu/lowry/ch8pt1.html
    .. [2] "Chi-squared test", http://en.wikipedia.org/wiki/Chi-squared_test
    
    Examples
    --------
    When just `f_obs` is given, it is assumed that the expected frequencies
    are uniform and given by the mean of the observed frequencies.
    
    >>> from scipy.stats import chisquare
    >>> chisquare([16, 18, 16, 14, 12, 12])
    (2.0, 0.84914503608460956)
    
    With `f_exp` the expected frequencies can be given.
    
    >>> chisquare([16, 18, 16, 14, 12, 12], f_exp=[16, 16, 16, 16, 16, 8])
    (3.5, 0.62338762774958223)
    
    When `f_obs` is 2-D, by default the test is applied to each column.
    
    >>> obs = np.array([[16, 18, 16, 14, 12, 12], [32, 24, 16, 28, 20, 24]]).T
    >>> obs.shape
    (6, 2)
    >>> chisquare(obs)
    (array([ 2.        ,  6.66666667]), array([ 0.84914504,  0.24663415]))
    
    By setting ``axis=None``, the test is applied to all data in the array,
    which is equivalent to applying the test to the flattened array.
    
    >>> chisquare(obs, axis=None)
    (23.31034482758621, 0.015975692534127565)
    >>> chisquare(obs.ravel())
    (23.31034482758621, 0.015975692534127565)
    
    `ddof` is the change to make to the default degrees of freedom.
    
    >>> chisquare([16, 18, 16, 14, 12, 12], ddof=1)
    (2.0, 0.73575888234288467)
    
    The calculation of the p-values is done by broadcasting the
    chi-squared statistic with `ddof`.
    
    >>> chisquare([16, 18, 16, 14, 12, 12], ddof=[0,1,2])
    (2.0, array([ 0.84914504,  0.73575888,  0.5724067 ]))
    
    `f_obs` and `f_exp` are also broadcast.  In the following, `f_obs` has
    shape (6,) and `f_exp` has shape (2, 6), so the result of broadcasting
    `f_obs` and `f_exp` has shape (2, 6).  To compute the desired chi-squared
    statistics, we use ``axis=1``:
    
    >>> chisquare([16, 18, 16, 14, 12, 12],
    ...           f_exp=[[16, 16, 16, 16, 16, 8], [8, 20, 20, 16, 12, 12]],
    ...           axis=1)
    (array([ 3.5 ,  9.25]), array([ 0.62338763,  0.09949846]))

In [146]:
 
ExpectedFrequency = np.array([Expectaion_11, Expectaion_21,Expectaion_12,Expectaion_22])
ObservedFrequency = np.array([4,4,2,10])
st.chisquare(ObservedFrequency, f_exp =ExpectedFrequency, ddof=2) # 自由度の計算に使うddofの値に注意
Out[146]:
Power_divergenceResult(statistic=2.5396825396825395, pvalue=0.1110171055121861)
 
<H2><A NAME="RS04:sampleSize">サンプルサイズの影響</A></h2>
<div class="content">
標本における連関の大きさが全く同じであっても、サンプルサイズが異なると検定の結果が変わることがある<BR>
サンプルサイズが大きくなると、有意になりやすい---統計的仮説検定一般にいえる性質

サンプルサイズの影響

標本における連関の大きさが全く同じであっても、サンプルサイズが異なると検定の結果が変わることがある
サンプルサイズが大きくなると、有意になりやすい---統計的仮説検定一般にいえる性質

In [160]:
 
import pandas as pd
data = { 'Mastered': [16,12], 'NotMastered':[4,8]}   # ある科目の履修 vs 未履修
df = pd.DataFrame(data)
df.index=['Humanities','Technicals']  # 文系 vs 理系
df # クロス集計表
Out[160]:
Mastered NotMastered
Humanities 16 4
Technicals 12 8
In [165]:
 
Exp_11 = sum(df['Mastered'])*sum(df.loc['Humanities',:])/40.0
Exp_12 = sum(df['NotMastered'])*sum(df.loc['Humanities',:])/40.0
Exp_21 = sum(df['Mastered'])*sum(df.loc['Technicals',:])/40.0
Exp_22 = sum(df['NotMastered'])*sum(df.loc['Technicals',:])/40.0
In [166]:
 
ExpectedFrequency = np.array([Exp_11, Exp_12,Exp_21,Exp_22])
ObservedFrequency = np.array([16,4,12,8])
st.chisquare(ObservedFrequency, f_exp =ExpectedFrequency, ddof=2) # 自由度の計算に使うddofの値に注意
# p値=0.17なので帰無仮説は棄却されない → 連関なし
Out[166]:
Power_divergenceResult(statistic=1.9047619047619047, pvalue=0.16754627748861739)
In [167]:
 
data10 = { 'Mastered': [160,120], 'NotMastered':[40,80]}   # ある科目の履修 vs 未履修 --- 前の10倍
df = pd.DataFrame(data10)
df.index=['Humanities','Technicals']  # 文系 vs 理系
df # クロス集計表
Out[167]:
Mastered NotMastered
Humanities 160 40
Technicals 120 80
In [168]:
 
Exp_11 = sum(df['Mastered'])*sum(df.loc['Humanities',:])/400.0
Exp_12 = sum(df['NotMastered'])*sum(df.loc['Humanities',:])/400.0
Exp_21 = sum(df['Mastered'])*sum(df.loc['Technicals',:])/400.0
Exp_22 = sum(df['NotMastered'])*sum(df.loc['Technicals',:])/400.0
In [169]:
 
ExpectedFrequency = np.array([Exp_11, Exp_12,Exp_21,Exp_22])
ObservedFrequency = np.array([160,40,120,80])
st.chisquare(ObservedFrequency, f_exp =ExpectedFrequency, ddof=2) # 自由度の計算に使うddofの値に注意
# p値=0.000013なので帰無仮説は棄却される→連関あり
Out[169]:
Power_divergenceResult(statistic=19.047619047619047, pvalue=1.2749674921097076e-05)
 
---
<H2><A NAME="Functions">関数のまとめ</A></H2>
注:  numpyをnp, np.randomを random、 matplotlib.pyplotをplt、pandasをpd、scipy.statsをstと略記する
|目的 | 関数名とモジュール       | 使い方|
|:---|:---|:---|
|指定された範囲からランダム抽出|random.choice(配列,個数)| random.choice(range(10),5)|
|標準正規分布で下側確率に対応する確率分布関数の値 | st.norm.ppf(p) | st.norm.ppf(0.025) # $Prob(Z < q)=0.025$となるqの値|
|標準正規分布で下側確率(p値)を求める|st.norm.cdf(z)|st.norm.cdf(1.96) # $Prob(Z < 1.96)$の値(p値)|
|t分布で下側確率に対応する確率分布関数の値 | st.t.ppf(p,自由度) | st.t.ppf(0.025,19) # 自由度19のt分布で$Prob(Z < q)=0.025$となるqの値|
|t分布で下側確率(p値)を求める|st.t.cdf(z,df)|st.t.cdf(1.96,19) # 自由度19のt分布で$Prob(Z < 1.96)$の値(p値)|
|t検定を行う|ttest_1samp(データ,$\mu$)|ttest_1samp(np.array([13,14,7,12,10,6,8,15,4,14,9,6,8,8,12,15]),12.0) # $\mu=12$の検定 |
|無相関検定を行う|st.pearsonr(データ1,データ2)|st.pearsonr(StatTest1, StatTest2) # 出力の第一要素は標本相関係数、第二要素は両側検定によるp値|
|カイ二乗分布の確率密度関数|st.distributions.chi2.pdf(x,自由度) | plt.plot(x,st.distributions.chi2.pdf(x,3)) \# 自由度3の$\chi^2$分布関数の描画|
|カイ二乗分布で上側確率に対応する値を求める|st.distributions.chi2.ppf(p,自由度)|st.distributions.chi2.ppf(0.95,2) # 自由度2で$Prob(Z < q)=0.95$となるq値を求める|
|カイ二乗分布で上側確率を求める|1-st.distributions.chi2.cdf(z,自由度)|1-st.distributions.chi2.cdf(3.5,1) # 自由度1のカイ2乗分布で$Prob(Z \geq 3.5)$となる確率|
|カイ二乗検定を行う(独立性の検定)|st.chisquare(観測度数リスト, f_exp =期待度数リスト, ddof=n) #n=観測個数-1-自由度 |st.chisquare(ObservedFrequency, f_exp =ExpectedFrequency, ddof=2) |

関数のまとめ

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

目的 関数名とモジュール 使い方
指定された範囲からランダム抽出 random.choice(配列,個数) random.choice(range(10),5)
標準正規分布で下側確率に対応する確率分布関数の値 st.norm.ppf(p) st.norm.ppf(0.025) # Prob(Z<q)=0.025となるqの値
標準正規分布で下側確率(p値)を求める st.norm.cdf(z) st.norm.cdf(1.96) # Prob(Z<1.96)の値(p値)
t分布で下側確率に対応する確率分布関数の値 st.t.ppf(p,自由度) st.t.ppf(0.025,19) # 自由度19のt分布でProb(Z<q)=0.025となるqの値
t分布で下側確率(p値)を求める st.t.cdf(z,df) st.t.cdf(1.96,19) # 自由度19のt分布でProb(Z<1.96)の値(p値)
t検定を行う ttest_1samp(データ,μ) ttest_1samp(np.array([13,14,7,12,10,6,8,15,4,14,9,6,8,8,12,15]),12.0) # μ=12の検定
無相関検定を行う st.pearsonr(データ1,データ2) st.pearsonr(StatTest1, StatTest2) # 出力の第一要素は標本相関係数、第二要素は両側検定によるp値
カイ二乗分布の確率密度関数 st.distributions.chi2.pdf(x,自由度) plt.plot(x,st.distributions.chi2.pdf(x,3)) # 自由度3のχ2分布関数の描画
カイ二乗分布で上側確率に対応する値を求める st.distributions.chi2.ppf(p,自由度) st.distributions.chi2.ppf(0.95,2) # 自由度2でProb(Z<q)=0.95となるq値を求める
カイ二乗分布で上側確率を求める 1-st.distributions.chi2.cdf(z,自由度) 1-st.distributions.chi2.cdf(3.5,1) # 自由度1のカイ2乗分布でProb(Z3.5)となる確率
カイ二乗検定を行う(独立性の検定) st.chisquare(観測度数リスト, f_exp =期待度数リスト, ddof=n) #n=観測個数-1-自由度 st.chisquare(ObservedFrequency, f_exp =ExpectedFrequency, ddof=2)
 
<H2><A NAME="Exercises">演習問題4</A></H2>
<h4>演習問題4-1</H4>
次のデータ(単位はcm)は、平均170cmの正規分布に従う20歳男性の母集団からの無作為抽出と考えてよいかどうかを検定せよ。

演習問題4

演習問題4-1

次のデータ(単位はcm)は、平均170cmの正規分布に従う20歳男性の母集団からの無作為抽出と考えてよいかどうかを検定せよ。

In [171]:
 
import numpy as np
Height = np.array([165,150,170,168,159,170,167,178,155,159,161,162,166,171,155,160,168,172,155,167])
 
<h4>演習問題4-2</H4>
以下に示すデータにおいて、勉強時間(StudyHours)と定期試験の成績(ExamResult)の相関係数の無相関検定を行え。

演習問題4-2

以下に示すデータにおいて、勉強時間(StudyHours)と定期試験の成績(ExamResult)の相関係数の無相関検定を行え。

In [170]:
 
import numpy as np
StudyHours = np.array([1, 3, 10, 12, 6, 3, 8, 4, 1, 5])
ExamResult = np.array([20, 40, 100, 80, 50, 50, 70, 50, 10, 60])
 
<h4>演習問題4-3</H4>
先の演習問題4-2のデータに対し、ピアソンの相関係数とスピアマンの順位相関係数を求め、
さらに無相関検定も行え。<BR>

演習問題4-3

先の演習問題4-2のデータに対し、ピアソンの相関係数とスピアマンの順位相関係数を求め、 さらに無相関検定も行え。

In [174]:
 
import scipy.stats as st
help(st.spearmanr)
Help on function spearmanr in module scipy.stats.stats:

spearmanr(a, b=None, axis=0, nan_policy='propagate')
    Calculates a Spearman rank-order correlation coefficient and the p-value
    to test for non-correlation.
    
    The Spearman correlation is a nonparametric measure of the monotonicity
    of the relationship between two datasets. Unlike the Pearson correlation,
    the Spearman correlation does not assume that both datasets are normally
    distributed. Like other correlation coefficients, this one varies
    between -1 and +1 with 0 implying no correlation. Correlations of -1 or
    +1 imply an exact monotonic relationship. Positive correlations imply that
    as x increases, so does y. Negative correlations imply that as x
    increases, y decreases.
    
    The p-value roughly indicates the probability of an uncorrelated system
    producing datasets that have a Spearman correlation at least as extreme
    as the one computed from these datasets. The p-values are not entirely
    reliable but are probably reasonable for datasets larger than 500 or so.
    
    Parameters
    ----------
    a, b : 1D or 2D array_like, b is optional
        One or two 1-D or 2-D arrays containing multiple variables and
        observations. When these are 1-D, each represents a vector of
        observations of a single variable. For the behavior in the 2-D case,
        see under ``axis``, below.
        Both arrays need to have the same length in the ``axis`` dimension.
    axis : int or None, optional
        If axis=0 (default), then each column represents a variable, with
        observations in the rows. If axis=1, the relationship is transposed:
        each row represents a variable, while the columns contain observations.
        If axis=None, then both arrays will be raveled.
    nan_policy : {'propagate', 'raise', 'omit'}, optional
        Defines how to handle when input contains nan. 'propagate' returns nan,
        'raise' throws an error, 'omit' performs the calculations ignoring nan
        values. Default is 'propagate'.
    
    Returns
    -------
    correlation : float or ndarray (2-D square)
        Spearman correlation matrix or correlation coefficient (if only 2
        variables are given as parameters. Correlation matrix is square with
        length equal to total number of variables (columns or rows) in a and b
        combined.
    pvalue : float
        The two-sided p-value for a hypothesis test whose null hypothesis is
        that two sets of data are uncorrelated, has same dimension as rho.
    
    Notes
    -----
    Changes in scipy 0.8.0: rewrite to add tie-handling, and axis.
    
    References
    ----------
    
    .. [1] Zwillinger, D. and Kokoska, S. (2000). CRC Standard
       Probability and Statistics Tables and Formulae. Chapman & Hall: New
       York. 2000.
       Section  14.7
    
    Examples
    --------
    >>> from scipy import stats
    >>> stats.spearmanr([1,2,3,4,5], [5,6,7,8,7])
    (0.82078268166812329, 0.088587005313543798)
    >>> np.random.seed(1234321)
    >>> x2n = np.random.randn(100, 2)
    >>> y2n = np.random.randn(100, 2)
    >>> stats.spearmanr(x2n)
    (0.059969996999699973, 0.55338590803773591)
    >>> stats.spearmanr(x2n[:,0], x2n[:,1])
    (0.059969996999699973, 0.55338590803773591)
    >>> rho, pval = stats.spearmanr(x2n, y2n)
    >>> rho
    array([[ 1.        ,  0.05997   ,  0.18569457,  0.06258626],
           [ 0.05997   ,  1.        ,  0.110003  ,  0.02534653],
           [ 0.18569457,  0.110003  ,  1.        ,  0.03488749],
           [ 0.06258626,  0.02534653,  0.03488749,  1.        ]])
    >>> pval
    array([[ 0.        ,  0.55338591,  0.06435364,  0.53617935],
           [ 0.55338591,  0.        ,  0.27592895,  0.80234077],
           [ 0.06435364,  0.27592895,  0.        ,  0.73039992],
           [ 0.53617935,  0.80234077,  0.73039992,  0.        ]])
    >>> rho, pval = stats.spearmanr(x2n.T, y2n.T, axis=1)
    >>> rho
    array([[ 1.        ,  0.05997   ,  0.18569457,  0.06258626],
           [ 0.05997   ,  1.        ,  0.110003  ,  0.02534653],
           [ 0.18569457,  0.110003  ,  1.        ,  0.03488749],
           [ 0.06258626,  0.02534653,  0.03488749,  1.        ]])
    >>> stats.spearmanr(x2n, y2n, axis=None)
    (0.10816770419260482, 0.1273562188027364)
    >>> stats.spearmanr(x2n.ravel(), y2n.ravel())
    (0.10816770419260482, 0.1273562188027364)
    
    >>> xint = np.random.randint(10, size=(100, 2))
    >>> stats.spearmanr(xint)
    (0.052760927029710199, 0.60213045837062351)

In [ ]:
 
 
<h4>演習問題4-4</H4>
以下に示す演習問題2-2のデータに対し、カイ二乗検定を行え。

演習問題4-4

以下に示す演習問題2-2のデータに対し、カイ二乗検定を行え。

In [179]:
 
import numpy as np
FoodTendency = np.array(["洋食","和食","和食","洋食","和食","洋食","洋食","和食","洋食","洋食","和食",
                         "洋食","和食","洋食","和食","和食","洋食","洋食","和食","和食"])
                         
TasteTendency = np.array(["甘党","辛党","甘党","甘党","辛党","辛党","辛党","辛党","甘党","甘党","甘党",
                         "甘党","辛党","辛党","甘党","辛党","辛党","甘党","辛党","辛党"])
 
<h4>演習問題4-5</h4>
次のそれぞれのデータについて無相関検定を行え。

演習問題4-5

次のそれぞれのデータについて無相関検定を行え。

In [ ]:
 
#5-1
import numpy as np
Kokugo = np.array([60,40,30,70,55])
Shakai = np.array([80,25,35,70,50])
In [ ]:
 
#5-2 --- 単純に(5-1)のデータを2回繰り返したもの
Kokugo = np.array([60,40,30,70,55,60,40,30,70,55])
Shakai = np.array([80,25,35,70,50,80,25,35,70,50])
 
<h4>演習問題4-6</H4>
<A HREF="SampleData/badmington.csv" target="_blank">badmington.csv</A>は区切り記号がコンマのCSVのファイルであり、
バドミントンのラケットの重量xと硬度yの表</A>(出典: 内田(2010)「すぐに使えるRによる統計解析とグラフの応用」東京図書)が収められている。このデータをデータフレームとして読み込み、硬度(y)と重量(x)の相関係数を算出し、無相関の検定を行え。
[参考] 区切り記号がコンマのcsvファイルを読み込み、その内容をデータフレームとして取り込むには、pandasモジュールの<B>read_csv</B>関数を用いる。

演習問題4-6

badmington.csvは区切り記号がコンマのCSVのファイルであり、 バドミントンのラケットの重量xと硬度yの表(出典: 内田(2010)「すぐに使えるRによる統計解析とグラフの応用」東京図書)が収められている。このデータをデータフレームとして読み込み、硬度(y)と重量(x)の相関係数を算出し、無相関の検定を行え。

[参考] 区切り記号がコンマのcsvファイルを読み込み、その内容をデータフレームとして取り込むには、pandasモジュールのread_csv関数を用いる。

In [ ]: