rpy2 を使って R の wilcox.exact を使う（scipy の マンホイットニーU検定 の実装が 少サンプル数非対応 なので、代わりに）

# 追記（2020/2/25）

マンホイットニーのU検定よりBrunner-Munzel検定の方が良いそうです。

そのため、下の``wilcox.exact``を``brunnermunzel.test``（片群のサンプルサイズが10程度以下なら``brunnermunzel.permutation.test``）に置き換えて利用すると良いはずです。

**参考**

奥村晴彦 先生の解説：
[Brunner-Munzel検定](https://oku.edu.mie-u.ac.jp/~okumura/stat/brunner-munzel.html)

他の検定（t検定、Welchのt検定、U検定）との比較：
[マイナーだけど最強の統計的検定 Brunner-Munzel 検定 - ほくそ笑む](https://hoxo-m.hatenablog.com/entry/20150217/p1)


# 経緯

Jupyter notebook でマンホイットニーのU検定をしようと思いました。U検定を実装しているライブラリを探してみた所、``scipy``で実装されていました。具体的には``scipy.stats.mannwhitneyu``という名前で実装されています。

このライブラリを利用しようと思い、引数等を確認するために[公式Doc](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.mannwhitneyu.html)を見た所、``Notes``に次の様な注記がありました。

```
Use only when the number of observation in each sample is > 20
```

訳すと、「2群それぞれのサンプルサイズが 20超 の場合のみ利用可能」となります。

なぜこのような条件が付されているのか調べた所、どうやらこのライブラリは「標準正規分布による近似をして検定する方法」のみを実装している事が分かりました。

この事は[ソースコード](https://github.com/scipy/scipy/blob/v1.4.1/scipy/stats/stats.py#L6334-L6429)と[Wikipedia](https://ja.wikipedia.org/wiki/%E3%83%9E%E3%83%B3%E3%83%BB%E3%83%9B%E3%82%A4%E3%83%83%E3%83%88%E3%83%8B%E3%83%BC%E3%81%AEU%E6%A4%9C%E5%AE%9A#%E4%BD%BF%E7%94%A8%E6%96%B9%E6%B3%95)を見て分かりました。

そして、[scipy の Github](https://github.com/scipy/scipy/pull/4933) を調べると、どうやら片群のサンプルサイズが20以下の場合の実装はプルリク中のようでした。

さて、私は片群のサンプルサイズが20以下の場合の検定を行いたかったので、別のライブラリを探す必要が生まれました。

そこで、Python ライブラリをざっと調べて見たのですが、良さそうなものがパッとは見つかりませんでした。

どうしようか考えた所、「Rの実装を使えばいいのでは」と思いました。

そこで、RでU検定を実装したライブラリを探した所、``wilcox.exact``というライブラリが見つかりました。

このライブラリは、[データ科学便覧](https://data-science.gr.jp/implementation/ist_r_wilcoxon_rank_sum_test.html)で紹介されています。

[データ科学便覧](https://data-science.gr.jp/implementation/ist_r_wilcoxon_rank_sum_test.html)の解説によると、サンプルサイズによって自動で「標準正規分布による近似」をするか「p値を正確に計算」するかを決定してくれる事が分かりました。そしてこの機能により、scipyの様なサンプルサイズの制約はない事が分かりました。

この様に、良さそうなライブラリがRで見つかったので、Pythonカーネル上でRを使えるようにするライブラリ``rpy2``を使って、``wilcox.exact``を使う事にしました。

# 方法

## 1. rpy2を install します

```terminal:terminal
pip install rpy2
```

## 2. exactRankTests を install します

```bash:terminal
# R を起動します
$ R
```

```R:R（teminal上）
# exactRankTests を install
> install.packages("exactRankTests", repos="http://cran.ism.ac.jp/")

# 終了
> q()
```

## 3. Jupyter notebook 上で次のコードを実行して動作確認をします

```python:Input（Jupyter　notebook）
# %R が使えるようにする
%load_ext rpy2.ipython

# U検定が入っているライブラリをロード
%R library(exactRankTests)

# Python の配列を引数にとる U検定メソッド を作成
def u_test(x, y):
    
    # ndarray に変換
    # ndarrayなら 「%R -i」 した時に R で 「c()」 した時の形と同じになる
    x = np.array(x)
    y = np.array(y)

    # R側に x と y を読み込む
    %R -i x
    %R -i y
    
    # U検定
    %R output = wilcox.exact(x=x,y=y,paired=F)
    
    # Python 側に出力
    %R -o output
    
    # ListVector型で出力されるので、dictionary 型に変換
    result = dict(zip(output.names, map(list,list(output))))
    
    return result

# 試してみる
x = np.array([0.15165101, 0.47173771, 0.46735185, 0.67045667, 0.17489188])
y = x * 10
u_test(x, y)
```

```:Output（Jupyter　notebook）
{'statistic': [0.0],
 'pointprob': [0.003968253968253968],
 'p.value': [0.007936507936507936],
 'null.value': [0.0],
 'alternative': ['two.sided'],
 'method': ['Exact Wilcoxon rank sum test'],
 'data.name': ['x and y']}
```

``p.value``がp値です。

# 参考

**rpy2の使い方**

- [Jupyterで1つのノートブック内にRとPythonを共存させる - Qiita](https://qiita.com/tanakarian/items/609b550c4b18e5a4b925#%E3%81%BE%E3%81%A8%E3%82%81)