JavaScript での Mann-Whitney U検定 の実装

# Reactive stat 内部の Mann-Whitney U検定の実装 です

Mann-Whitney U検定は、2つの独立した標本が同じ母集団から抽出されたか、あるいは同じ分布を持つ母集団から抽出されたかを判断するためのノンパラメトリックな検定方法です。

この記事では、Mann-Whitney U検定の概要と、JavaScriptでの実装方法を説明します。

## はじめに

弊社では **[ブラウザだけで使える無料統計ソフト Reactive stat](https://www.emuyn.net/stats)** を提供しています。
信頼性の高い R で統計解析し、その結果を AI が解説します！
その背景には、統計に苦労している医療者の助けになりたい、という気持ちがあります。
最終的な統計解析は R で行うのですが、レスポンスとサーバー負荷の改善のため、一部の統計計算はブラウザ内で行っています。
そのうち、汎用性の高い部分を共有させていただきたいと思います。
できるだけ R の出力と整合性を持つように javascript をインプリメントしてあるので、参考になれば幸いです。

## 使用したライブラリ

このコードでは、jStatライブラリを使用しています。jStatは、JavaScriptで統計計算を行うための便利なライブラリです。
jStatの詳細については、[公式ドキュメント](https://github.com/jstat/jstat) を参照してください。
全関数の一覧は [jStat v1.9.3 Documentation](https://jstat.github.io/all.html) を参照してください。

## Mann-Whitney U検定の概要

Mann-Whitney U検定は以下のような場合に利用します:
- 2つの独立した標本があり、それらが同じ母集団から抽出されたかどうかを検証したい場合
- データが正規分布に従っていない、またはそのような仮定が不適切な場合
- 比較されるグループのサイズが小さい場合や、データに外れ値が含まれる場合

**検定の手順は以下の通りです**:

- 両方のグループの全データ点を合わせて順位付けします
- 各グループでの順位の合計を計算します
- $U_1$と$U_2$の値を計算し、最小値を選択します
- $U$値を使用してp値を計算します

**数式**:
- $U = \min(U_1, U_2)$

ここで、
- $U_1 = R_1 - \frac{n_1(n_1 + 1)}{2}$
- $U_2 = R_2 - \frac{n_2(n_2 + 1)}{2}$
- $n_1$ と $n_2$ はそれぞれ1つ目と2つ目のグループのサンプルサイズです
- $R_1$ と $R_2$ はそれぞれ1つ目と2つ目のグループの順位の合計です

    Mann-WhitneyのU検定では、正規近似を使用する場合には、連続性補正（continuity correction）を適用することが一般的です。
    
    $\text{連続性補正値} = \frac{1}{2}$
    
    これを使用して、標準化されたU値（Z）が次のように計算されます：
    
    $Z = \frac{U - \text{平均U} + \text{連続性補正値}}{\text{標準偏差U}}$
    
    連続性補正を適用することで、正確な分布に近い結果を得ることができます。

**注意点**:

Mann-Whitney U検定はノンパラメトリック検定であるため、データの分布が同じであるという仮定のもとで行われます。
大きなサンプルサイズの場合には正規近似を使用することができますが、小さなサンプルサイズの場合には正確な分布を使用する必要があります。

## コード

ライブラリ込みの html にしてありますので、
[jsfiddle](https://jsfiddle.net/) などにコピペして簡単に試せます。

```html:html
<!DOCTYPE html>
<html>
<head>
  <meta charset="utf-8">
  <title>Mann-Whitney U Test</title>
  <script src="https://cdn.jsdelivr.net/npm/jstat@latest/dist/jstat.min.js"></script>
</head>
<body>
  <h1>Mann-Whitney U Test Example</h1>
  <p>Check the console for the result.</p>
  
  <script>
    // Mann-Whitney のU検定
    // 正規近似を使用
    function mannWhitneyU(groupA, groupB) {
        // 全データのランクを計算
        const combined = groupA.concat(groupB);
        const ranks = jStat.rank(combined);
        const n1 = groupA.length;
        const n2 = groupB.length;

        // ランク和を計算
        const sum1 = ranks.slice(0, n1).reduce((a, b) => a + b, 0);
        const sum2 = ranks.slice(n1).reduce((a, b) => a + b, 0);

        const U1 = sum1 - (n1 * (n1 + 1)) / 2;
        const U2 = sum2 - (n2 * (n2 + 1)) / 2;
        const U = Math.min(U1, U2);

        const meanU = n1 * n2 / 2;
        const sdU = Math.sqrt((n1 * n2 * (n1 + n2 + 1)) / 12);

        // 連続性補正の適用
        const continuityCorrection = 0.5;
        const Z = (U - meanU + continuityCorrection) / sdU;

        return 2 * (1 - jStat.normal.cdf(Math.abs(Z), 0, 1)); // 両側検定のp値
    }

    // サンプルデータ
    const group1 = [31, 32, 43, 54, 65, 72, 80, 90, 92];
    const group2 = [22, 26, 27, 38, 49, 60, 65, 78, 84];

    const pValue = mannWhitneyU(group1, group2);
    console.log("Mann–Whitney U Test p-value: ", pValue);
  </script>
</body>
</html>
```

## 検証用 R コード

[rdrr](https://rdrr.io/snippets/) にコピペして簡単に試せます。

ここでは、Mann–Whitney U 検定を実行するために wilcox.test() 関数を使用します。
この関数は、Frank Wilcoxon によって開発された Wilcoxon 順位和検定を実装しており、その特別なケースとして Mann–Whitney U 検定が含まれています。
wilcox.test() 関数は Mann–Whitney U 検定を正確に実装しているため、結果は Mann–Whitney U 検定と同等になります。

```R:R
# Mann–Whitney U Test

# サンプルデータ
group1 <- c(31, 32, 43, 54, 65, 72, 80, 90, 92)
group2 <- c(22, 26, 27, 38, 49, 60, 65, 78, 84)

# Mann–Whitney U 検定を実行
result <- wilcox.test(group1, group2)

# 結果の表示
print(result)
```

## 注意点

- このコードは、jStatライブラリに依存しています。jStatを読み込まないと動作しません
- 浮動小数点数の計算では、わずかな誤差が生じる可能性があります
- このコードでは、常に正規近似を使用し、連続補正を適用しています。
  - 正確な分布を使用するコードは含まれていません
  - R でのデフォルトの動作は、サンプルサイズが大きい場合（両方のグループの合計が50以上）、正規近似が使用され、連続性補正が適用されます。そうでない場合、正確な方法が使用され、連続性補正は適用されません
  - そのため、得られる結果は、R の wilcox.test() とわずかに異なっています

## 最後に

コードの内容はよく吟味し、R との食い違いが極力ないように気を付けていますが、もし間違いに気づかれた場合には指摘していただけますとありがたいです。

ただ、Chat GPT に読ませてその出力を検証もせず「問題があるので修正します」として貼り付けられると、他の方を混乱させてしまうのと、弊社の製品に対する信頼を損ねることにもつながりかねませんので、どうかそのようなことのないようにお願いいたします。

検証のための R のコードもつけてありますので、ご活用いただけますと幸いです。