この規格ページの目次
7
Z 9031 : 2012
がある。そのための代表的な方法を5.4.2に規定する。
5.4.2 3項GFSR法
3項GFSR(Generalized Feedback Shift Register)法は,次の漸化式によってwビット二進整数を生成し,
疑似乱数として用いるものである。
Xn p Xn q Xn (n ,1 ,3,2)
ここに,p>qは整数定数,Xnはwビット二進整数である。この生成法のパラメータは(p, q, w)であり,
シードはX1,···,Xpである。多項式tp+tq+1が特性多項式となる。
注記1 この生成法は,非常に高速で,周期も長くできる。表2の特性多項式を使用すれば,周期は
289−1から29 689−1までの間となる。
注記2 シードとしてX1,X2,···,Xpのp個の値を与える必要があるが,生成される数列の乱数性は
これらのシードに強く依存するので,ユーザが勝手に与えるのは危険である。附属書Bに与
えられているような初期化ルーチンを使用することが望ましい。この手法は,自己相関関数
及びk次元均等分布性(附属書JA参照)を考慮して設計されている。
注記3 pを超えた個数の連続する疑似乱数の組を見ると,各ビットに現れる0の数及び1の数は対
称的に分布しない。特に,ランダムウォーク,パーコレーションなど,過去の長い(p以上
の)履歴に依存する高次元分布が問題となるシミュレーションに使うのは危険である。
――――― [JIS Z 9031 pdf 11] ―――――
8
Z 9031 : 2012
表2−3項原始特性多項式
p q
89 38
127 1,7,15,30,63
521 32,48,158,168
607 105,147,273
1279 216,418
2281 715,915,1029
3217 67,576
4423 271,369,370,649,1393,1419,2098
9689 84,471,1836,2444,4187
注記 qは,特性多項式の非零項のべき(冪)を表す。
5.4.3 5項GFSR法
これは,5項の特性多項式を使用するものであり,漸化式
Xn p Xn q1 Xn q2 Xn q3 Xn (n ,1 ,2 ,3)
でwビットの二進整数の列を生成するものである。3項GFSR法に比べ,1,0の偏りは小さくなっている。
最大周期は2p−1である。(p, q1, q2, q3, w) がパラメータであり,X1,···,Xpがシードである。最大周期を
与えるp,q1,q2,q3の組合せの例を表3に示す。
表3−5項原始特性多項式
P q1 q2 q3
89 20 40 69
107 31 57 82
127 22 63 83
521 86 197 447
607 167 307 461
1279 339 630 988
2203 585 1197 1656
2281 577 1109 1709
3217 809 1621 2381
4253 1093 2254 3297
4423 1171 2273 3299
9689 2799 5463 7712
注記 q1, q2, q3は,非零項のべき(冪)を表す。
5.4.4 結合トーズワース(Combined Tausworthe)法
0, 1数列x0,x1,x2,x3,···を漸化式
xn p xn q xn (mod )2 (n ,1,0 ,3,2)
によって定義される3項M系列とする。この0, 1数列から連続するw個をつなげたものをw桁の二進整
数とみなす。最初にw桁の二進整数x0x1···xw−1が得られ,次にt個ずらしてw桁の二進整数xtxt+1···xt+p−1
を得る。これを繰り返すとwビット整数の列
Xn xnt w 1 (n
xntxnt 1 ,1,0 ,3,2)
が得られる。ただし,tは2p−1と互いに素な自然数,wはpを超えない自然数とする。この系列を (p, q, t)
をパラメータとするトーズワース系列という。シードは,x0,···,xp−1である。シードは,どれかは1で
――――― [JIS Z 9031 pdf 12] ―――――
9
Z 9031 : 2012
あるように選ぶ。この系列の周期も2p−1となる。
注記1 二つの自然数が互いに素とは,それらの最大公約数が1であることである。
例 特性多項式として,原始3項式t4+t+1を選んだとする。上記の漸化式でいえば,p=4,q=1と選
んだことになる。もしシードとして(x0,x1,x2,x3)=(1,1,1,1)が漸化式に与えられたとすると,生成さ
れるM系列は1,1,1,1, 0,0,0,1, 0,0,1,1, 1,1,1,0, ···であり,その周期は24−1=15である。例えば,15
と互いに素であるt=4を選び,w=4とすると,トーズワース法のパラメータは(4,1,4)で表される。
生成される系列[{Xn}]は,次のとおりとなる。
X0=x0x1x2x3=1111(=15)
X1=x4x5x6x7=0001(=1)
X2=x8x9x10x11=0011(=3)
X3=x12x13x14x0=0101(=5)
X4=x1x2x3x4=1110(=14)
X5=x5x6x7x8=0010(=2)
···
このトーズワース法で生成された系列は,10進記法で15, 1, 3, 5, 14, 2, 6, 11, 12, 4, 13, 7, 8, 9, 10, 15, 1, 3,
···となり,周期は24−1=15である。
結合トーズワース法は上記のトーズワース系列を複数生成し,排他的論理和で合成して疑似乱数を生成
する手法である。例えばJ個のトーズワース系列[{Xn(j)}], j=1,2,···,Jを用意したとする。これらの系列は全
て同じワード長wをもつとする。結合トーズワース法とは,これらの系列の2進表現のビットごとの排他
的論理和
Xn=Xn(1) Xn(2) ··· Xn(J) (n=0,1,2,···)
によって,一つの系列[{Xn}]を生成する方法である。
結合トーズワース法のパラメータ,シードはそれぞれのトーズワース法のパラメータ,シードを合わせ
たものである。
注記2 この生成法は,高速で多次元均等分布性(附属書JA参照)のよいものが作れる。附属書B
に添付されたコードの周期は (231−1) (229−1) (228−1) (ほぼ288)であり,メモリも3ワー
ドしか消費しない。また,他に多くの組合せが参考文献の[10]に提案されている。
5.4.5 メルセンヌツイスター(Mersenne Twister)法
Xnをwビット二進整数とする。メルセンヌツイスター法は,整数定数p,q,r及びwビット二進整数定
数aを用い,次の漸化式によってwビット二進整数列Xnを生成し,
u l
Xn p Xn q (Xn|Xn 1 ) A(n ,1 ,2 ,3
) (1)
u l
それを後述の変換によってwビット二進整数疑似乱数列とするものである。ここで, (XnX
| n 1 ) はXnの上
位w−rビット及びXn+1の下位rビットをこの順に連接して得られる二進整数を表す。Aは,aによって決
まる0及び1からなるw×w行列であり,
shiftright (X) (X の最下位ビットが 0 の場合)
XA
shiftright ( X) a(X の最下位ビットが 1 の場合)
によって積が与えられるものである(ただし,ここではXをw次元の0−1横ベクトルとみなしている)。
――――― [JIS Z 9031 pdf 13] ―――――
10
Z 9031 : 2012
ここで,shiftrightは1ビットの右シフトを表す。
上位ビットの乱数性を改善するため,Xnに次の変換を施す。
y←Xn
y←y (y>>u)
y←y ((y<<s) &b)
y←y ((y<<t) &c)
y←y (y>>l)
最後に得られたyを疑似乱数として出力する。ここに,(y>>u) はyの右uビットシフト,(y<<u) は左u
ビットシフトを表す。b,cは上位ビットの乱数性を改善するための定数ビットマスクである。
このアルゴリズムのパラメータは,(p, q, r, w, a, u, s, t, l, b, c) である。シードはX1の上位w−rビット及
び
X2,···,Xpである。シードとしては,全てが0であるものは選ばない。
注記1 消費されるメモリは,pワードで周期は2pw−r−1となり,前述のGFSR法より効率がよい。
注記2 この生成法は,高速で非常に長い周期と優れた多次元均等分布性(附属書JA参照)をもっ
ている。
6 各種の分布に従う乱数の生成
6.1 一般
整数形の一様乱数Xを使っていろいろな分布に従う乱数Yを生成する方法を示す。目標とする分布の分
布関数をF (y) で表し,それが連続分布の場合は確率密度関数をf (y) で,また,離散分布の場合は確率質
量関数(Y=yとなる確率)をp (y) で表すことにする。
6.2 標準一様分布
6.2.1 確率密度関数
,1 0 ≦ y≦ 1
f ( y)
,0 その他
6.2.2 乱数生成法
箇条5に定める方法で生成される整数形一様乱数Xの取る値の最大値をm−1とすると,
Y X /m
とする。
注記1 Xの取る値が離散的なため,Yの取る値も離散的となり,本当の連続分布とはならないこと
に注意する。
注記2 Yの値が1.0となることは決してない。Yの値が0.0となるのは,X=0の場合に限られる。こ
れが起こり得るのは,線形合同法で表1のパラメータを使う場合は,1行目を使うときだけ
である。M系列乱数の場合には,どの生成法でも起こり得る。
注記3 標準一様分布に従う乱数を標準一様乱数と呼び,U, U1, U2などで表す。このように表した場
合,これらは互いに独立であるものとする。
注記4 区間 [a, a+b]上の一様分布に従う乱数Yは,標準一様乱数Uを次の式で変換して生成する。
Y=bU+a
――――― [JIS Z 9031 pdf 14] ―――――
11
Z 9031 : 2012
6.3 標準ベータ分布
6.3.1 確率密度関数
d 1
yc 1 1y
, 0 ≦ y≦ 1
fy d
B c,
0,その他
1
ここで, B c, d xc 11 x
d 1
dxはベータ関数で,パラメータc及びdは正である。
0
6.3.2 ジェーンク(Jhnk)の方法
6.2に定める方法で,互いに独立な標準一様乱数U1及びU2を生成し,次の手順で標準ベータ乱数Yを生
成する。
~
1° Y U11 cU21 d
とする。
c~
2°もし不等式 Y~≦1が成り立てば,
Y U11 Y を採用して終了。
3°上記の不等式が成り立たなければ,U1及びU2を生成し直して1°に戻る。
6.3.3 チェン(Cheng)の方法
6.2に定める方法で,互いに独立な0でない標準一様乱数U1及びU2を生成し,次の手順で標準ベータ乱
数Yを生成する。
[準備]
0°次の式でqを定義する。
1
min c, d , min c, d
q 2cd c d
, その他
c d 2
[生成]
1°
1 U1
V , W cexp V
q1 U1
とする。
2°もし不等式
c d
c d ln c qV ln 4 ≧ ln U12U2
d W
が成り立てば,
Y=W/(d+W)を採用して終了。
3°上記の不等式が成り立たなければ,0でない標準一様乱数U1及びU2を生成し,1°に戻る。
注記1 max(c,d)≦1の場合はJhnkの方法,そうでない場合はチェンの方法を使うことが望ましい。
注記2 [a,a+b]上の一般のベータ分布に従う乱数Y' は,標準ベータ乱数Yを次の式で変換して生成
できる。
Y'=b Y+a
――――― [JIS Z 9031 pdf 15] ―――――
次のページ PDF 16
JIS Z 9031:2012の引用国際規格 ISO 一覧
- ISO 28640:2010(MOD)
JIS Z 9031:2012の国際規格 ICS 分類一覧
- 03 : サービス.経営組織,管理及び品質.行政.運輸.社会学. > 03.120 : 品質 > 03.120.30 : 統計的方法の応用
JIS Z 9031:2012の関連規格と引用規格一覧
- 規格番号
- 規格名称
- JISX0001:1994
- 情報処理用語―基本用語
- JISZ8101-1:2015
- 統計―用語及び記号―第1部:一般統計用語及び確率で用いられる用語
- JISZ8101-2:2015
- 統計―用語及び記号―第2部:統計の応用
- JISZ8101-3:1999
- 統計―用語と記号―第3部:実験計画法
- JISZ8201:1981
- 数学記号