QuantLib 金融計算——數學工具之隨機數發生器
阿新 • • 發佈：2018-12-16
目錄
如果未做特別說明，文中的程式都是 Python3 程式碼。
QuantLib 金融計算——數學工具之隨機數發生器
載入模組
import QuantLib as ql
import scipy
print(ql.__version__)1.12概述
隨機模擬通常從產生均勻分佈的隨機數開始。假設 \(X \sim U [0, 1]\)
均勻分佈的隨機數發生器主要分兩種：
偽隨機數
quantlib-python 提供了以下三種均勻分佈的（偽）隨機數發生器：
KnuthUniformRng，高德納（Knuth）演算法LecuyerUniformRng，L'Ecuyer 演算法MersenneTwisterUniformRng，著名的梅森旋轉（Mersenne-Twister）演算法
隨機數發生器的建構函式，
Rng(seed)其中
seed，整數，預設值是 0，作為種子用於初始化相應的確定性序列；
隨機數發生器的成員函式：
next()：返回一個SampleNumber物件，作為模擬的結果。
r = rng.next()
v = r.value(r)使用者通過反覆呼叫成員函式 next() 獲得一連串的隨機數，需要注意的是 r 的型別是 SampleNumber，需要呼叫 value()
例子 1，
def testingRandomNumbers1():
seed = 1
unifMt = ql.MersenneTwisterUniformRng(seed)
unifLec = ql.LecuyerUniformRng(seed)
unifKnuth = ql.KnuthUniformRng(seed)
print('{0:<25}{1:<25}{2:<25}'.format(
'Mersenne Twister', 'Lecuyer', 'Knut'))
for i in range(10):
print('{0:<25}{1:<25}{2:<25}'.format(
unifMt.next().value(),
unifLec.next().value(),
unifKnuth.next().value()))
testingRandomNumbers1()Mersenne Twister Lecuyer Knut
0.41702199855353683 0.2853808990946861 0.4788952510312594
0.9971848082495853 0.2533581892659171 0.7694635535665499
0.7203244894044474 0.09346853100919404 0.47721285286866455
0.9325573613168672 0.6084968907396475 0.15752737762851
0.00011438119690865278 0.90342026007861 0.6065713927733087 正態分佈（偽）隨機數
隨機模擬中最常見的分佈是正態分佈，quantlib-python 提供的正態分佈隨機數發生器有 4 類：
CentralLimitABCGaussianRngBoxMullerABCGaussianRngMoroInvCumulativeABCGaussianRngInvCumulativeABCGaussianRng
其中 ABC 特指一種均勻隨機數發生器。
具體來講 4 類發生器分為 12 種：
CentralLimitLecuyerGaussianRngCentralLimitKnuthGaussianRngCentralLimitMersenneTwisterGaussianRngBoxMullerLecuyerGaussianRngBoxMullerKnuthGaussianRngBoxMullerMersenneTwisterGaussianRngMoroInvCumulativeLecuyerGaussianRngMoroInvCumulativeKnuthGaussianRngMoroInvCumulativeMersenneTwisterGaussianRngInvCumulativeLecuyerGaussianRngInvCumulativeKnuthGaussianRngInvCumulativeMersenneTwisterGaussianRng
隨機數發生器的建構函式：
rng = Rng(seed)
grng = Gaussianrng(rng)正態分佈隨機數發生器接受一個對應的均勻分佈隨機數發生器作為源，以 BoxMullerMersenneTwisterGaussianRng 為例，需要配置一個 MersenneTwisterUniformRng 物件作為隨機數的源，使用經典的 Box-Muller 演算法得到正態分佈隨機數。
例子 2，
def testingRandomNumbers2():
seed = 12324
unifMt = ql.MersenneTwisterUniformRng(seed)
bmGauss = ql.BoxMullerMersenneTwisterGaussianRng(unifMt)
for i in range(5):
print(bmGauss.next().value())
testingRandomNumbers2()-1.1756781173398896
0.14110041851886157
1.569582906805544
-0.026736779238941934
-0.8220676600472409擬隨機數
相較於之前描述的“偽”隨機數，隨機模擬中另一類重要的隨機數成為“擬”隨機數，也稱為低偏差序列。因為收斂性更好，擬隨機數通常用於高維隨機變數的模擬。quantlib-python 提供的擬隨機數有兩類，
HaltonRsg: Halton 序列SobolRsg: Sobol 序列
HaltonRsg
HaltonRsg 的建構函式，
HaltonRsg(dimensionality,
seed,
randomStart,
randomShift)其中，
dimensionality：整數，設定維度；seed，整數，預設值是 0，作為種子用於初始化相應的確定性序列；randomStart：布林值，預設是True，是否隨機開始；randomShift：布林值，預設是False，是否隨機平移。
HaltonRsg 的成員函式，
nextSequence()：返回一個SampleRealVector物件，作為模擬的結果；lastSequence()：返回一個SampleRealVector物件，作為上一個模擬的結果；dimension()：返回維度。
SobolRsg
SobolRsg 的建構函式，
SobolRsg(dimensionality,
seed,
directionIntegers=Jaeckel)其中，
dimensionality：整數，設定維度；seed，整數，預設值是 0，作為種子用於初始化相應的確定性序列；directionIntegers，quantlib-python 的內建變數，預設值是SobolRsg.Jaeckel，用於 Sobol 序列的初始化。
SobolRsg 的成員函式，
nextSequence()：返回一個SampleRealVector物件，作為模擬的結果；lastSequence()：返回一個SampleRealVector物件，作為上一個模擬的結果；dimension()：返回維度。skipTo(n)：n是整數，跳轉到抽樣結果的第 n 個維度；nextInt32Sequence()：返回一個IntVector物件。
例子 3，
def testingRandomNumbers4():
dim = 5
haltonGen = ql.HaltonRsg(dim)
sobolGen = ql.SobolRsg(dim)
sampleHalton = haltonGen.nextSequence().value()
sampleSobol = sobolGen.nextSequence().value()
print('{0:<25}{1:<25}'.format(
'Halton', 'Sobol'))
for i in range(dim):
print('{0:<25}{1:<25}'.format(
sampleHalton[i],
sampleSobol[i]))
testingRandomNumbers4()Halton Sobol
0.04081786540336907 0.5
0.8535710143553551 0.5
0.69400573329408 0.5
0.818105927979147 0.5
0.878826694887864 0.5 兩類隨機數的收斂性比較
最後用一個例子比較兩類隨機數的收斂性，分別產生正態分佈的偽隨機數和擬隨機數，計算分佈的四個統計指標：
- 均值（理論值等於 0.0）；
- 方差（理論值等於 1.0）；
- 偏度（理論值等於 0.0）；
- 超額峰度（理論值等於 0.0）。
def testingRandomNumbers5():
sobolGen = ql.SobolRsg(1)
seed = 12324
unifMt = ql.MersenneTwisterUniformRng(seed)
bmGauss = ql.BoxMullerMersenneTwisterGaussianRng(unifMt)
boxMullerStat = ql.IncrementalStatistics()
sobolStat = ql.IncrementalStatistics()
invGauss = ql.MoroInverseCumulativeNormal()
numSim = 10000
for j in range(numSim):
boxMullerStat.add(bmGauss.next().value())
currSobolNum = sobolGen.nextSequence().value()[0]
sobolStat.add(invGauss(currSobolNum))
stats = {
"BoxMuller Mean:": boxMullerStat.mean(),
"Sobol Mean:": sobolStat.mean(),
"BoxMuller Var:": boxMullerStat.variance(),
"Sobol Var:": sobolStat.variance(),
"BoxMuller Skew:": boxMullerStat.skewness(),
"Sobol Skew:": sobolStat.skewness(),
"BoxMuller Kurtosis:": boxMullerStat.kurtosis(),
"Sobol Kurtosis:": sobolStat.kurtosis()}
for k, v in stats.items():
print('{0:>25}{1:>30}'.format(k, v))
testingRandomNumbers5() BoxMuller Mean: 0.005966482725988245
Sobol Mean: -0.0002364019095203635
BoxMuller Var: 1.0166044844467006
Sobol Var: 0.9986010126883317
BoxMuller Skew: 0.02100635339070779
Sobol Skew: -7.740573185322994e-05
BoxMuller Kurtosis: -0.0340476839897507
Sobol Kurtosis: -0.020768126049145776直觀上看 Sobol 序列的結果更加接近理論值，這證明使用擬隨機數做模擬的收斂速度更好。