すでに知っていることを使う
コールをシミュレートするとき、満期の株価 S_T も同時にシミュレートしており、その真の平均 S0 * exp(r * T) は正確に分かっています。標本の S_T の平均が少しずれて、たとえば高すぎたとします。コールのペイオフと S_T は連動するので、コールの標本平均もおそらく高すぎます。コントロール変量法はこれを使い、知りたい量の推定値から、既知の量で観測された誤差の何倍かを引きます。
import numpy as np
rng = np.random.default_rng(1)
z = rng.standard_normal(2000)
s_t = 100.0 * np.exp((0.05 - 0.5 * 0.2 ** 2) + 0.2 * z) # terminal prices
call = np.maximum(s_t - 100.0, 0.0) # the payoff we want
known = 100.0 * np.exp(0.05) # E[S_T], exactly
print(round(float(s_t.mean()), 3), round(known, 3)) # the sample mean of S_T vs its true mean
print(round(float(s_t.mean() - known), 3)) # the observed error of the control最適な倍率
コントロールの誤差のどれだけを引くべきでしょうか。payoff - b * control の分散を最小にする倍率は b = cov(payoff, control) / var(control) で、まさにペイオフをコントロールに回帰したときの傾きです（CAPMのトピックのベータと同じ和）。補正後の推定値は mean(payoff) - b * (mean(control) - known_mean) です。引いている項の期待値がゼロなので、不偏です。
import numpy as np
def control_coefficient(payoff, control):
dev = control - control.mean()
return float(np.sum(dev * (payoff - payoff.mean())) / np.sum(dev ** 2))
rng = np.random.default_rng(1)
z = rng.standard_normal(2000)
s_t = 100.0 * np.exp((0.05 - 0.5 * 0.2 ** 2) + 0.2 * z)
call = np.maximum(s_t - 100.0, 0.0)
known = 100.0 * np.exp(0.05)
b = control_coefficient(call, s_t)
controlled = call.mean() - b * (s_t.mean() - known)
print(round(b, 4), round(float(call.mean()), 4), round(float(controlled), 4)) # the plain mean and the corrected one# A payoff that is exactly proportional to the control is estimated without any error
import numpy as np
x = np.random.default_rng(2).standard_normal(50) + 5.0 # the control, known mean 5.0
y = 2.0 * x # the payoff
dev = x - x.mean()
b = np.sum(dev * (y - y.mean())) / np.sum(dev ** 2)
print(round(float(b), 6), round(float(y.mean() - b * (x.mean() - 5.0)), 6)) # 2.0 and exactly 10.0どれだけ誤差が減るか
調整後の標本 payoff - b * control がノイズを担うので、補正後の推定の標準誤差はその標準偏差（ddof=1）を sqrt(n) で割ったものです。ペイオフとコントロールの相関を rho とすると、分散は 1 - rho^2 倍に減ります。相関0.9で81%、0.99で98%が除かれます。コントロールに追加のシミュレーションは要らず、数行増えるだけです。
import numpy as np
def control_variate_error(payoff, control):
dev = control - control.mean()
b = np.sum(dev * (payoff - payoff.mean())) / np.sum(dev ** 2)
return float((payoff - b * control).std(ddof=1) / np.sqrt(len(payoff)))
rng = np.random.default_rng(1)
z = rng.standard_normal(20000)
s_t = 100.0 * np.exp((0.05 - 0.5 * 0.2 ** 2) + 0.2 * z)
for strike in (80.0, 100.0, 120.0):
call = np.maximum(s_t - strike, 0.0)
plain = float(call.std(ddof=1) / np.sqrt(len(call)))
print(strike, round(plain, 4), round(control_variate_error(call, s_t), 4), round(plain / control_variate_error(call, s_t), 2))良いコントロールは相関が高く、既知である両方の条件が重要です。コントロールはペイオフと連動しなければならず（そうでないと
b は0に近く何も得られません）、正確な平均が既知でなければなりません（そうでないと補正は当て推量です）。アジアンの算術平均に適したコントロールは、価格に公式がある幾何平均や、満期価格そのものです。複数のコントロールを同時に使うこともでき、その場合の係数は重回帰から得られます。- コントロール変量は平均が既知のシミュレーション量。推定値からその観測誤差の何倍かを引く：
mean(Y) - b * (mean(X) - E[X]) - 最適な倍率は回帰の傾き
b = cov(Y, X) / var(X)。誤差はY - b * Xの標準偏差をsqrt(n)で割ったもの - 分散は
1 - rho^2倍に減る。良いコントロールはペイオフと強く相関し、平均が正確に既知で、追加のシミュレーションも要らない
次はリターンのためのscikit-learnパイプラインです。リターンデータで予測モデルを当てはめて評価します。
ラボ：採点される内容
control_coefficientで最適な倍率を選ぶ25点control_variate_estimateで推定値を補正する25点control_variate_errorで誤差を小さくする25点error_reductionで改善の大きさを数える25点