ハイスクールPython

10. Python を使った数値シミュレーション

目次

この講座で使用する Google Colab の URL

10. Python を使った数値シミュレーション (Google Colab)

演習課題

Ex.10. Python を使った数値シミュレーション (Google Colab)

この講座で使用する Python, Jupyter Notebook のファイルと実行環境

Lesson 10: high-school-python-code (GitHub)

10.1 シミュレーションとは?やり方とその意義

シミュレーションとは、現実の事象をコンピュータで再現することです。例えば、サイコロを 1000 回振ったらどうなるか、雨の日が続く確率はどのくらいかなど、実際に試すのが難しい場合でもコンピュータで簡単に実験できます。

シミュレーションの基本例

import random

# サイコロを100回振るシミュレーション
results = []
for _ in range(100):
    dice = random.randint(1, 6)
    results.append(dice)

# 結果の集計
for number in range(1, 7):
    count = results.count(number)
    percentage = count / len(results) * 100
    print(f"{number}の目: {count}回 ({percentage:.1f}%)")

10.2 乱数の基礎と random モジュールの使い方

random モジュールを使って、様々な種類の乱数を生成できます:

random モジュールの基本的な使い方

import random

# 整数の乱数(1から6まで)
dice = random.randint(1, 6)

# 0から1までの小数
probability = random.random()

# 範囲を指定した小数
temperature = random.uniform(20.0, 30.0)

# リストからランダムな要素を選択
colors = ["赤", "青", "緑", "黄"]
chosen = random.choice(colors)

# リストをシャッフル
cards = list(range(1, 53))
random.shuffle(cards)

10.3 確率分布とは(一様分布・二項分布・正規分布)

代表的な確率分布をシミュレーションで体験します:

様々な確率分布のシミュレーション

import random
import matplotlib.pyplot as plt
import numpy as np

# 一様分布(サイコロのような、全ての値が同じ確率で出る)
uniform_numbers = [random.uniform(0, 1) for _ in range(1000)]

# 二項分布(コイン投げを10回試行する実験を1000回)
def coin_flip_experiment():
    flips = [random.choice(['表', '裏']) for _ in range(10)]
    return flips.count('表')

binomial_results = [coin_flip_experiment() for _ in range(1000)]

# 正規分布(身長や体重のような、平均値を中心に分布する)
normal_numbers = np.random.normal(170, 5, 1000)  # 平均170、標準偏差5

# ヒストグラムの描画
plt.figure(figsize=(15, 5))

plt.subplot(1, 3, 1)
plt.hist(uniform_numbers, bins=30)
plt.title('一様分布')

plt.subplot(1, 3, 2)
plt.hist(binomial_results, bins=11)
plt.title('二項分布')

plt.subplot(1, 3, 3)
plt.hist(normal_numbers, bins=30)
plt.title('正規分布')

plt.tight_layout()
plt.show()

10.4 大数の法則:試行回数を増やすと平均値が一定に近づく

コイン投げを例に、大数の法則を確認します:

大数の法則のシミュレーション

import random
import matplotlib.pyplot as plt

def coin_flip_simulation(trials):
    heads = 0
    heads_ratio = []

    for i in range(1, trials + 1):
        if random.random() < 0.5:  # コインを投げる
            heads += 1
        heads_ratio.append(heads / i)  # 表の割合を記録

    return heads_ratio

# 異なる試行回数でシミュレーション
trials_1000 = coin_flip_simulation(1000)
trials_10000 = coin_flip_simulation(10000)

# グラフの描画
plt.figure(figsize=(12, 6))
plt.plot(trials_1000, label='1000回試行')
plt.plot(trials_10000, label='10000回試行')
plt.axhline(y=0.5, color='r', linestyle='--', label='理論値(0.5)')
plt.xlabel('試行回数')
plt.ylabel('表の出る確率')
plt.title('コイン投げシミュレーション:大数の法則の確認')
plt.legend()
plt.grid(True)
plt.show()

10.5 モンテカルロ法:乱数を使って円周率を求めてみよう

モンテカルロ法を使って円周率を計算します:

モンテカルロ法で円周率を計算

import random
import matplotlib.pyplot as plt

def estimate_pi(points):
    inside_circle = 0
    x_inside, y_inside = [], []
    x_outside, y_outside = [], []

    for _ in range(points):
        x = random.uniform(-1, 1)
        y = random.uniform(-1, 1)

        # 点が単位円の中にあるか判定
        if x*x + y*y <= 1:
            inside_circle += 1
            x_inside.append(x)
            y_inside.append(y)
        else:
            x_outside.append(x)
            y_outside.append(y)

    # 円周率の計算
    pi_estimate = 4 * inside_circle / points

    return pi_estimate, x_inside, y_inside, x_outside, y_outside

# シミュレーション実行
points = 10000
pi_estimate, x_in, y_in, x_out, y_out = estimate_pi(points)

# 結果の表示
print(f"計算された円周率: {pi_estimate}")
print(f"実際の円周率との誤差: {abs(pi_estimate - 3.14159265359)}")

# 点のプロット
plt.figure(figsize=(8, 8))
plt.scatter(x_in, y_in, c='blue', s=1, alpha=0.6, label='円の内側')
plt.scatter(x_out, y_out, c='red', s=1, alpha=0.6, label='円の外側')
plt.axis('equal')
plt.title(f'モンテカルロ法による円周率の計算\n推定値: {pi_estimate:.6f}')
plt.legend()
plt.grid(True)
plt.show()

10.6 色々なシミュレーションをしてみよう

実践的なシミュレーション例を見てみましょう:

じゃんけんのシミュレーション

import random
from collections import Counter

def play_janken():
    hands = ['グー', 'チョキ', 'パー']
    rules = {
        ('グー', 'チョキ'): 'グー',
        ('チョキ', 'パー'): 'チョキ',
        ('パー', 'グー'): 'パー',
    }

    player1 = random.choice(hands)
    player2 = random.choice(hands)

    if player1 == player2:
        return '引き分け'

    # 勝敗判定
    if (player1, player2) in rules:
        if rules[(player1, player2)] == player1:
            return 'プレイヤー1の勝ち'
    else:
        if rules[(player2, player1)] == player2:
            return 'プレイヤー2の勝ち'

# 1000回のじゃんけんをシミュレーション
results = []
for _ in range(1000):
    result = play_janken()
    results.append(result)

# 結果の集計
counter = Counter(results)
for outcome, count in counter.items():
    percentage = count / len(results) * 100
    print(f"{outcome}: {count}回 ({percentage:.1f}%)")

10.7 グラフで結果を可視化して考察する

シミュレーション結果を視覚的に表現します:

サイコロの期待値のシミュレーション

import random
import matplotlib.pyplot as plt
import numpy as np

def dice_simulation(trials):
    # サイコロを振って平均値を記録
    averages = []
    running_sum = 0

    for i in range(1, trials + 1):
        roll = random.randint(1, 6)
        running_sum += roll
        average = running_sum / i
        averages.append(average)

    return averages

# 複数回のシミュレーション
num_simulations = 5
trials = 1000

plt.figure(figsize=(12, 6))

# 各シミュレーションを実行してプロット
for i in range(num_simulations):
    averages = dice_simulation(trials)
    plt.plot(averages, alpha=0.5, label=f'シミュレーション{i+1}')

# 理論値(3.5)の線を追加
plt.axhline(y=3.5, color='r', linestyle='--', label='理論値(3.5)')

plt.xlabel('試行回数')
plt.ylabel('平均値')
plt.title('サイコロの平均値の収束')
plt.legend()
plt.grid(True)
plt.show()

Ex.1 コインを 100 回投げて表と裏の回数を数えるプログラム

コイン投げシミュレーション

import random
from collections import Counter
import matplotlib.pyplot as plt

def coin_flip_experiment(flips):
    # コインを指定回数投げる
    results = []
    for _ in range(flips):
        result = 'タ' if random.random() < 0.5 else 'オ'
        results.append(result)

    # 結果を集計
    count = Counter(results)

    # 結果を表示
    print(f"\n{flips}回コインを投げた結果:")
    for face, num in count.items():
        percentage = (num / flips) * 100
        print(f"{face}: {num}回 ({percentage:.1f}%)")

    # グラフの作成
    plt.figure(figsize=(8, 6))
    plt.bar(count.keys(), count.values())
    plt.title(f'コイン投げの結果 ({flips}回)')
    plt.ylabel('回数')
    plt.show()

    return count

# プログラムの実行
result = coin_flip_experiment(100)

# チャレンジ:
# - 複数回実験を行い、結果の平均を取る
# - 試行回数を変えて実験する
# - アニメーションで結果をリアルタイム表示する

Ex.2 サイコロを繰り返し振って、各目の出る確率が 1/6 に近づく様子をグラフ化するプログラム

サイコロの確率収束シミュレーション

import random
import matplotlib.pyplot as plt
import numpy as np

def dice_probability_simulation(trials):
    # 結果を記録する配列
    counts = {i: 0 for i in range(1, 7)}
    probabilities = {i: [] for i in range(1, 7)}

    # サイコロを振るシミュレーション
    for i in range(1, trials + 1):
        roll = random.randint(1, 6)
        counts[roll] += 1

        # 各目の確率を計算して記録
        for num in range(1, 7):
            prob = counts[num] / i
            probabilities[num].append(prob)

    # グラフの描画
    plt.figure(figsize=(12, 6))

    # 各目の確率の推移をプロット
    for num in range(1, 7):
        plt.plot(probabilities[num],
                label=f'{num}の目',
                alpha=0.7)

    # 理論値(1/6)の線を追加
    plt.axhline(y=1/6, color='r', linestyle='--',
                label='理論値(1/6)')

    plt.xlabel('試行回数')
    plt.ylabel('確率')
    plt.title('サイコロの確率収束シミュレーション')
    plt.legend()
    plt.grid(True)
    plt.show()

    # 最終的な確率を表示
    print("\n最終的な確率:")
    for num in range(1, 7):
        final_prob = probabilities[num][-1]
        print(f"{num}の目: {final_prob:.3f}")

# プログラムの実行
dice_probability_simulation(1000)

# チャレンジ:
# - 複数のサイコロを同時に投げる
# - 目の合計の分布を調べる
# - 特定の目が出る確率を予測する

まとめ

この章で学んだこと:

  • シミュレーションの基本的な考え方
  • random モジュールの使い方
  • 様々な確率分布の特徴
  • 大数の法則の実践的な理解
  • モンテカルロ法による数値計算
  • データの可視化と結果の考察方法