Reactive stat 内部の Bartlett検定の実装 です
Bartlett検定は、3つ以上のグループの分散が等しいかどうかを検定する方法の一つです。
分散分析のようないくつかの統計的検定では、群やサンプル間で分散が等しいことを前提とします。Bartlett検定は、その仮定を検証するために使用されます。
この記事では、Bartlett検定の概要と、JavaScriptでの実装方法を説明します。
はじめに
弊社では ブラウザだけで使える無料統計ソフト Reactive stat を提供しています。
信頼性の高い R で統計解析し、その結果を AI が解説します！
その背景には、統計に苦労している医療者の助けになりたい、という気持ちがあります。
最終的な統計解析は R で行うのですが、レスポンスとサーバー負荷の改善のため、一部の統計計算はブラウザ内で行っています。
そのうち、汎用性の高い部分を共有させていただきたいと思います。
できるだけ R の出力と整合性を持つように javascript をインプリメントしてあるので、参考になれば幸いです。
使用したライブラリ
このコードでは、jStatライブラリを使用しています。jStatは、JavaScriptで統計計算を行うための便利なライブラリです。
jStatの詳細については、公式ドキュメント を参照してください。
全関数の一覧は jStat v1.9.3 Documentation を参照してください。
Bartlett検定の概要
仮説
$H_0$: すべての母集団の分散は等しい。
$H_a$: 少なくとも2つの分散は異なる。
検定統計量$K^2の計算:
$K^2 = \frac{(n - k) \ln(s_p^2) - \sum_{i=1}^{k} (n_i - 1) \ln(s_i^2)}{1 + \frac{1}{3(k-1)}(\sum_{i=1}^{k}\frac{1}{n_i-1} - \frac{1}{n-k})}$
ここで：
- $k$ はグループの数,
- $n$ は合計サンプルサイズ,
- $n_i$ は $i$ 番目のグループのサンプルサイズ,
- $s_p^2$ はプールされた分散,
- $s_p^2 = \frac{\sum_{i=1}^{k} (n_i - 1)s_i^2}{n - k}$
- $s_i^2$ は $i$ 番目のグループの分散.
帰無仮説のもとで検定統計量 $K^2$ は、自由度が $k-1$ のカイ二乗分布に近似的に従います。したがって、p値はカイ二乗分布を用いて求めることができ、このp値が予め定められた有意水準（しばしば $\alpha = 0.05$ ）よりも小さい場合、等分散の帰無仮説は棄却されます。
仮定と注意
Bartlett検定は以下の仮定を置いています:
- 各データ群は独立しています。
- 各データは正規分布しています。
Bartlett検定は、正規性からの逸脱に対して敏感です。サンプルが非正規分布の場合、Levene検定を検討してください。
コード
ライブラリ込みの html にしてありますので、
jsfiddle などにコピペして簡単に試せます。
<!DOCTYPE html>
<html>
<head>
<meta charset="utf-8">
<title>Bartlett's Test</title>
<script src="https://cdn.jsdelivr.net/npm/jstat@latest/dist/jstat.min.js"></script>
</head>
<body>
<h1>Bartlett's Test Example</h1>
<p>Check the console for the result.</p>
<script>
// Bartlettの検定
function bartlettTest(...groups) {
const k = groups.length; // グループ数
// 各グループのサンプルサイズと分散を計算
const groupSizes = groups.map(group => group.length);
const variances = groups.map(group => {
const mean = group.reduce((acc, val) => acc + val) / group.length;
return group.reduce((acc, val) => acc + Math.pow(val - mean, 2), 0) / (group.length - 1);
});
// 全サンプルサイズを計算
const totalSize = groupSizes.reduce((acc, val) => acc + val); // 全サンプル数
// プールされた分散を計算
const pooledVariance = groups.reduce((acc, group, idx) => {
return acc + (groupSizes[idx] - 1) * variances[idx];
}, 0) / (totalSize - k);
// 修正項の計算
const correction = 1 + (1 / (3 * (k - 1))) * (groupSizes.reduce((acc, groupSize) => {
return acc + 1 / (groupSize - 1);
}, 0) - 1 / (totalSize - k));
// K^2統計量の計算
const K2 = ((totalSize - k) * Math.log(pooledVariance) - groupSizes.reduce((acc, groupSize, idx) => {
return acc + (groupSize - 1) * Math.log(variances[idx]);
}, 0)) / correction;
// p値の計算
const pValue = 1 - jStat.chisquare.cdf(K2, k - 1);
return pValue;
}
// サンプルデータ
const group1 = [85, 90, 92, 87, 99];
const group2 = [80, 82, 78, 83, 79];
const group3 = [84, 88, 91, 93, 86];
const pValue = bartlettTest(group1, group2, group3);
console.log("Bartlett's Test p-value: ", pValue);
</script>
</body>
</html>
検証用 R コード
rdrr にコピペして簡単に試せます。
# Bartlett's Test
# サンプルデータ
group1 <- c(85, 90, 92, 87, 99)
group2 <- c(80, 82, 78, 83, 79)
group3 <- c(84, 88, 91, 93, 86)
# Bartlett検定を実行
result <- bartlett.test(list(group1, group2, group3))
# 結果の表示
print(result)
注意点
- このコードは、jStatライブラリに依存しています。jStatを読み込まないと動作しません
- 浮動小数点数の計算では、わずかな誤差が生じる可能性があります
- Bartlett検定は、データが正規分布に従っていることを前提としています。データが正規分布に従っていない場合、結果の解釈には注意が必要です
- 群の数が2つの場合、F検定を使用することをお勧めします。Bartlett検定は3つ以上の群の比較に適しています
最後に
コードの内容はよく吟味し、R との食い違いが極力ないように気を付けていますが、もし間違いに気づかれた場合には指摘していただけますとありがたいです。
ただ、Chat GPT に読ませてその出力を検証もせず「問題があるので修正します」として貼り付けられると、他の方を混乱させてしまうのと、弊社の製品に対する信頼を損ねることにもつながりかねませんので、どうかそのようなことのないようにお願いいたします。
検証のための R のコードもつけてありますので、ご活用いただけますと幸いです。