轉換密度拒絕 (TDR)#
必要條件：T-凹 PDF，dPDF
可選：模式，中心
速度
設定：慢
採樣：快
TDR 是一種接受/拒絕方法，它利用轉換密度的凹性來自動建構帽函數和擠壓函數。這類 PDF 稱為 T-凹。目前已實作以下轉換
除了 PDF 之外，它還需要 PDF 相對於 x（即變數）的導數。這些函數必須作為 Python 物件的方法存在，然後可以將其傳遞給生成器以實例化它們的物件。已實作的變體使用與帽函數成比例的擠壓函數 ([1])。
下面顯示了使用此方法的範例
>>> import numpy as np
>>> from scipy.stats.sampling import TransformedDensityRejection
>>> from scipy.stats import norm
>>>
>>> class StandardNormal:
... def pdf(self, x):
... # note that the normalization constant is not required
... return np.exp(-0.5 * x*x)
... def dpdf(self, x):
... return -x * np.exp(-0.5 * x*x)
...
>>> dist = StandardNormal()
>>>
>>> urng = np.random.default_rng()
>>> rng = TransformedDensityRejection(dist, random_state=urng)
>>> rng.rvs()
-1.526829048388144
在上面的範例中，我們使用了 TDR 方法從標準常態分佈中採樣。請注意，在計算 PDF 時，我們可以捨棄正規化常數。這通常有助於加速採樣階段。另外，請注意 PDF 不需要向量化。它應該接受並傳回一個純量。
也可以使用 ppf_hat 方法直接評估帽分佈的 CDF 的反函數。
>>> rng.ppf_hat(0.5)
-0.00018050266342362759
>>> norm.ppf(0.5)
0.0
>>> u = np.linspace(0, 1, num=10)
>>> rng.ppf_hat(u)
array([ -inf, -1.22227372, -0.7656556 , -0.43135731, -0.14002921,
0.13966423, 0.43096141, 0.76517113, 1.22185606, inf])
>>> norm.ppf(u)
array([ -inf, -1.22064035, -0.76470967, -0.4307273 , -0.1397103 ,
0.1397103 , 0.4307273 , 0.76470967, 1.22064035, inf])
除了 PPF 方法之外，還可以存取其他屬性，以查看生成器與給定分佈的擬合程度。這些是
‘squeeze_hat_ratio’：生成器的（擠壓函數下方的面積）/（帽函數下方的面積）。它是一個介於 0 和 1 之間的數字。越接近 1 表示帽函數和擠壓函數緊密地包絡分佈，並且產生樣本所需的 PDF 評估次數越少。密度的預期評估次數以每個樣本
(1/squeeze_hat_ratio) - 1為界。預設情況下，它保持在 0.99 以上，但可以透過傳遞max_squeeze_hat_ratio參數來更改。‘hat_area’：生成器的帽函數下方的面積。
‘squeeze_area’：生成器的擠壓函數下方的面積。
>>> rng.squeeze_hat_ratio 0.9947024204884917 >>> rng.hat_area 2.510253139791547 >>> rng.squeeze_area 2.4969548741894876 >>> rng.squeeze_hat_ratio == rng.squeeze_area / rng.hat_area True
可以透過傳遞 domain 參數來截斷分佈
>>> urng = np.random.default_rng()
>>> rng = TransformedDensityRejection(dist, domain=[0, 1], random_state=urng)
>>> rng.rvs(10)
array([0.05452512, 0.97251362, 0.49955877, 0.82789729, 0.33048885,
0.55558548, 0.23168323, 0.13423275, 0.73176575, 0.35739799])
如果未指定 domain，則使用 dist 物件的 support 方法來確定 domain
>>> class StandardNormal:
... def pdf(self, x):
... return np.exp(-0.5 * x*x)
... def dpdf(self, x):
... return -x * np.exp(-0.5 * x*x)
... def support(self):
... return -np.inf, np.inf
...
>>> dist = StandardNormal()
>>>
>>> urng = np.random.default_rng()
>>> rng = TransformedDensityRejection(dist, random_state=urng)
>>> rng.rvs(10)
array([-1.52682905, 2.06206883, 0.15205036, 1.11587367, -0.30775562,
0.29879802, -0.61858268, -1.01049115, 0.78853694, -0.23060766])
如果 dist 物件未提供 support 方法，則 domain 假定為 (-np.inf, np.inf)。
若要增加 squeeze_hat_ratio，請傳遞 max_squeeze_hat_ratio
>>> dist = StandardNormal()
>>> rng = TransformedDensityRejection(dist, max_squeeze_hat_ratio=0.999,
... random_state=urng)
>>> rng.squeeze_hat_ratio
0.999364900465214
讓我們看看這如何影響對分佈的 PDF 方法的回呼
>>> from copy import copy
>>> class StandardNormal:
... def __init__(self):
... self.callbacks = 0
... def pdf(self, x):
... self.callbacks += 1
... return np.exp(-0.5 * x*x)
... def dpdf(self, x):
... return -x * np.exp(-0.5 * x*x)
...
>>> dist1 = StandardNormal()
>>> urng1 = np.random.default_rng()
>>> urng2 = copy(urng1)
>>> rng1 = TransformedDensityRejection(dist1, random_state=urng1)
>>> dist1.callbacks # evaluations during setup
139
>>> dist1.callbacks = 0 # don't consider evaluations during setup
>>> rvs = rng1.rvs(100000)
>>> dist1.callbacks # evaluations during sampling
527
>>> dist2 = StandardNormal()
>>> # use the same stream of uniform random numbers
>>> rng2 = TransformedDensityRejection(dist2, max_squeeze_hat_ratio=0.999,
... random_state=urng2)
>>> dist2.callbacks # evaluations during setup
467
>>> dist2.callbacks = 0 # don't consider evaluations during setup
>>> rvs = rng2.rvs(100000)
>>> dist2.callbacks # evaluations during sampling
84 # may vary
如我們所見，當我們增加 squeeze_hat_ratio 時，採樣期間所需的 PDF 評估次數要少得多。PPF-hat 函數也更準確
>>> abs(norm.ppf(0.975) - rng1.ppf_hat(0.975))
0.0027054565421578136
>>> abs(norm.ppf(0.975) - rng2.ppf_hat(0.975))
0.00047824084476300044
但是，請注意，這是以增加設定期間的 PDF 評估次數為代價的。
對於模式不接近 0 的密度，建議透過傳遞 mode 或 center 參數來設定分佈的模式或中心。後者是模式或分佈平均值的近似位置。此位置提供有關 PDF 主要部分的一些資訊，並用於避免數值問題。
>>> # mode = 0 for our distribution
>>> # if exact mode is not available, pass 'center' parameter instead
>>> rng = TransformedDensityRejection(dist, mode=0.)
預設情況下，此方法使用 30 個建構點來建構帽函數。可以透過傳遞 construction_points 參數來更改此設定，該參數可以是建構點陣列，也可以是表示要使用的建構點數量的整數。
>>> rng = TransformedDensityRejection(dist,
... construction_points=[-5, 0, 5])
此方法接受許多其他設定參數。請參閱文件以取得完整清單。有關參數和方法的更多資訊，請參閱 UNU.RAN 使用者手冊的第 5.3.16 節。