scipy.stats: 点双列相関係数 pointbiserialr

# scipy.stats: 点双列相関係数 pointbiserialr

特別な名前がついているが，ピアソンの積率相関係数において，片方の変数が 2 値データの場合ということである。

`pointbiserialr` は point biserial correlation coefficient r で，訳すと，点双列相関係数ということである。

2 値変数は連続変数なので（知らない人も多いかもしれないが），当たり前なのだが，その昔，計算環境が劣悪だった頃は，特別な場合に簡単な計算式で計算できるということで誕生したものかな。

普通に（数値に変換しなければならないかもしれないが） `scipy.stats.pearsonr` を使えばよい。


```python
import numpy as np

a = np.array([0, 0, 0, 1, 1, 1, 1])
b = np.arange(7)
```

## 1． 点双列相関係数

`pointbiserialr(x, y)`

戻り値はピアソンの積率相関係数と無相関検定（$H_0:\rho=0$）の両側検定の $p$ 値である。


```python
from scipy.stats import pointbiserialr

pointbiserialr(a, b)
```




    PointbiserialrResult(correlation=0.8660254037844387, pvalue=0.011724811003954626)



ちなみに，相関係数が $r$，サンプルサイズが $n$ のとき無相関検定の両側検定の $p$ 値は以下のように求められる。


```python
from scipy.stats import t

def pvalue(r, n):
    return t.sf(abs(r) * np.sqrt(n - 2) / np.sqrt(1 - r**2), n-2) * 2 if r**2 != 0 else 1

pvalue(0.8660254037844387, len(a))
```




    0.011724811003954614



## 2.　ピアソンの積率相関係数
    
`pearsonr(x, y)`

戻り値はピアソンの積率相関係数と無相関検定（$H_0: \rho = 0$）の両側検定の $p$ 値である。


```python
from scipy.stats import pearsonr

pearsonr(a, b)
```




    (0.8660254037844387, 0.011724811003954626)



## 3.　numpy の `corrcoef`

numpy には `corrcoef` がある。

`corrcoef(x, y=None, rowvar=True, bias=<no value>, ddof=<no value>, *, dtype=None)`


```python
np.corrcoef(a, b)[0,1]
```




    0.8660254037844385



ちなみに，`np.corrcoef()` は相関係数行列を求めることができる。無相関検定の結果は得られない。


```python
x = [[2,1,2,3,4], [3,2,1,2,3], [5,4,3,4,5], [2,1,3,4,5]]
r = np.corrcoef(x)
print(r)
```

    [[1.         0.41931393 0.41931393 0.97072534]
     [0.41931393 1.         1.         0.18898224]
     [0.41931393 1.         1.         0.18898224]
     [0.97072534 0.18898224 0.18898224 1.        ]]


以下のような関数で $p$ 値の行列を求めることができる（前出の `pvalue()` を使用する）。


```python
def r(x):
    k, n = np.array(x).shape
    p = np.ones((k, k))
    r = np.corrcoef(x)
    for i in range(k):
        for j in range(i+1, k):
            p[i, j] = p[j, i] = pvalue(r[i, j], n)
    return p

r(x)
```




    array([[1.00000000e+00, 4.82198952e-01, 4.82198952e-01, 5.98625570e-03],
           [4.82198952e-01, 1.00000000e+00, 1.40426542e-24, 7.60820376e-01],
           [4.82198952e-01, 1.40426542e-24, 1.00000000e+00, 7.60820376e-01],
           [5.98625570e-03, 7.60820376e-01, 7.60820376e-01, 1.00000000e+00]])



# 4.　φ 係数

ファイ係数は，2変数がともに 2 値データのときのピアソンの積率相関係数の絶対値に該当する。

さらに，2×2 分割表の場合には $\phi$ 係数は `scipy.stats.contingency` のデフォルトで得られる クラメールの $V$ に一致する。


```python
import numpy as np
np.random.seed(1234)
from scipy.stats import uniform
x = uniform.rvs(size=100) > 0.5
y = uniform.rvs(size=100) > 0.5
import pandas as pd
df = pd.DataFrame({'x': x, 'y': y})
tbl = pd.crosstab(df['x'], df['y'])
print(tbl)
```

    y      False  True
    x                 
    False     17    26
    True      27    30



```python
from scipy.stats.contingency import association

association(tbl)
```




    0.07812845233116676




```python
from scipy.stats import pearsonr

pearsonr(x, y)
```




    (-0.07812845233116661, 0.43973315252059714)



なお，$\phi$ 係数に適切な符号を付けてピアソンの積率相関係数と同じにするには，分割表の 4 つのマス目の数値を a, b, c, d としたとき，$a\cdot d - b \cdot c$ の符号をつければよい。

例の場合だと `17*30 - 26*27 = -192` なので，`0.07812845233116676` の符号は `-` である。

同じ分割表 `tbl` に基づく検定は，相関係数の検定の他に独立性の検定（いわゆる $\chi^2$ 検定）もある。同じ結果になるだろうか？

独立性の検定は `scipy.stats.chi2_contingency` で行える。
連続性の補正は `correction` で指定できる（デフォルトで　　`True`）。
戻り値は，$\chi^2$ 値，$p$ 値，自由度，期待値行列である。


```python
from scipy.stats import chi2_contingency
chi2_contingency(tbl)
```




    (0.3338817445304601,
     0.5633822668015476,
     1,
     array([[18.92, 24.08],
            [25.08, 31.92]]))




```python
chi2_contingency(tbl, correction=False)
```




    (0.6104055063663398,
     0.43463517068036106,
     1,
     array([[18.92, 24.08],
            [25.08, 31.92]]))



同じ $p$ 値になるわけではないということか。