NumPy乱数配列の完全マスター!様々な種類の乱数生成テクニック集

 

NumPyの乱数生成機能は、データサイエンス、機械学習、統計解析において欠かせないツールです。本記事では、様々な種類の乱数配列の生成方法を、基本から応用まで豊富なサンプルコードとともに詳しく解説します。

NumPy乱数生成の基本構造

NumPyには旧式のrandom関数新式のGeneratorの2つの方式があります。

import numpy as np

# 旧式(NumPy 1.17以前の方式)
old_random = np.random.rand(5)

# 新式(推奨)
rng = np.random.default_rng(seed=42)
new_random = rng.random(5)

print("旧式:", old_random)
print("新式:", new_random)

新式Generatorを推奨する理由:

  • より高性能
  • 予測可能な結果
  • 並列処理に安全
  • より豊富な分布

1. 基本的な乱数配列の生成

1-1. 一様乱数(0~1)

# [0, 1)の一様乱数
rng = np.random.default_rng(42)

# 1次元配列
uniform_1d = rng.random(5)
print("1次元:", uniform_1d)

# 2次元配列
uniform_2d = rng.random((3, 4))
print("2次元形状:", uniform_2d.shape)

# 3次元配列
uniform_3d = rng.random((2, 3, 4))
print("3次元形状:", uniform_3d.shape)

1-2. 指定範囲の一様乱数

# [low, high)の一様乱数
low, high = 10, 50
uniform_range = rng.uniform(low, high, 8)
print(f"[{low}, {high})の乱数:", uniform_range)

# 整数の一様乱数
int_range = rng.integers(1, 100, 10)
print("整数乱数:", int_range)

# 浮動小数点の範囲指定
float_range = rng.uniform(-5.0, 5.0, 6)
print("浮動小数点:", float_range)

1-3. 正規分布(ガウス分布)

# 標準正規分布(平均0、標準偏差1)
normal_std = rng.standard_normal(8)
print("標準正規分布:", normal_std)

# 平均と標準偏差を指定
mean, std = 100, 15
normal_custom = rng.normal(mean, std, 10)
print(f"正規分布(μ={mean}, σ={std}):", normal_custom)

# 2次元の正規分布
normal_2d = rng.normal(0, 1, (3, 5))
print("2次元正規分布形状:", normal_2d.shape)

2. 様々な確率分布からの乱数生成

2-1. 指数分布

# 指数分布(λ=1/scale)
scale = 2.0
exponential = rng.exponential(scale, 8)
print(f"指数分布(scale={scale}):", exponential)

# 複数のスケールパラメータ
scales = [1, 2, 5]
for s in scales:
    exp_sample = rng.exponential(s, 3)
    print(f"scale={s}: {exp_sample}")

2-2. ポアソン分布

# ポアソン分布(離散分布)
lam = 3.0  # 平均発生率
poisson = rng.poisson(lam, 10)
print(f"ポアソン分布(λ={lam}):", poisson)

# 異なる発生率での比較
lambdas = [1, 5, 10]
for l in lambdas:
    pois_sample = rng.poisson(l, 5)
    print(f"λ={l}: {pois_sample}")

2-3. 二項分布

# 二項分布(n回試行、確率p)
n, p = 20, 0.3
binomial = rng.binomial(n, p, 8)
print(f"二項分布(n={n}, p={p}):", binomial)

# コイン投げのシミュレーション
coin_flips = rng.binomial(1, 0.5, 100)  # 100回コイン投げ
heads_ratio = coin_flips.mean()
print(f"表の割合: {heads_ratio:.2f}")

2-4. ガンマ分布

# ガンマ分布
shape, scale = 2.0, 1.0
gamma = rng.gamma(shape, scale, 8)
print(f"ガンマ分布(shape={shape}):", gamma)

# カイ二乗分布(ガンマ分布の特殊ケース)
df = 5  # 自由度
chi_square = rng.chisquare(df, 8)
print(f"カイ二乗分布(df={df}):", chi_square)

2-5. ベータ分布

# ベータ分布([0,1]区間)
alpha, beta = 2.0, 5.0
beta_dist = rng.beta(alpha, beta, 8)
print(f"ベータ分布(α={alpha}, β={beta}):", beta_dist)

# 一様分布(ベータ分布の特殊ケース)
uniform_beta = rng.beta(1, 1, 5)
print("一様分布(ベータ版):", uniform_beta)

3. 多変量分布の乱数生成

3-1. 多変量正規分布

# 2次元多変量正規分布
mean = [0, 0]
cov = [[1, 0.5], [0.5, 2]]  # 共分散行列
multivariate = rng.multivariate_normal(mean, cov, 5)
print("多変量正規分布:")
print(multivariate)

# 相関のない場合
uncorr_cov = [[1, 0], [0, 1]]
uncorr_mv = rng.multivariate_normal([0, 0], uncorr_cov, 3)
print("無相関多変量正規分布:")
print(uncorr_mv)

3-2. ディリクレ分布

# ディリクレ分布(確率の分布)
alpha = [1, 2, 3]
dirichlet = rng.dirichlet(alpha, 5)
print("ディリクレ分布:")
print(dirichlet)

# 各行の和が1であることを確認
row_sums = dirichlet.sum(axis=1)
print("行の和:", row_sums)

4. 離散的な選択とサンプリング

4-1. 配列からのランダムサンプリング

# 配列からの復元抽出
data = np.array([10, 20, 30, 40, 50])
sample_with_replace = rng.choice(data, 8, replace=True)
print("復元抽出:", sample_with_replace)

# 非復元抽出
sample_without_replace = rng.choice(data, 3, replace=False)
print("非復元抽出:", sample_without_replace)

# 確率重み付きサンプリング
weights = [0.1, 0.2, 0.3, 0.3, 0.1]
weighted_sample = rng.choice(data, 6, p=weights, replace=True)
print("重み付きサンプリング:", weighted_sample)

4-2. インデックスのランダム選択

# インデックスのシャッフル
indices = np.arange(10)
shuffled = indices.copy()
rng.shuffle(shuffled)
print("シャッフル前:", indices)
print("シャッフル後:", shuffled)

# ランダムな順列
permutation = rng.permutation(8)
print("ランダム順列:", permutation)

# 配列の順列
arr = np.array(['A', 'B', 'C', 'D', 'E'])
perm_arr = rng.permutation(arr)
print("配列の順列:", perm_arr)

5. 特殊な乱数パターンの生成

5-1. 条件付き乱数生成

def generate_conditional_random(n, condition_func, max_attempts=1000):
    """条件を満たす乱数のみを生成"""
    rng = np.random.default_rng()
    result = []
    attempts = 0
    
    while len(result) < n and attempts < max_attempts:
        candidate = rng.standard_normal()
        if condition_func(candidate):
            result.append(candidate)
        attempts += 1
    
    return np.array(result)

# 正の値のみの正規乱数
positive_normal = generate_conditional_random(8, lambda x: x > 0)
print("正の正規乱数:", positive_normal)

# 絶対値が1以上の乱数
large_normal = generate_conditional_random(5, lambda x: abs(x) >= 1)
print("絶対値≥1の乱数:", large_normal)

5-2. 相関のある乱数系列

def generate_correlated_series(n, correlation=0.7):
    """相関のある時系列乱数を生成"""
    rng = np.random.default_rng(42)
    series = np.zeros(n)
    series[0] = rng.standard_normal()
    
    for i in range(1, n):
        series[i] = correlation * series[i-1] + \
                   np.sqrt(1 - correlation**2) * rng.standard_normal()
    
    return series

# 相関のある時系列
corr_series = generate_correlated_series(20, 0.8)
print("相関時系列(最初の5個):", corr_series[:5])

# 実際の相関計算
actual_corr = np.corrcoef(corr_series[:-1], corr_series[1:])[0,1]
print(f"実際の相関係数: {actual_corr:.3f}")

5-3. ノイズの重畳

def add_noise(signal, noise_level=0.1, noise_type='gaussian'):
    """信号にノイズを追加"""
    rng = np.random.default_rng()
    
    if noise_type == 'gaussian':
        noise = rng.normal(0, noise_level, signal.shape)
    elif noise_type == 'uniform':
        noise = rng.uniform(-noise_level, noise_level, signal.shape)
    elif noise_type == 'salt_pepper':
        noise = np.zeros_like(signal)
        # 塩胡椒ノイズ
        salt = rng.random(signal.shape) < noise_level/2
        pepper = rng.random(signal.shape) < noise_level/2
        noise[salt] = 1
        noise[pepper] = -1
    
    return signal + noise

# 正弦波にノイズを追加
t = np.linspace(0, 2*np.pi, 50)
clean_signal = np.sin(t)
noisy_signal = add_noise(clean_signal, 0.2, 'gaussian')

print("ノイズ追加前(最初の3点):", clean_signal[:3])
print("ノイズ追加後(最初の3点):", noisy_signal[:3])

6. 実践的な応用例

6-1. モンテカルロシミュレーション

def monte_carlo_pi(n_samples):
    """モンテカルロ法でπを推定"""
    rng = np.random.default_rng(42)
    
    # 単位正方形内のランダム点
    x = rng.uniform(-1, 1, n_samples)
    y = rng.uniform(-1, 1, n_samples)
    
    # 単位円内の点を数える
    inside_circle = (x**2 + y**2) <= 1
    pi_estimate = 4 * inside_circle.sum() / n_samples
    
    return pi_estimate

# π推定の精度確認
for n in [1000, 10000, 100000]:
    pi_est = monte_carlo_pi(n)
    error = abs(pi_est - np.pi)
    print(f"n={n:6d}: π≈{pi_est:.4f}, 誤差={error:.4f}")

6-2. ブートストラップサンプリング

def bootstrap_confidence_interval(data, statistic_func, n_bootstrap=1000, confidence=0.95):
    """ブートストラップ信頼区間"""
    rng = np.random.default_rng(42)
    n = len(data)
    
    bootstrap_stats = []
    for _ in range(n_bootstrap):
        # 復元抽出でブートストラップサンプル
        bootstrap_sample = rng.choice(data, n, replace=True)
        stat = statistic_func(bootstrap_sample)
        bootstrap_stats.append(stat)
    
    bootstrap_stats = np.array(bootstrap_stats)
    
    # 信頼区間の計算
    alpha = 1 - confidence
    lower = np.percentile(bootstrap_stats, 100 * alpha/2)
    upper = np.percentile(bootstrap_stats, 100 * (1 - alpha/2))
    
    return lower, upper, bootstrap_stats

# 平均の信頼区間
data = rng.normal(50, 10, 100)
ci_lower, ci_upper, boot_means = bootstrap_confidence_interval(
    data, np.mean, confidence=0.95
)

print(f"データ平均: {data.mean():.2f}")
print(f"95%信頼区間: [{ci_lower:.2f}, {ci_upper:.2f}]")

6-3. 異常値の生成

def generate_data_with_outliers(n_normal, n_outliers, outlier_factor=5):
    """異常値を含むデータセットを生成"""
    rng = np.random.default_rng(42)
    
    # 正常データ
    normal_data = rng.normal(0, 1, n_normal)
    
    # 異常値
    outlier_scale = outlier_factor
    outliers = rng.choice([-1, 1], n_outliers) * \
               rng.uniform(outlier_scale, outlier_scale*2, n_outliers)
    
    # 結合
    all_data = np.concatenate([normal_data, outliers])
    labels = np.concatenate([np.zeros(n_normal), np.ones(n_outliers)])
    
    # シャッフル
    indices = rng.permutation(len(all_data))
    
    return all_data[indices], labels[indices]

# 異常値付きデータセット
data, labels = generate_data_with_outliers(95, 5)
outlier_indices = np.where(labels == 1)[0]

print(f"データサイズ: {len(data)}")
print(f"異常値の位置: {outlier_indices}")
print(f"異常値: {data[outlier_indices]}")

7. 再現性とシード管理

7-1. 基本的なシード設定

# グローバルシード設定(旧式)
np.random.seed(42)
old_way = np.random.rand(3)

# Generatorでのシード設定(推奨)
rng1 = np.random.default_rng(42)
rng2 = np.random.default_rng(42)

new_way1 = rng1.random(3)
new_way2 = rng2.random(3)

print("旧式:", old_way)
print("新式1:", new_way1)
print("新式2:", new_way2)
print("新式同じ結果:", np.array_equal(new_way1, new_way2))

7-2. 複数のランダムストリーム

def create_multiple_streams(master_seed, n_streams):
    """独立した複数の乱数ストリームを作成"""
    # マスターシードから子シードを生成
    master_rng = np.random.default_rng(master_seed)
    seeds = master_rng.integers(0, 2**32, n_streams)
    
    # 独立したGenerator群
    streams = [np.random.default_rng(seed) for seed in seeds]
    
    return streams

# 3つの独立ストリーム
streams = create_multiple_streams(42, 3)

for i, stream in enumerate(streams):
    sample = stream.normal(0, 1, 3)
    print(f"ストリーム{i}: {sample}")

7-3. 状態の保存と復元

# 乱数生成器の状態を保存
rng = np.random.default_rng(42)
initial_state = rng.bit_generator.state

# いくつか乱数を生成
sample1 = rng.random(3)
print("最初のサンプル:", sample1)

# 状態を復元
rng.bit_generator.state = initial_state
sample2 = rng.random(3)
print("復元後のサンプル:", sample2)
print("同じ結果:", np.array_equal(sample1, sample2))

8. パフォーマンス最適化

8-1. 大規模配列の効率的生成

import time

def benchmark_random_generation(size):
    """乱数生成のベンチマーク"""
    rng = np.random.default_rng()
    
    # 単精度 vs 倍精度
    start = time.time()
    float32_array = rng.random(size, dtype=np.float32)
    time_float32 = time.time() - start
    
    start = time.time()
    float64_array = rng.random(size, dtype=np.float64)
    time_float64 = time.time() - start
    
    print(f"サイズ: {size:,}")
    print(f"float32: {time_float32:.4f}秒")
    print(f"float64: {time_float64:.4f}秒")
    print(f"メモリ比: {float32_array.nbytes / float64_array.nbytes:.1f}")

# ベンチマーク実行
benchmark_random_generation(1_000_000)

8-2. チャンク処理による大規模データ生成

def generate_large_random_chunked(total_size, chunk_size, distribution_func):
    """チャンク単位での大規模乱数生成"""
    rng = np.random.default_rng()
    
    result = np.empty(total_size)
    processed = 0
    
    while processed < total_size:
        current_chunk = min(chunk_size, total_size - processed)
        chunk = distribution_func(rng, current_chunk)
        result[processed:processed + current_chunk] = chunk
        processed += current_chunk
        
        if processed % (chunk_size * 10) == 0:
            print(f"進捗: {processed:,} / {total_size:,}")
    
    return result

# 大規模正規分布データの生成
def normal_generator(rng, size):
    return rng.normal(0, 1, size)

# 使用例(小さいサイズでデモ)
large_data = generate_large_random_chunked(
    1000, 200, normal_generator
)
print(f"生成完了: {len(large_data):,}個の要素")

9. よくある問題と解決方法

9-1. シード設定の落とし穴

# 問題: グローバルシードの干渉
def problematic_function():
    np.random.seed(42)  # グローバル状態を変更
    return np.random.rand(3)

def better_function(seed=None):
    rng = np.random.default_rng(seed)  # 独立したGenerator
    return rng.random(3)

# テスト
print("問題のある方法:")
result1 = problematic_function()
result2 = problematic_function()
print("同じ結果:", np.array_equal(result1, result2))

print("\n改善された方法:")
result3 = better_function(42)
result4 = better_function(42)
print("同じ結果:", np.array_equal(result3, result4))

9-2. 分布パラメータの検証

def safe_normal_generation(mean, std, size, min_std=1e-6):
    """安全な正規分布生成(パラメータ検証付き)"""
    rng = np.random.default_rng()
    
    # パラメータ検証
    if std <= 0:
        raise ValueError(f"標準偏差は正の値である必要があります: {std}")
    
    if std < min_std:
        print(f"警告: 標準偏差が非常に小さいです: {std}")
    
    if not isinstance(size, (int, tuple)):
        raise TypeError("sizeは整数またはタプルである必要があります")
    
    return rng.normal(mean, std, size)

# 使用例
try:
    safe_result = safe_normal_generation(0, 1, 5)
    print("正常生成:", safe_result)
    
    # エラーケース
    error_result = safe_normal_generation(0, -1, 5)
except ValueError as e:
    print(f"エラー捕捉: {e}")

まとめ

NumPyの乱数生成機能を使いこなすことで、様々な統計・機械学習タスクが効率的に実行できます。

重要なポイント:

  • 新式Generatordefault_rng())を使用
  • 適切な分布を選択(一様、正規、指数など)
  • 再現性のためのシード管理
  • 大規模データでのメモリ・速度最適化
  • パラメータ検証でエラーを防止

これらの技術を活用して、データサイエンスプロジェクトを成功に導きましょう!

■プロンプトだけでオリジナルアプリを開発・公開してみた!!

■AI時代の第一歩!「AI駆動開発コース」はじめました!

テックジム東京本校で先行開始。

■テックジム東京本校

「武田塾」のプログラミング版といえば「テックジム」。
講義動画なし、教科書なし。「進捗管理とコーチング」で効率学習。
より早く、より安く、しかも対面型のプログラミングスクールです。

<短期講習>5日で5万円の「Pythonミニキャンプ」開催中。

<月1開催>放送作家による映像ディレクター養成講座

<オンライン無料>ゼロから始めるPython爆速講座