【C++】Boost.Mathのfisher_fでフィッシャーF分布のPDF・CDF・逆累積分布を簡単に計算
C++のBoost.Mathライブラリにあるboost::math::fisher_fは、分子自由度df1と分母自由度df2を指定するだけでフィッシャーのF分布のPDFや累積分布、逆累積分布を計算できます。
pdf(dist,x)で確率密度、cdf(dist,x)で累積、quantile(dist,p)でpに対応するx値が得られ、分散分析などに便利です。
Boost.Math fisher_fクラスの基本構造
ヘッダーファイルと名前空間
Boostライブラリの数学関数を利用するためには、必要なヘッダーファイルをインクルードします。
フィッシャーのF分布を扱うためには、<boost/math/distributions/fisher_f.hpp>をインクルードします。
このヘッダーファイルには、boost::math名前空間内に定義されているfisher_fクラスが含まれています。
#include <boost/math/distributions/fisher_f.hpp>このインクルードにより、boost::math::fisher_fクラスや、そのメンバ関数にアクセスできるようになります。
また、Boostの数学関数はboost::math名前空間に属しているため、コード内ではboost::mathを明示的に指定するか、usingディレクティブを用いて名前空間を省略することも可能です。
using namespace boost::math;ただし、名前空間の汚染を避けるために、必要な部分だけを明示的に指定するのが一般的です。
コンストラクタと内部パラメータ
boost::math::fisher_fクラスは、フィッシャーのF分布を表現するためのクラスです。
このクラスのインスタンスを作成するには、2つの自由度パラメータを指定します。
boost::math::fisher_f dist(df1, df2);ここで、df1は分子の自由度(正の実数)、df2は分母の自由度(正の実数)を表します。
これらの値は、分布の形状を決定し、統計的な検定や確率計算に用いられます。
自由度パラメータの詳細
df1(分子自由度)
分子側のカイ二乗分布の自由度を示し、通常は正の実数値を取ります。
値が大きくなるほど、分布は正規分布に近づきます。
df2(分母自由度)
分母側のカイ二乗分布の自由度を示し、こちらも正の実数値です。
df2が大きくなると、分布の尾部が短くなり、より集中した形状になります。
インスタンスの作成例
// 自由度を設定
double df1 = 5.0; // 分子自由度
double df2 = 10.0; // 分母自由度
// フィッシャーのF分布のインスタンスを作成
boost::math::fisher_f dist(df1, df2);パラメータの検証と例外処理
自由度の値は正の実数でなければなりません。
もし0または負の値を指定した場合、boost::math::domain_error例外がスローされるため、例外処理を適切に行う必要があります。
try {
boost::math::fisher_f dist(-3.0, 10.0); // 不正な値
} catch (const boost::math::domain_error& e) {
std::cerr << "エラー: " << e.what() << std::endl;
}このように、fisher_fクラスはシンプルなインターフェースで、自由度の設定とともに分布のインスタンスを作成し、その後の確率計算や統計的操作に利用します。
フィッシャーF分布オブジェクトの生成
df1(分子自由度)とdf2(分母自由度)の指定
boost::math::fisher_fクラスのインスタンスを作成する際には、まず分子自由度df1と分母自由度df2を明示的に指定します。
これらのパラメータは、分布の形状を決定し、確率計算の基礎となるため、正確な値を設定することが重要です。
// 例：分子自由度と分母自由度を設定
double df1 = 5.0; // 分子自由度
double df2 = 10.0; // 分母自由度
// 分布オブジェクトの生成
boost::math::fisher_f dist(df1, df2);この例では、df1が5、df2が10に設定されており、これらの値に基づいて分布の確率密度や累積確率を計算できます。
正の実数チェック
自由度の値は、常に正の実数でなければなりません。
もし0または負の値を指定した場合、boost::math::domain_error例外がスローされます。
これを防ぐために、インスタンス生成前に値の妥当性を検証することが推奨されます。
// 自由度の値を検証
if (df1 <= 0 || df2 <= 0) {
throw std::invalid_argument("自由度は正の実数でなければなりません。");
}また、例外処理を用いて安全にインスタンスを作成することも可能です。
try {
boost::math::fisher_f dist(df1, df2);
} catch (const boost::math::domain_error& e) {
std::cerr << "エラー: " << e.what() << std::endl;
}動的パラメータ調整
一度作成したfisher_fオブジェクトのパラメータdf1やdf2は、コンストラクタ呼び出し後に変更できません。
もし異なる自由度の値で再度分布を定義したい場合は、新たにインスタンスを生成し直す必要があります。
ただし、分布のパラメータを動的に調整したい場合は、次のように新しいインスタンスを作成します。
// パラメータを変更したい場合は、新たにインスタンスを作成
boost::math::fisher_f new_dist(8.0, 15.0);この方法により、異なる自由度の設定を持つ複数の分布を同時に扱うことが可能です。
コピー・ムーブ操作のサポート
boost::math::fisher_fクラスは、コピーコンストラクタとムーブコンストラクタをサポートしています。
これにより、インスタンスの複製や効率的な移動が可能です。
// コピーによるインスタンスの複製
boost::math::fisher_f dist_copy = dist;
// ムーブによるインスタンスの移動(C++11以降)
boost::math::fisher_f dist_moved = std::move(dist);これらの操作は、分布の状態を保持したまま、複数の場所で同じパラメータの分布を使いたい場合に便利です。
ただし、boost::math::fisher_fは基本的に値型であり、内部に動的なリソースを持たないため、コピーやムーブの操作は非常に高速です。
確率密度関数(PDF)の計算
pdf(dist, x)の呼び出し方
boost::math::fisher_fクラスのインスタンスに対して、特定の値xにおける確率密度関数(PDF)の値を計算するには、pdf関数を用います。
これは、次のように呼び出します。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
// 自由度を設定
double df1 = 5.0;
double df2 = 10.0;
// 分布のインスタンスを作成
boost::math::fisher_f dist(df1, df2);
// 計算したい点
double x = 2.5;
// PDFの計算
double pdf_value = boost::math::pdf(dist, x);
// 結果の出力
std::cout << "x = " << x << "におけるPDFの値: " << pdf_value << std::endl;
return 0;
}x = 2.5におけるPDFの値: 0.0935948この例では、distはboost::math::fisher_fのインスタンス、xは計算したい点です。
boost::math::pdf関数は、第一引数に分布のインスタンス、第二引数に点xを取り、対応する確率密度の値を返します。
引数と戻り値の型
- 引数:
dist:const boost::math::fisher_f&型の分布オブジェクトx:double型の点の値
- 戻り値:
double型の確率密度値
この関数は、xの値に対して分布の確率密度を計算し、その結果を返します。
計算結果は、xの周辺の確率密度を示し、値が高いほど、その点付近に確率質量が集中していることを意味します。
計算例と結果の解釈
例えば、df1=5、df2=10の分布に対して、x=2.5のときのPDF値を計算します。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
boost::math::fisher_f dist(5, 10);
double x = 2.5;
double pdf_value = boost::math::pdf(dist, x);
std::cout << "PDF at x = " << x << ": " << pdf_value << std::endl;
return 0;
}このプログラムを実行すると、次のような出力が得られます。
PDF at x = 2.5: 0.0935948この値は、x=2.5の点における確率密度が約0.093であることを示しています。
確率密度は確率そのものではなく、単位区間あたりの確率の密度を表すため、値が高いほど、その点付近に確率質量が集中していると解釈します。
例えば、xの値を変化させて複数の点でPDFを計算し、その分布の形状を把握することも可能です。
double points[] = {1.0, 2.0, 3.0, 4.0, 5.0};
for (double point : points) {
double val = boost::math::pdf(dist, point);
std::cout << "x = " << point << "のPDF値: " << val << std::endl;
}x = 1のPDF値: 0.49548
x = 2のPDF値: 0.162006
x = 3のPDF値: 0.0558267
x = 4のPDF値: 0.0218973
x = 5のPDF値: 0.00963064これにより、分布のピークや尾部の挙動を視覚的に理解できます。
複数サンプルへの適用
複数の点に対して一括で確率密度を計算したい場合、forループやSTLのアルゴリズムを活用できます。
forループによる一括計算
#include <vector>
#include <iostream>
#include <boost/math/distributions/fisher_f.hpp>
int main() {
boost::math::fisher_f dist(5, 10);
std::vector<double> x_values = {1.0, 1.5, 2.0, 2.5, 3.0};
for (size_t i = 0; i < x_values.size(); ++i) {
double x = x_values[i];
double pdf_val = boost::math::pdf(dist, x);
std::cout << "x = " << x << "のPDF値: " << pdf_val << std::endl;
}
return 0;
}
x = 1のPDF値: 0.49548
x = 1.5のPDF値: 0.286459
x = 2のPDF値: 0.162006
x = 2.5のPDF値: 0.0935948
x = 3のPDF値: 0.0558267この例では、x_valuesの各値に対してpdfを計算し、結果を逐次出力します。
範囲ベースfor／STLアルゴリズムの活用
C++11以降の範囲ベースfor文を用いると、コードがより簡潔になります。
#include <vector>
#include <iostream>
#include <boost/math/distributions/fisher_f.hpp>
int main() {
boost::math::fisher_f dist(5, 10);
std::vector<double> x_values = {1.0, 1.5, 2.0, 2.5, 3.0};
for (double x : x_values) {
std::cout << "x = " << x << "のPDF値: " << boost::math::pdf(dist, x) << std::endl;
}
return 0;
}x = 1のPDF値: 0.49548
x = 1.5のPDF値: 0.286459
x = 2のPDF値: 0.162006
x = 2.5のPDF値: 0.0935948
x = 3のPDF値: 0.0558267また、std::transformを用いて結果を別のコンテナに格納することも可能です。
#include <vector>
#include <algorithm>
#include <iostream>
#include <boost/math/distributions/fisher_f.hpp>
int main() {
boost::math::fisher_f dist(5, 10);
std::vector<double> x_values = {1.0, 1.5, 2.0, 2.5, 3.0};
std::vector<double> pdf_results(x_values.size());
std::transform(x_values.begin(), x_values.end(), pdf_results.begin(),
[&dist](double x) { return boost::math::pdf(dist, x); });
for (size_t i = 0; i < x_values.size(); ++i) {
std::cout << "x = " << x_values[i] << "のPDF値: " << pdf_results[i] << std::endl;
}
return 0;
}x = 1のPDF値: 0.49548
x = 1.5のPDF値: 0.286459
x = 2のPDF値: 0.162006
x = 2.5のPDF値: 0.0935948
x = 3のPDF値: 0.0558267このように、複数の点に対して効率的に確率密度を計算し、分布の形状を詳細に把握することが可能です。
累積分布関数(CDF)の計算
cdf(dist, x)の使い方
boost::math::fisher_fクラスのインスタンスに対して、ある点xまでの確率を求めるには、cdf関数を使用します。
これは、次のように呼び出します。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
// 自由度を設定
double df1 = 5.0;
double df2 = 10.0;
// 分布のインスタンスを作成
boost::math::fisher_f dist(df1, df2);
// 計算したい点
double x = 2.5;
// 累積分布関数の計算
double cdf_value = boost::math::cdf(dist, x);
// 結果の出力
std::cout << "x = " << x << "までの累積確率: " << cdf_value << std::endl;
return 0;
}x = 2.5までの累積確率: 0.897998この例では、distはboost::math::fisher_fのインスタンス、xは計算したい点です。
boost::math::cdf関数は、x以下の値をとる確率を返します。
累積確率の意味
cdf(dist, x)が返す値は、x以下の範囲における確率の総和です。
つまり、確率密度関数(PDF)をxから無限大まで積分した結果に相当します。
これにより、ある閾値xを超える確率や、特定の範囲内にデータが収まる確率を計算できます。
例えば、x=2.5のときのcdf値が0.85であれば、全体の85%の確率がx=2.5以下の値に収まることを意味します。
上側確率との関係
cdfは、ある点xまでの確率を表しますが、逆にxより大きい値の確率も重要です。
これを計算するには、boost::math::complement関数を用います。
#include <boost/math/distributions/complement.hpp>
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
double df1 = 5.0, df2 = 10.0;
boost::math::fisher_f dist(df1, df2);
double x = 2.5;
// 上側確率 P(X > x) を得る
double upper_prob = boost::math::cdf(boost::math::complement(dist, x));
std::cout << "x = " << x << " より大きい確率: " << upper_prob << std::endl;
return 0;
}
この関数は、xより大きい値をとる確率を返します。
以下のように修正すると正確です：
cdf(boost::math::complement(dist, x)) は計算上1.0 - boost::math::cdf(dist, x)と等しくなります。
信頼水準(α)の例
信頼区間や統計的検定において、上側確率は重要な役割を果たします。
例えば、信頼水準α=0.05に対応する上側確率は0.05です。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
double df1 = 5.0;
double df2 = 10.0;
boost::math::fisher_f dist(df1, df2);
double alpha = 0.05;
// 信頼水準に対応する上側確率の閾値
double critical_value = boost::math::quantile(complement(dist, alpha));
std::cout << "信頼水準 " << (1 - alpha) * 100 << "% の閾値: " << critical_value << std::endl;
return 0;
}信頼水準 95% の閾値: 3.32583この例では、quantile(complement(dist, alpha))を用いて、α=0.05に対応する上側確率の臨界値を求めています。
これにより、観測値がこの閾値を超える確率が5%未満となるような閾値を得ることができます。
このように、cdfとcdf_complementは、統計的な判断や信頼区間の設定において非常に重要な役割を果たします。
逆累積分布関数(Quantile)の取得
quantile(dist, p)の用法
boost::math::fisher_fクラスのインスタンスに対して、確率pに対応する分位点(逆累積分布関数の値)を求めるには、quantile関数を使用します。
これは次のように呼び出します。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
// 自由度を設定
double df1 = 5.0;
double df2 = 10.0;
// 分布のインスタンスを作成
boost::math::fisher_f dist(df1, df2);
// 確率p
double p = 0.95;
// pに対応する分位点を計算
double q_value = boost::math::quantile(dist, p);
// 結果の出力
std::cout << "確率 p = " << p << " に対応する分位点: " << q_value << std::endl;
return 0;
}確率 p = 0.95 に対応する分位点: 3.32583この例では、p=0.95に対して、分布の上側5%点(95%信頼区間の上限)を求めています。
p値の範囲チェック
quantile関数に渡すpは、0から1の範囲内の値でなければなりません。
範囲外の値を指定した場合、boost::math::domain_error例外がスローされるため、事前に値の妥当性を確認することが重要です。
double p = 1.2; // 不正な値
if (p < 0.0 || p > 1.0) {
throw std::invalid_argument("pは0から1の範囲内で指定してください。");
}または、例外処理を用いて安全に呼び出すことも推奨されます。
try {
double q_value = boost::math::quantile(dist, p);
} catch (const boost::math::domain_error& e) {
std::cerr << "エラー: " << e.what() << std::endl;
}下側・上側の区別
quantileは、pに対応する下側確率の分位点を返します。
例えば、p=0.05の場合は、分布の下側5%点を示します。
一方、上側確率に対応する分位点を求めたい場合は、pを1 - pに変換してquantileに渡します。
double p_upper = 0.05; // 上側5%
double p_lower = 1.0 - p_upper; // 下側95%
double upper_quantile = boost::math::quantile(dist, p_upper);
double lower_quantile = boost::math::quantile(dist, p_lower);このようにして、信頼区間の境界値や臨界値を計算できます。
信頼区間の境界値計算
信頼区間の境界値は、quantile関数を用いて計算します。
例えば、95%信頼区間の下限と上限を求める場合、次のようにします。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
double df1 = 5.0;
double df2 = 10.0;
boost::math::fisher_f dist(df1, df2);
double alpha = 0.05; // 5%の有意水準
double lower_bound = boost::math::quantile(dist, alpha / 2);
double upper_bound = boost::math::quantile(dist, 1 - alpha / 2);
std::cout << "95%信頼区間の下限: " << lower_bound << std::endl;
std::cout << "95%信頼区間の上限: " << upper_bound << std::endl;
return 0;
}95%信頼区間の下限: 0.151077
95%信頼区間の上限: 4.23609この例では、alpha/2と1 - alpha/2を用いて、両側の臨界値を計算しています。
これにより、観測値がこの範囲内に収まる確率が95%となる信頼区間の境界値を得ることができます。
また、信頼区間の計算は、統計的仮説検定やパラメータ推定の際に頻繁に用いられ、quantile関数を使うことで簡潔に求めることが可能です。
例外処理とエラー対策
domain_errorの捕捉
boost::mathライブラリの関数を使用する際には、入力値やパラメータが不正な場合にstd::domain_error例外がスローされることがあります。
これを適切に捕捉し、プログラムのクラッシュを防ぐためには、try -catch構造を用いることが重要です。
#include <boost/math/distributions/fisher_f.hpp> // F 分布
#include <exception> // std::exception
#include <iostream>
int main() {
double df1 = -5.0; // 不正な自由度
double df2 = 10.0; // 正常な自由度
try {
// F 分布オブジェクトを生成（ここで自由度が負だと例外が投げられる）
boost::math::fisher_f dist(df1, df2);
double x = 2.0;
// 確率密度関数を評価
double result = boost::math::pdf(dist, x);
std::cout << "PDF: " << result << std::endl;
} catch (const std::domain_error& e) {
// 例外をキャッチしてエラーメッセージを表示
std::cerr << "例外が発生しました: " << e.what() << std::endl;
}
return 0;
}
例外が発生しました: Error in function fisher_f_distribution<double>::fisher_f_distribution: Degrees of freedom argument is -5, but must be > 0 !この例では、df1に不正な値(負の値)を設定しているため、std::domain_error例外がスローされます。
tryブロック内で例外が発生すると、catchブロックに制御が移り、エラーメッセージを出力します。
不正パラメータ時の対応策
パラメータの入力値が不正な場合に備え、事前に値の妥当性を検証することが推奨されます。
具体的には、自由度の値が正の実数であるかどうかを確認し、不正な場合は例外を投げたり、エラー処理を行ったりします。
#include <iostream>
#include <stdexcept>
#include <boost/math/distributions/fisher_f.hpp>
void validate_parameters(double df1, double df2) {
if (df1 <= 0 || df2 <= 0) {
throw std::invalid_argument("自由度は正の実数でなければなりません。");
}
}
int main() {
double df1 = -3.0; // 不正な値
double df2 = 10.0;
try {
validate_parameters(df1, df2);
boost::math::fisher_f dist(df1, df2);
double x = 2.0;
double pdf_value = boost::math::pdf(dist, x);
std::cout << "PDF: " << pdf_value << std::endl;
} catch (const std::invalid_argument& e) {
std::cerr << "パラメータエラー: " << e.what() << std::endl;
} catch (const std::domain_error& e) {
std::cerr << "ドメインエラー: " << e.what() << std::endl;
}
return 0;
}パラメータエラー: 自由度は正の実数でなければなりません。この例では、validate_parameters関数で自由度の値を検証し、不正な場合はstd::invalid_argument例外を投げます。
これにより、パラメータの誤りを早期に検出し、適切なエラー処理を行うことが可能です。
また、ユーザからの入力や外部データを扱う場合は、入力値の範囲や型を厳密に検証し、不正な値に対して適切なエラー通知や修正処理を行うことが、堅牢なプログラム設計のポイントとなります。
ランダムサンプリングとの組み合わせ
Boost.Randomとの統合
boost::randomライブラリを用いることで、boost::math::fisher_f分布からの乱数サンプルを生成できます。
boost::randomのエンジンと分布オブジェクトを組み合わせることで、擬似乱数を用いたシミュレーションやモンテカルロ法などの確率的計算を効率的に行えます。
まず、必要なヘッダーファイルをインクルードします。
#include <boost/random.hpp>
#include <iostream>
#include <vector>次に、乱数エンジンと分布オブジェクトを作成し、サンプルを生成します。
#include <boost/random.hpp>
#include <iostream>
#include <vector>
int main() {
// 乱数エンジンの作成（メルセンヌ・ツイスタ）
// 乱数エンジンのシード値を指定することで、同じ乱数列を再現可能
boost::random::mt19937 rng(12345);
// fisher_f 分布（分子自由度 5, 分母自由度 10）
boost::random::fisher_f_distribution<double> dist_gen(5.0, 10.0);
// 100 個のサンプルを生成
std::vector<double> samples;
samples.reserve(100);
for (int i = 0; i < 100; ++i) {
double sample = dist_gen(rng);
samples.push_back(sample);
}
// 生成したサンプルの先頭を 5 件出力
for (size_t i = 0; i < 5; ++i) {
std::cout << "サンプル[" << i << "] = " << samples[i] << std::endl;
}
return 0;
}
サンプル[0] = 0.243988
サンプル[1] = 0.897873
サンプル[2] = 0.749119
サンプル[3] = 1.20324
サンプル[4] = 3.01069この例では、boost::random::fisher_f_distributionを用いて、rngからフィッシャーF分布に従う乱数を生成しています。
samplesに格納された値は、実際の分布に従った乱数サンプルです。
サンプル取得例
std::vectorに格納したサンプルを用いて、統計的解析やヒストグラム作成、確率推定などに利用できます。
#include <boost/random/fisher_f_distribution.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <iostream>
#include <numeric> // std::accumulate
#include <vector>
int main() {
// 乱数エンジンと Fisher F 分布 (自由度 5, 10)
boost::random::mt19937 rng(12345); // 乱数エンジンの初期化
boost::random::fisher_f_distribution<double> dist(5.0, 10.0);
const int sample_size = 1000;
std::vector<double> samples;
samples.reserve(sample_size);
// サンプルの生成
for (int i = 0; i < sample_size; ++i) {
samples.push_back(dist(rng));
}
// 平均値を計算
double sum = std::accumulate(samples.begin(), samples.end(), 0.0);
double mean = sum / sample_size;
std::cout << "サンプルの平均値: " << mean << std::endl;
return 0;
}
サンプルの平均値: 1.28782この例では、1000個のサンプルを生成し、その平均値を計算しています。
これにより、分布の特性をシミュレーションによって把握できます。
また、サンプルのヒストグラムを作成したり、分布の確率密度関数と比較したりすることで、理論値と実測値の整合性を検証することも可能です。
このように、boost::randomとboost::mathを組み合わせることで、確率分布からの乱数生成と統計解析を効率的に行えるため、シミュレーションやモデリングの幅が広がります。
パフォーマンス最適化
計算コストの把握
boost::mathの関数は高精度な計算を提供しますが、その分計算コストも考慮する必要があります。
特に、大量の確率計算やシミュレーションを行う場合、パフォーマンスの最適化は重要です。
関数呼び出しオーバーヘッド
pdfやcdf、quantileといった関数は、内部で複雑な数値計算を行います。
これらの関数呼び出しには一定のオーバーヘッドが伴い、多数回呼び出すとパフォーマンスに影響します。
例えば、ループ内で何度も同じ分布に対して関数を呼び出す場合、呼び出し回数を減らす工夫が必要です。
事前に計算結果をキャッシュしたり、必要な値だけを計算したりすることで、効率的な実装が可能です。
// 例：ループ内での関数呼び出しの最適化
double x_values[] = {1.0, 2.0, 3.0, 4.0, 5.0};
double cached_pdf = boost::math::pdf(dist, x_values[0]); // 最初に一度だけ計算
for (size_t i = 1; i < sizeof(x_values)/sizeof(x_values[0]); ++i) {
double pdf_value = boost::math::pdf(dist, x_values[i]);
// 追加処理
}精度と速度のトレードオフ
高精度な計算は時間がかかる場合があります。
boost::mathはデフォルトで高い精度を提供しますが、場合によっては計算速度を優先したいこともあります。
boost::mathの関数には、精度を調整できるオプションは基本的にありませんが、必要に応じて近似手法や低精度の計算を選択することも検討できます。
例えば、特定の範囲内での近似式を用いることで、計算時間を短縮できる場合があります。
並列処理への応用
大量の確率計算やサンプル生成を行う場合、並列処理を活用してパフォーマンスを向上させることが可能です。
OpenMP／std::threadの活用
OpenMPを用いると、簡潔に並列化が可能です。
例えば、複数のx値に対してpdfを計算する場合、次のように書きます。
#include <omp.h>
#include <vector>
#include <iostream>
#include <boost/math/distributions/fisher_f.hpp>
int main() {
boost::math::fisher_f dist(5, 10);
std::vector<double> x_values = {1.0, 1.5, 2.0, 2.5, 3.0};
std::vector<double> results(x_values.size());
#pragma omp parallel for
for (size_t i = 0; i < x_values.size(); ++i) {
results[i] = boost::math::pdf(dist, x_values[i]);
}
for (size_t i = 0; i < results.size(); ++i) {
std::cout << "x = " << x_values[i] << "のPDF: " << results[i] << std::endl;
}
return 0;
}x = 1のPDF: 0.49548
x = 1.5のPDF: 0.286459
x = 2のPDF: 0.162006
x = 2.5のPDF: 0.0935948
x = 3のPDF: 0.0558267このコードでは、#pragma omp parallel forディレクティブにより、ループの各反復が複数のスレッドで並列に実行され、計算時間を短縮します。
std::threadを用いる場合は、スレッドを明示的に管理し、各スレッドに計算を割り当てる必要があります。
#include <thread>
#include <vector>
#include <iostream>
#include <boost/math/distributions/fisher_f.hpp>
void compute_pdf(const boost::math::fisher_f& dist, double x, double& result) {
result = boost::math::pdf(dist, x);
}
int main() {
boost::math::fisher_f dist(5, 10);
std::vector<double> x_values = {1.0, 1.5, 2.0, 2.5, 3.0};
std::vector<double> results(x_values.size());
std::vector<std::thread> threads;
for (size_t i = 0; i < x_values.size(); ++i) {
threads.emplace_back(compute_pdf, std::ref(dist), x_values[i], std::ref(results[i]));
}
for (auto& t : threads) {
t.join();
}
for (size_t i = 0; i < results.size(); ++i) {
std::cout << "x = " << x_values[i] << "のPDF: " << results[i] << std::endl;
}
return 0;
}x = 1のPDF: 0.49548
x = 1.5のPDF: 0.286459
x = 2のPDF: 0.162006
x = 2.5のPDF: 0.0935948
x = 3のPDF: 0.0558267この例では、各スレッドがcompute_pdf関数を呼び出し、並列に計算を行います。
スレッドの管理や同期に注意が必要ですが、大規模な計算では効果的です。
パフォーマンス最適化は、計算コストを理解し、適切な並列化や工夫を行うことで、より効率的なプログラムを作成できます。
実践的な利用シナリオ
F検定への応用
F検定は、二つの母集団の分散が等しいかどうかを検定するために広く用いられます。
boost::math::fisher_f分布は、この検定の統計量の分布として直接利用されます。
データ前処理と統計量算出
まず、二つのデータセットの分散を計算します。
次に、それぞれのサンプルサイズと分散を用いてF値を算出します。
#include <boost/math/distributions/complement.hpp>
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
#include <numeric>
#include <vector>
// 不偏分散を計算する関数
double variance(const std::vector<double>& data) {
double mean = std::accumulate(data.begin(), data.end(), 0.0) / data.size();
double sum_sq = 0.0;
for (double val : data) {
sum_sq += (val - mean) * (val - mean);
}
return sum_sq / (data.size() - 1);
}
int main() {
// サンプルデータ
std::vector<double> data1 = {2.3, 2.5, 2.7, 2.6, 2.4};
std::vector<double> data2 = {3.1, 3.3, 3.2, 3.4, 3.5};
// 各データの不偏分散を計算
double var1 = variance(data1);
double var2 = variance(data2);
// 標本サイズ
int n1 = static_cast<int>(data1.size());
int n2 = static_cast<int>(data2.size());
// F値 = 分散の比
double F = var1 / var2;
// 自由度
double df1 = n1 - 1;
double df2 = n2 - 1;
// F分布オブジェクト（ブースト）
boost::math::fisher_f dist(df1, df2);
// 片側下側確率
double lower_tail = boost::math::cdf(dist, F);
// 片側上側確率
double upper_tail = boost::math::cdf(dist, F);
// 両側検定の p 値
double p_value = 2.0 * std::min(lower_tail, upper_tail);
// 結果出力
std::cout << "F値: " << F << std::endl;
std::cout << "p値: " << p_value << std::endl;
return 0;
}
この例では、二つのデータセットの分散を計算し、F値を求め、その値に基づいてp値を算出しています。
p値が有意水準(例：0.05)未満であれば、分散の等質性は棄却されます。
モデル評価での活用
回帰分析やその他の統計モデルの評価においても、F分布は重要な役割を果たします。
特に、分散分析(ANOVA)では、モデルの説明変数が従属変数に与える影響の有意性を検定します。
分散分析結果の解釈
分散分析の結果から得られるF値は、モデルの説明力を示す指標です。
boost::math::fisher_fを用いて、F値の有意性を検定します。
#include <boost/math/distributions/fisher_f.hpp>
#include <iostream>
int main() {
double F_value = 4.67;
int df_between = 3;
int df_within = 20;
// F 分布を作成
boost::math::fisher_f dist(df_between, df_within);
// 補完分布関数（上側確率）を直接計算
double p_value = boost::math::cdf(boost::math::complement(
dist, F_value)); // ← complement(dist, F_value) 形式で渡す
// :contentReference[oaicite:0]{index=0}
std::cout << "F値: " << F_value << "\n";
std::cout << "p値: " << p_value << "\n";
if (p_value < 0.05) {
std::cout << "有意差あり (帰無仮説棄却)\n";
} else {
std::cout << "有意差なし (帰無仮説採択)\n";
}
return 0;
}
F値: 4.67
p値: 0.0124682
有意差あり (帰無仮説棄却)この例では、F値と自由度からp値を計算し、有意水準と比較しています。
p値が0.05未満であれば、グループ間に有意な差があると判断します。
このように、boost::math::fisher_fは、実務の統計解析や研究において、検定やモデル評価のための重要なツールとして広く活用されます。
トラブルシューティング
コンパイルエラーの原因例
ライブラリ未リンク／ヘッダー忘れ
boost::mathやboost::randomの関数を使用しているにもかかわらず、コンパイル時にエラーが発生する場合、最も一般的な原因は必要なライブラリやヘッダーファイルが正しくインクルードまたはリンクされていないことです。
主な原因と対策:
- ヘッダーファイルの未インクルード
必要なヘッダーファイルを忘れていると、関数やクラスが未定義となりエラーになります。
例として、<boost/math/distributions/fisher_f.hpp>を忘れると、boost::math::fisher_fが未定義となります。
- ライブラリの未リンク
Boostライブラリはヘッダオンリーの部分もありますが、一部のコンパイル環境やビルドシステムでは、ライブラリのリンク設定が必要です。
特に、Boost.RandomやBoost.Threadなどを使用する場合は、リンクオプションに-lboost_threadや-lboost_systemを追加します。
例：g++のコンパイルコマンド
g++ -std=c++17 main.cpp -lboost_system -lboost_thread -o my_program- Boostのインストールパスが標準でない場合は、
-Iや-Lオプションでパスを指定します - Boostのバージョンに応じて、リンクすべきライブラリ名が異なる場合があります
ランタイムエラーの対処
ドメインエラー／数値オーバーフロー
実行時にエラーが発生する場合、最も多い原因は入力値の不正や計算結果の数値的な問題です。
主な原因と対策:
- ドメインエラー
boost::mathの関数に不正なパラメータを渡すと、std::domain_error例外がスローされます。
例えば、自由度や確率pが範囲外の場合です。
- パラメータの範囲を事前に検証します
- 例外処理を用いてエラーを捕捉し、適切なメッセージを出力または処理を中断します
try {
double p = 1.2; // 不正な値
double q = boost::math::quantile(dist, p);
} catch (const std::domain_error& e) {
std::cerr << "ドメインエラー: " << e.what() << std::endl;
}- 数値オーバーフローやアンダーフロー
大きな値や極端に小さな値を扱うと、計算結果が浮動小数点の範囲を超えてしまい、infやNaNになることがあります。
- 入力値の範囲を制限し、極端な値を避けます
- 計算途中の値を逐次検査し、異常値が出たら処理を中断または修正します
- 必要に応じて、数値的に安定なアルゴリズムや近似式を採用します
double safe_pdf(const boost::math::fisher_f& dist, double x) {
if (x < 0) {
// F分布の定義域外
throw std::domain_error("xは0以上の値でなければなりません。");
}
double result = boost::math::pdf(dist, x);
if (std::isnan(result) || std::isinf(result)) {
// 計算結果が不正
throw std::runtime_error("計算結果が不正です。");
}
return result;
}エラーの原因を理解し、事前の入力検証と例外処理を徹底することで、プログラムの堅牢性を高めることができます。
特に、数値計算においては、入力値の範囲や計算結果の妥当性を常に確認することが重要です。
まとめ
この記事では、Boostライブラリを用いたフィッシャーF分布の基本操作や計算方法、例外処理のポイント、並列処理の活用例などを解説しました。
確率密度関数や累積分布関数、逆関数の計算方法や、実務での応用例も紹介しています。
これにより、統計解析やシミュレーションの効率化と信頼性向上に役立てられます。
![[C++] boost::jsonでのJSON書き込み方法](https://af-e.net/wp-content/uploads/2024/09/thumbnail-36107.png)