Pythonで指示確率変数から関数Fの期待値を計算するプログラムの実装方法

問題の概要

2つの整数 k と n が与えられているとします。1 から n までの自然数 {1, 2, ..., n} のランダムな順列 p₁, p₂, ..., p_n を考え、次の値 F を計算します。

F = (X₂ + ... + X_n−1)^k

ここで X_i は指示確率変数(インジケーター確率変数)であり、以下のいずれかの条件が成り立つときに 1、それ以外のときに 0 を取ります。

• p_i−1 < p_i > p_i+1(局所的な極大となる場合)
• p_i−1 > p_i < p_i+1(局所的な極小となる場合)

このとき、F の期待値 E[F] を求めることが目標です。

たとえば、入力が k = 1、n = 1000 の場合、出力は 1996/3 になります。

解法のアプローチ

この問題は、k の値ごとに事前に導出された多項式(閉形式)を利用することで、効率的に答えを求められます。全体の流れは次の通りです。

1. 関数 exp_factor() を定義します。引数として n と k を受け取ります。
2. k の値に応じて、対応する多項式の「分子」と「分母」のペアを返します。
3. メイン処理では、分子と分母の最大公約数(GCD)で約分し、既約分数の形で結果を返します。

exp_factor() の戻り値一覧

• k = 1 の場合:(2×(n−2), 3) を返す
• k = 2 の場合:(40n² − 144n + 131, 90) を返す
• k = 3 の場合:(280n³ − 1344n² + 2063n − 1038, 945) を返す
• k = 4 の場合:(2800n⁴ − 15680n³ + 28844n² − 19288n + 4263, 14175) を返す
• k = 5 の場合:(12320n⁵ − 73920n⁴ + 130328n³ − 29568n² − 64150n − 5124, 93555) を返す
• 上記以外の場合:1.0 を返す

メイン処理の手順

• M := n − 2
• p := 2.0 / 3
• q := 1 − p
• (num, den) := exp_factor(n, k)
• g := gcd(num, den)
• 分数 (num/g) / (den/g) を文字列として返す

なぜ p = 2/3 なのか?

連続する 3 つの要素 p_i−1, p_i, p_i+1 は互いに異なる値を持つため、その大小関係の並び方は 3! = 6 通りあり、すべて同じ確率で現れます。このうち、中央の要素が極大または極小になる(山または谷の形になる)のは 6 通りのうち 4 通りです。したがって、P(X_i = 1) = 4/6 = 2/3 となります。これが p = 2/3 という設定の由来です。

実装例(Python)

以下のコードで実際の動作を確認できます。

from math import gcd

def exp_factor(n, k):
    if k == 1:
        return (2*(n-2), 3)
    elif k == 2:
        return (40*n**2 - 144*n + 131, 90)
    elif k == 3:
        return (280*n**3 - 1344*n**2 + 2063*n - 1038, 945)
    elif k == 4:
        return (2800*n**4 - 15680*n**3 + 28844*n**2 - 19288*n + 4263, 14175)
    elif k == 5:
        return (12320*n**5 - 73920*n**4 + 130328*n**3 - 29568*n**2 - 64150*n - 5124, 93555)
    return 1.0

def solve(k, n):
    M = n - 2
    p = 2.0 / 3
    q = 1 - p

    num, den = exp_factor(n, k)
    g = gcd(num, den)
    return str(int(num/g)) + '/' + str(int(den/g))

k = 1
n = 1000
print(solve(k, n))

入力

1, 1000

出力

1996/3

まとめ

このプログラムでは、順列の中で各位置が局所的な極大・極小になる事象を指示確率変数で表し、その k 乗和の期待値を k ごとの閉形式多項式で計算しています。GCD による約分を行うことで、結果を常に既約分数として出力できる点がポイントです。