【C++】Boost.Mathでガンマ分布を扱う：PDF・CDF計算と乱数生成

#include <boost/math/distributions/gamma.hpp>

using namespace boost::math;

using boost::math::pdf;
using boost::math::cdf;
using boost::math::gamma_distribution;

boost::math::gamma_distribution<double> dist(α, β);

#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 形状母数を2.0、尺度母数を3.0に設定
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    // 作成した分布の情報を出力
    std::cout << "ガンマ分布のパラメータ:" << std::endl;
    std::cout << "形状母数 (α): 2.0" << std::endl;
    std::cout << "尺度母数 (β): 3.0" << std::endl;
    return 0;
}

ガンマ分布のパラメータ:
形状母数 (α): 2.0
尺度母数 (β): 3.0

#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 形状母数2.0、尺度母数3.0のガンマ分布を生成
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    // 計算したい点を設定
    double x = 4.0;
    // 確率密度関数の値を計算
    double density = boost::math::pdf(dist, x);
    // 結果を出力
    std::cout << "x = " << x << "におけるPDFの値: " << density << std::endl;
    return 0;
}

x = 4におけるPDFの値: 0.117154

#include <boost/math/distributions/gamma.hpp>
#include <iostream>
#include <vector>
int main() {
    // 形状母数3.0、尺度母数2.0のガンマ分布を生成
    boost::math::gamma_distribution<> dist(3.0, 2.0);
    // 評価したい点のリスト
    std::vector<double> points = {1.0, 2.0, 3.0, 4.0, 5.0};
    // 各点でのPDF値を計算し出力
    for (double x : points) {
        double density = boost::math::pdf(dist, x);
        std::cout << "x = " << x << " のときのPDF値: " << density << std::endl;
    }
    return 0;
}

x = 1 のときのPDF値: 0.0379082
x = 2 のときのPDF値: 0.0919699
x = 3 のときのPDF値: 0.125511
x = 4 のときのPDF値: 0.135335
x = 5 のときのPDF値: 0.128258

#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 形状母数2.0、尺度母数3.0のガンマ分布を生成
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    // 評価したい点
    double x = 4.0;
    // 累積分布関数の値を計算
    double cdf_value = boost::math::cdf(dist, x);
    // 結果を出力
    std::cout << "x = " << x << "におけるCDFの値: " << cdf_value << std::endl;
    return 0;
}

x = 4におけるCDFの値: 0.38494

#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 形状母数2.0、尺度母数3.0のガンマ分布を生成
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    // 例えば、0.75の確率に対応する値を求める
    double probability = 0.75;
    // 逆CDF(パーセンタイル)を計算
    double quantile_value = boost::math::quantile(dist, probability);
    // 結果を出力
    std::cout << "確率 " << probability << " に対応する分位点: " << quantile_value << std::endl;
    return 0;
}

確率 0.75 に対応する分位点: 8.0779

#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include <boost/random/variate_generator.hpp>
#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 乱数エンジンの生成
    boost::random::mt19937 engine;
    // 乱数エンジンを用いた一様分布の例(必要に応じて他の分布も使用可能)
    boost::random::uniform_real_distribution<> dist_uniform(0.0, 1.0);
    // 一様乱数を生成
    double uniform_sample = dist_uniform(engine);
    std::cout << "一様乱数: " << uniform_sample << std::endl;
    return 0;
}

一様乱数: 0.814724

#include <boost/math/distributions/gamma.hpp>
#include <iostream>

int main() {
    try {
        // 不正なパラメータ(負の値)を設定
        boost::math::gamma_distribution<> dist(-1.0, 2.0);
    } catch (const std::domain_error& e) {
        std::cerr << "ドメインエラー: " << e.what() << std::endl;
    }
    return 0;
}

double alpha = 2.0;
double beta = 3.0;
if (alpha <= 0 || beta <= 0) {
    std::cerr << "パラメータは正の値でなければなりません。" << std::endl;
    // 必要に応じて処理を中断またはデフォルト値に設定
}

try {
    boost::math::gamma_distribution<> dist(alpha, beta);
    double result = boost::math::pdf(dist, x);
    // 以降の処理
} catch (const boost::math::domain_error& e) {
    std::cerr << "パラメータまたは評価点が不正です: " << e.what() << std::endl;
    // 例外発生時の処理(例：デフォルト値の設定やエラー通知)
}

#include <boost/random/gamma_distribution.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <iostream>
#include <random>
#include <vector>

int main() {
    // 乱数エンジンの初期化
    std::random_device rd;
    boost::random::mt19937 engine(rd());

    // ガンマ分布のパラメータ設定
    double shape = 2.0;
    double scale = 3.0;

    // gamma_distributionを作成
    boost::random::gamma_distribution<> gamma_dist(shape, scale);

    // サンプル数
    const int sample_size = 10000;
    double sum = 0.0;
    double sum_sq = 0.0;

    // サンプルを生成し、平均と分散を計算
    for (int i = 0; i < sample_size; ++i) {
        double sample_value = gamma_dist(engine); // 直接エンジンに呼び出し
        sum += sample_value;
        sum_sq += sample_value * sample_value;
    }

    double mean = sum / sample_size;
    double variance = (sum_sq / sample_size) - (mean * mean);

    std::cout << "サンプル平均: " << mean << std::endl;
    std::cout << "サンプル分散: " << variance << std::endl;

    // 確率の推定例：値が5未満の確率
    int count_below_5 = 0;
    for (int i = 0; i < sample_size; ++i) {
        double sample_value = gamma_dist(engine);
        if (sample_value < 5.0) {
            ++count_below_5;
        }
    }

    double prob_below_5 = static_cast<double>(count_below_5) / sample_size;
    std::cout << "値が5未満の確率の推定: " << prob_below_5 << std::endl;

    return 0;
}

サンプル平均: 6.01233
サンプル分散: 17.6753
値が5未満の確率の推定: 0.4997

#include <boost/random/gamma_distribution.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/normal_distribution.hpp>
#include <iostream>
#include <random> // C++11
int main() {
    const int sample_size = 1000;
    double shape = 2.0;
    double scale = 2.0;

    // === C++11 <random> を使う場合 ===
    {
        std::mt19937 engine(12345); // 任意のシード
        std::gamma_distribution<double> gamma_dist(shape, scale);
        std::normal_distribution<double> normal_dist(0.0, 1.0);

        double sum = 0;
        for (int i = 0; i < sample_size; ++i) {
            double g = gamma_dist(engine);
            double n = normal_dist(engine);
            sum += (g + n);
        }
        std::cout << "[std::random] 平均値: " << (sum / sample_size) << "\n";
    }

    // === Boost.Random を使う場合 ===
    {
        boost::random::mt19937 engine(12345);
        boost::random::gamma_distribution<double> gamma_dist(shape, scale);
        boost::random::normal_distribution<double> normal_dist(0.0, 1.0);

        double sum = 0;
        for (int i = 0; i < sample_size; ++i) {
            double g = gamma_dist(engine);
            double n = normal_dist(engine);
            sum += (g + n);
        }
        std::cout << "[Boost.Random] 平均値: " << (sum / sample_size)
                  << std::endl;
    }

    return 0;
}

#include <boost/random/mersenne_twister.hpp>
#include <boost/random/variate_generator.hpp>
#include <boost/math/distributions/gamma.hpp>
#include <iostream>
int main() {
    // 分布オブジェクトを一度だけ作成
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    // 乱数エンジンの作成
    boost::random::mt19937 engine;
    // variate_generatorも一度だけ作成
    boost::variate_generator<boost::random::mt19937&, boost::math::gamma_distribution<> > sampler(engine, dist);
    // 何度もサンプリング
    for (int i = 0; i < 10000; ++i) {
        double sample = sampler();
        // サンプルの処理
    }
    return 0;
}

#include <boost/thread.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/variate_generator.hpp>
#include <boost/math/distributions/gamma.hpp>
#include <vector>
void thread_func() {
    // スレッドごとにエンジンを作成
    boost::random::mt19937 engine;
    boost::math::gamma_distribution<> dist(2.0, 3.0);
    boost::variate_generator<boost::random::mt19937&, boost::math::gamma_distribution<> > sampler(engine, dist);
    // 乱数生成
    for (int i = 0; i < 1000; ++i) {
        double sample = sampler();
        // 乱数を使った処理
    }
}
int main() {
    const int thread_count = 4;
    std::vector<boost::thread> threads;
    for (int i = 0; i < thread_count; ++i) {
        threads.emplace_back(thread_func);
    }
    for (auto& t : threads) {
        t.join();
    }
    return 0;
}