C 言語でモンテカルロシミュレーションを用いた確率計算

CBeginner
オンラインで実践に進む

はじめに

この実験では、C プログラミングでモンテカルロシミュレーションを使用して確率を推定する方法を学びます。コイン投げのような単純なランダム実験から始め、複数の試行を実行して成功した結果の数をカウントします。最後に、成功回数と総試行回数で割ることで、推定確率を計算します。この実験は、データ分析や意思決定に不可欠な、確率と組み合わせ論の基礎概念を C を使って実践的に学ぶ導入です。

ランダム実験の定義

このステップでは、C 言語でモンテカルロシミュレーションを用いてランダム実験を定義する方法を学びます。ランダム実験とは、確率的手法を用いてシミュレーションできる、不確かな結果を持つプロセスです。

ランダム実験の理解

基本的なランダム実験であるコイン投げシミュレーションを示す簡単な C プログラムを作成しましょう。

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

#define NUM_TRIALS 1000

int main() {
    // 乱数生成器のシードを設定
    srand(time(NULL));

    // 表が出た回数
    int heads_count = 0;

    // コイン投げをシミュレート
    for (int i = 0; i < NUM_TRIALS; i++) {
        // 0 または 1 の乱数を生成
        int flip = rand() % 2;

        // 表をカウント
        if (flip == 0) {
            heads_count++;
        }
    }

    // 表が出る確率を計算
    double probability = (double)heads_count / NUM_TRIALS;

    printf("コイン投げ実験:\n");
    printf("総試行回数:%d\n", NUM_TRIALS);
    printf("表が出た回数:%d\n", heads_count);
    printf("表が出る確率の推定値:%.2f\n", probability);

    return 0;
}

出力例:

コイン投げ実験:
総試行回数: 1000
表が出た回数: 502
表が出る確率の推定値: 0.50

重要な概念の説明

  1. 乱数生成:

    • srand(time(NULL)) は乱数生成器のシードを設定します。
    • rand() % 2 は、0 または 1 を等しい確率で生成します。
  2. 実験設計:

    • コイン投げをランダム実験として定義します。
    • 複数の試行(ここでは 1000 回)を実行します。
    • 成功した結果(表)の数をカウントします。
  3. 確率の推定:

    • 確率 = (成功した結果の数) / (総試行回数)
    • この場合、表が出る確率は 0.5 に近い値になることが期待されます。

プログラムのコンパイルと実行

## ソースファイルを作成
nano ~/project/coin_flip_experiment.c

## プログラムをコンパイル
gcc ~/project/coin_flip_experiment.c -o ~/project/coin_flip_experiment

## 実験を実行
~/project/coin_flip_experiment

コンパイルと実行の出力例:

## コンパイル
gcc ~/project/coin_flip_experiment.c -o ~/project/coin_flip_experiment

## 実行
~/project/coin_flip_experiment
コイン投げ実験:
総試行回数:1000
表が出た回数:502
表が出る確率の推定値:0.50

複数回のランダム試行の実行と成功カウント

このステップでは、モンテカルロシミュレーションを拡張し、複数のランダム試行を実行して成功した結果を正確にカウントします。より複雑な確率実験、ランダムな点の生成を用いたπ(パイ)の推定を通して、これを示します。

π推定のためのモンテカルロ法

四分円を含む正方形内にランダムな点を生成することで、幾何学的アプローチを用いてπを推定します。

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>

#define NUM_TRIALS 100000

int main() {
    // 乱数生成器のシードを設定
    srand(time(NULL));

    // 総点と円内点のカウント
    int total_points = 0;
    int points_inside_circle = 0;

    // 複数回の試行を実行
    for (int i = 0; i < NUM_TRIALS; i++) {
        // 0 から 1 の間のランダムな x 座標と y 座標を生成
        double x = (double)rand() / RAND_MAX;
        double y = (double)rand() / RAND_MAX;

        // 点が四分円の内部にあるかどうかをチェック
        if (sqrt(x*x + y*y) <= 1.0) {
            points_inside_circle++;
        }
        total_points++;
    }

    // πの推定
    double pi_estimate = 4.0 * points_inside_circle / total_points;

    printf("π推定実験:\n");
    printf("総点:%d\n", total_points);
    printf("円内点:%d\n", points_inside_circle);
    printf("推定π: %.6f\n", pi_estimate);
    printf("実際のπ: %.6f\n", M_PI);

    return 0;
}

実験のコンパイルと実行

## ソースファイルを作成
nano ~/project/pi_estimation.c

## プログラムをコンパイル(mathライブラリのリンクに注意)
gcc ~/project/pi_estimation.c -o ~/project/pi_estimation -lm

## 実験を実行
~/project/pi_estimation

出力例:

π推定実験:
総点: 100000
円内点: 78540
推定π: 3.141600
実際のπ: 3.141593

重要な概念の説明

  1. 複数回の試行:

    • 多くのランダム試行(10 万回)を実行します。
    • 各試行で 1x1 の正方形内のランダムな点を生成します。
  2. 成功カウント:

    • 総点と四分円内の点を追跡します。
    • 成功は、点が円内に落ちることを意味します。
  3. 確率の推定:

    • 確率 = (成功した点の数) / (総点の数)
    • 四分円法のため、4 を掛けてπを推定します。

重要なテクニック

  • (double)rand() / RAND_MAX は 0 から 1 の間のランダムな小数を生成します。
  • sqrt(x*x + y*y) は原点からの距離を計算します。
  • 試行回数を増やすことで、推定の精度を高めることができます。

確率の推定 = 成功回数/試行回数

この最終ステップでは、より実践的なシナリオで、成功した結果の割合と総試行回数を分析することで、確率を計算する方法を示します。

サイコロ振り確率実験

2 個のサイコロを振って、特定の合計値(例:合計 7)を得る確率を計算します。

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

#define NUM_TRIALS 100000
#define TARGET_SUM 7

int roll_die() {
    return rand() % 6 + 1;
}

int main() {
    // 乱数生成器のシードを設定
    srand(time(NULL));

    // 総試行回数と成功試行回数のカウント
    int total_rolls = 0;
    int successful_rolls = 0;

    // 複数回の試行を実行
    for (int i = 0; i < NUM_TRIALS; i++) {
        // 2 個のサイコロを振る
        int die1 = roll_die();
        int die2 = roll_die();

        // 合計が目標値と一致するかどうかをチェック
        if (die1 + die2 == TARGET_SUM) {
            successful_rolls++;
        }
        total_rolls++;
    }

    // 確率を計算
    double probability = (double)successful_rolls / total_rolls;

    printf("サイコロ振り確率実験:\n");
    printf("目標合計:%d\n", TARGET_SUM);
    printf("総試行回数:%d\n", total_rolls);
    printf("成功試行回数:%d\n", successful_rolls);
    printf("推定確率:%.4f\n", probability);

    // 比較のための理論確率
    printf("理論確率:%.4f\n", 1.0/6);

    return 0;
}

実験のコンパイルと実行

## ソースファイルを作成
nano ~/project/dice_probability.c

## プログラムをコンパイル
gcc ~/project/dice_probability.c -o ~/project/dice_probability

## 実験を実行
~/project/dice_probability

出力例:

サイコロ振り確率実験:
目標合計: 7
総試行回数: 100000
成功試行回数: 16644
推定確率: 0.1664
理論確率: 0.1667

重要な概念の説明

  1. 確率の計算:

    • 確率 = (成功した結果の数) / (総試行回数)
    • この場合:成功試行回数 / 総試行回数
  2. モンテカルロシミュレーション:

    • 多くの試行回数(10 万回)により、正確な推定が可能になります。
    • シミュレーションされた確率は、理論確率と非常に近い値になります。
  3. ランダム性と精度:

    • srand(time(NULL)) は異なるランダムなシーケンスを保証します。
    • 試行回数を増やすことで、推定の精度を高めることができます。

確率の解釈

  • 推定確率(0.1664)は、理論確率(1/6 ≈ 0.1667)と非常に近いです。
  • モンテカルロ法が確率を推定する方法を示しています。

まとめ

この実験では、C 言語を用いてモンテカルロシミュレーションでランダムな実験を定義する方法を学びました。重要な概念を示すために、単純なコイン投げシミュレーションを作成しました。最初に、乱数生成器のシードを設定し、コイン投げをシミュレーションして、成功した結果(表)の数をカウントしました。次に、成功した結果の数と総試行回数で割ることで、表が出る確率を推定しました。出力結果から、公平なコイン投げの場合の期待値である 0.5 に近い推定確率が表示されました。