VBAの乱数を観察する

VBAの乱数は粗悪？

先日投稿したVBA記事にて、VBAの乱数の品質が悪いと書いたところ、これは使い方が悪いのだとご指摘をいただきました。その通りでした ...

よい機会ですので、VBAの乱数を観察してみましょう。あんまり数学が得意ではないので、細かいところはご愛嬌

VBAの乱数について

VBAのRnd関数は「線形合同法（LCG）」を用いた疑似乱数を生成します。乱数値を生成する種となる内部値をシードと呼び、シード列は漸化式で定義されます。VB5の場合は、

r_{n+1} = (r_n * 0xFD43FD + 0xC39EC3) mod 0x1000000

で定義され、VB6、VBAも恐らくこれに準じます。

この場合、0以上0x1000000未満の24bitの数を取得することができます。これは、戻り値であるSingle型のサイズです。

乱数列は循環します。すでに出た値が出るまでの生成回数を周期と呼びます。線形合同法では各パラメーターを、可能な限り周期が法（modの後ろの数）に近付くように設定します。

VBAのRandomize関数は、引数を与えずに実行した場合には、以前のRandomizeによって設定された値とTimer関数の戻り値から、次のRnd関数のシードを決定するために使用します。シードがグローバル変数としてどこかに置かれており、それを変更しているようなイメージです（？）

件の問題があるコードは、このRandomizeを高速に何度も呼んだことが不適切で、想定よりも周期が短くなっていたようです。

周期の測定

かんたんに思いつくところで、同じ値が出るまで乱数値を生成し、その回数で周期を測定してみましょう。今度こそ適切に書いてみます。

Option Explicit

Const NOT_FOUND As Long = -1

Sub Main()
    
    ' できるだけ速いVBAを書く一般的なコツとして
    ' ReDimは遅いので配列のサイズは宣言時に決めたほうがよいです。
    Dim arr(0 To 100000000) As Single
    Dim i As Long
    i = 0
    
    Dim startTime As Single
    startTime = Timer
    
    ' Randomizeはここで一回
    Randomize
    Do While True
        
        Dim rand As Single
        rand = Rnd
        
        ' IndexOfが遅いので、判定の頻度を下げます。
        If i Mod 1000000 = 0 Then
            Debug.Print i & vbTab & Timer - startTime & "s"
            
            Dim duplicate As Long
            duplicate = Array_IndexOf(arr, rand, 0, i)
            
            If duplicate <> NOT_FOUND Then
                Debug.Print "0x" & Hex(i - duplicate)
                Exit Do
            End If
        End If
        
        arr(i) = rand
        i = i + 1
    Loop
End Sub

' 配列なめて見てるのでめちゃ遅いです。
Function Array_IndexOf(ByRef arr() As Single, ByVal value As Single, ByVal indexFrom As Long, ByVal indexTo As Long) As Long
    Dim i As Long
    For i = indexFrom To indexTo
        If arr(i) = value Then
            Array_IndexOf = i
            Exit Function
        End If
    Next i
    Array_IndexOf = NOT_FOUND
End Function

出力: 0x1000000

なんとVBAの乱数の周期は、理論値である0x1000000となるように調整されていました！

前述の問題のあるコードの再現として、Randomizeをループの中に入れてみましょう。

' 前略
    Do While True
        Randomize
        Dim rand As Single
' 後略

0x500
0x7900
0x500
0x7F00
0x4D00
0x6800
0x4400
0x5100

周期が非常に短くなっていることが確認できます。

一様性の確認

周期の1/10回乱数値を生成した時、0~0.1, 0.1~0.2, ... 0.9~1の範囲の数値をそれぞれ何回生成したか調べて、偏りを見てみます。断じてカイ2乗検定をぱっと書けなかったわけではありません、断じて ...

Sub Main()

    Dim arr() As Single
    arr = Generate(&H100000 - 1)
    
    Dim i As Long
    For i = 0 To 9
    
        Dim min As Single
        Dim max As Single
        min = Round(i / 10, 2)
        max = Round(min + 0.1, 2)
        
        Dim count As Long
        count = UBound(Array_Filter(arr, min, max))
        
        Debug.Print Left(CStr(min) + "  ", 3) & "~" & max & vbTab & count & vbTab & Round(count / UBound(arr) * 100, 2) & "%"
    Next i
    
End Sub

Function Generate(ByVal size As Long) As Single()
    Dim arr() As Single
    ReDim arr(0 To size)
    Dim i As Long
    
    Randomize
    For i = LBound(arr) To UBound(arr)
        arr(i) = Rnd
    Next i
    Generate = arr
End Function

' min以上max未満の値のみを含む新しい配列を返します。
Function Array_Filter(ByRef arr() As Single, ByVal min As Single, ByVal max As Single) As Single()
    Dim ret() As Single
    ReDim ret(0 To &H100000 - 1)
    Dim curIndex As Long
    curIndex = 0

    Dim i As Long
    For i = LBound(arr) To UBound(arr)
        If min <= arr(i) And arr(i) < max Then
            ret(curIndex) = arr(i)
            curIndex = curIndex + 1
        End If
    Next i
    
    ReDim Preserve ret(0 To curIndex - 1)
    Array_Filter = ret
End Function

出力:
0  ~0.1 104739  9.99%
0.1~0.2 104624  9.98%
0.2~0.3 104568  9.97%
0.3~0.4 105094  10.02%
0.4~0.5 105060  10.02%
0.5~0.6 104945  10.01%
0.6~0.7 105447  10.06%
0.7~0.8 104644  9.98%
0.8~0.9 104663  9.98%
0.9~1   104782  9.99%

すべての階級が、ほとんど偏りなく出ていることがわかります。

おわりに（懺悔）

実際に観察した結果、VBAの乱数は周期および一様性の観点で、一般的な用途には十分なものでした。下位ビットを使用しないようにするなど、線形合同法の弱点に気を付けて、安心して使用することができますね。

デマ書いてすいませんでした ...