「その他」カテゴリーアーカイブ

HadoopからSparkへの移行

こんにちは。配信/インフラチームの佐々木です。今回は、adstirの集計システムをHadoopからSparkに移行したお話をしたいと思います。

以前の構成

adstirでは以前はAWS上のEMRでHadoopを使って集計していました。一つの集計のフローとしては以下のようになります。

Mapper(UNIQ集計)
↓
SQSにメッセージを送信
↓
Reducer(UNIQ集計)
↓
SQSにメッセージを送信
↓
Mapper(通常集計)
↓
SQSにメッセージを送信
↓
Reducer(通常集計)

このシステムでは以下のような問題点があり、それを解消するために移行を計画しました。

・Hadoopのパフォーマンスがあまり良くなく、集計に時間がかかり、かつ料金的なコストもかかる。

・MapperとReducerでアプリケーションがを分ける必要があり、UNIQ集計も別で必要になるため、煩雑な仕組みになる。

・かなり昔に開発したシステムの為、利用しているサービスや構成が古く、パフォーマンス・安定性等お世辞にも良いと言えない為、最新のアーキテクチャに刷新したい。

Sparkに移行すれば、UNIQ集計を分ける必要もなく、SQSでの通知も不要になります。つまり上記のフローが一つの処理で完結するわけで、システムが大幅に簡素化し、かつ時間短縮・コスト削減するのを見込み開発をしました。

移行する集計用のソフトウェアの検討

当初はAthenaに移行する想定でしたが、Athenaは分析用のサービスで、今回のようなバッチでの集計処理には向いていないことがわかり、Hadoopから移行しやすくかつパフォーマンスの良いSparkを採用することにしました。ですが手動での分析ではAthenaを使いたいとの要望があり、バッチ集計をSpark、分析用にAthenaを使い分けをするようにしました。

Sparkについて

Sparkとは一言で言えば「オンメモリで動くHadoop」の事で、AWSではEMRでのマネージドサービスを使うことが出来ます。利用できる言語はScala、Python、R、Java等ですが一般的にはScalaとPythonを使うことが多いようです。Spark自体はScalaで実装されています。

採用した言語とAPI

当初はPythonで実装したのですが、テストしてみたところHadoopとあまりパフォーマンスが変わりませんでした。原因はRDDでの処理がボトルネックになっている事でした。 Sparkで使える用意されているAPIはRDD、Dataframe、Datasetの3種類があります。APIの特徴は以下のページに詳しく書かれています。

https://yubessy.hatenablog.com/entry/2016/12/11/095915 https://databricks.com/blog/2015/02/17/introducing-dataframes-in-spark-for-large-scale-data-science.html

一見PythonでDataframeを使うのが良さそうなのですが、Dataframeは色々制限がある為、複雑なデータ構造を扱う場合はRDDの方が向いています。そこで全てScalaに書き換え、APIにRDDとDataframeを併用する実装にしました。またSparkはExecuter（分散）とDriver（非分散）で処理が分かれており、できる限りExecuter側で処理をすればパフォーマンスが向上しますので、そこを意識して実装する事が大事でした。

結果

集計時間とコストが1/3程度になり、期待する結果を得る事が出来ました。 AWSでは集計用のマネージドサービスが他にも多々あり、これからも追加されていくと思われますので、今後もキャッチアップして行きたいと思います。

今回は以上となります。

Windows Mixed RealityのadstirSDK研究開発について

MRサービス研究開発室の佐々木(竜)と申します。 MRサービス研究開発の一環としてHoloLensのアプリにおいて広告表示を行うSDKを開発しましたのでその紹介をさせていただきます。

前提知識

Windows Mixed Realityとは？

MRは、現在一般的な仮想現実技術のAR、VRを含んだ技術です。現実空間に仮想的なものを表示する部分から、周りのもの全てを仮想的に表示するところまでを技術範囲としています。

HoloLensとは？

Windows Mixed Realityに準じたAR端末で、端末内にPCが内蔵されており単体で動作します。視線方向の環境データを認識して現実物体に対してオクルージョンなどの処理をしながら表示することができます。

MixedRealityToolKit-Unityとは？

HoloLensのアプリケーションを作成するときに必要となる、HoloLensデバイスの入力などを容易に制御できるようにするスクリプト群です。これを用いることにより、Unity上で空間認識などを扱うアプリケーションを容易に開発することができます。

実際の動作動画

SDK設計

動作環境設定

開発環境 Unity 5.6.2f1
Visual Studio 2015

device
HoloLens Developer Edition
Windows10 Microsoft

Plugin
MixedRealityToolKit#Unity 5.6.2f1

Scene 今回はテストのため、広告を表示するのに最低限のオブジェクトのみを配置しています

機能要件

UnityのScene内に広告を表示します。
広告は実際にadstir経由で入札を行い、案件をとってきます。
HoloLens側の空間認識を利用し、実際に壁にポスターのように表示します。
人のタップしたジェスチャーを認識して、Edgeを起動LPを表示します。

動作フロー

①空間を認識する

MixedRealityToolKit内のSpatialMapping.prefabをScene内に配置することで、HoloLensで認識した空間データをもとに、空間をメッシュ化して空間に沿ったメッシュとメッシュコライダーを配置します。作成した空間が、アプリケーションの動作上十分な数に達したとき空間のメッシュへの変換を停止します。

void Start()
 {
      // Callbackとして、空間認識が終わったタイミングをハンドリングする
      SurfaceMeshesToPlanes.Instance.MakePlanesComplete += SurfaceMeshesToPlanes_MakePlanesComplete;
 }

void Update()
{
    if (!meshesProcessed)
    {
        // 設定した空間の認識時間が過ぎたかどうか判定
        if ((Time.time - SpatialMappingManager.Instance.StartTime) > scanTime)
        {  
            // 空間認識処理が走っているか
            if (SpatialMappingManager.Instance.IsObserverRunning())
            {
                // 空間認識処理を止める
                SpatialMappingManager.Instance.StopObserver();
            }
            // 認識したデータをもとにメッシュを作成する
            CreatePlanes();

            meshesProcessed = true;
        }
    }
}

void Start()

{

// Callbackとして、空間認識が終わったタイミングをハンドリングする

SurfaceMeshesToPlanes.Instance.MakePlanesComplete += SurfaceMeshesToPlanes_MakePlanesComplete;

}

void Update()

{

if (!meshesProcessed)

{

// 設定した空間の認識時間が過ぎたかどうか判定

if ((Time.time - SpatialMappingManager.Instance.StartTime) > scanTime)

{

// 空間認識処理が走っているか

if (SpatialMappingManager.Instance.IsObserverRunning())

{

// 空間認識処理を止める

SpatialMappingManager.Instance.StopObserver();

}

// 認識したデータをもとにメッシュを作成する

CreatePlanes();

meshesProcessed = true;

}

上記のようにコルーチンのStart処理で、マッピングの終了をハンドリングします。そして、あらかじめ設定した空間認識の時間を過ぎたときMapping処理を止めて、空間メッシュを作成します。空間認識の処理は重い処理となっているので、アプリケーションの要件に合わせて、空間スキャンの時間を設定します。

参照 : HoloToolkitのSpatialMappingを理解する

②広告データを取得する

①の動作と同時に広告リクエストを非同期で飛ばします。

OnEnable()
{
       // 広告の読み込みが終わったことを受け取るコールバックを登録
       HttpCommon.OnLoadAdDataJson += OnSetAdImage;
       // 広告リクエストする
       ShowAd()
}

 private void ShowAd()
 {
        AdstirMRAdNative ad = new AdstirMRAdNative("メディアID", 枠ID  );
        ad.AdLoadStart();
 }

OnEnable()

{

// 広告の読み込みが終わったことを受け取るコールバックを登録

HttpCommon.OnLoadAdDataJson += OnSetAdImage;

// 広告リクエストする

ShowAd()

}

private void ShowAd()

{

AdstirMRAdNative ad = new AdstirMRAdNative("メディアID", 枠ID );

ad.AdLoadStart();

}

広告をリクエストする前に、処理が終わったイベントを登録します。OnSetAdImage内では広告が読み込まれたフラグを立てる処理と、タップされたときのイベントをセットします。Unityの他のGameObjectから広告の読み込み終了を参照しやするするために、フラグとしてイベントの終了をUnity側で持っておきます。

Managed Pluginを作成するとき

UnityEditor側からDebugで動作させるときと、実際にHoloLensに配置して動作させるときでアプリケーションのrun timeが変わり、参照できるライブラリが変わってきます。今回のSDKでは、実際にUnityのSceneに配置する前提でした。なので、UnityEditor側からある程度動作を再現したかったので下記のようなプロジェクト構成をしました。

Pluginの実態はSharedの中に持っていて、UWP、UnityのフォルダはSharedファイルからのシンボリックリンクを張っています。また、コード内ではコンパイルを切り替えるマクロを使って、動作環境ごとのコードを書いています。

public sealed partial class AdstirMRAdNative
{
    private readonly String mediaId;
    private readonly int spotNo;

    public AdstirMRAdNative(String mediaId, int spotNo, Build_Types mode = Build_Types.RELEASE)
    {
        this.mediaId = mediaId;
        this.spotNo = spotNo;
        ConfigLoader.SetConfigSettings(mode);
    }

    public void AdLoadStart()
    {
#if NETFX_CORE
        UWPAdStirNativeAd Ad = new UWPAdStirNativeAd();
        Ad.getAd(mediaId, spotNo);
#else
        UnityAdStirNativeAd Ad = new UnityAdStirNativeAd();
        Ad.getAd(mediaId, spotNo);
#endif
    }
}

public sealed partial class AdstirMRAdNative

{

private readonly String mediaId;

private readonly int spotNo;

public AdstirMRAdNative(String mediaId, int spotNo, Build_Types mode = Build_Types.RELEASE)

{

this.mediaId = mediaId;

this.spotNo = spotNo;

ConfigLoader.SetConfigSettings(mode);

}

public void AdLoadStart()

{

#if NETFX_CORE

UWPAdStirNativeAd Ad = new UWPAdStirNativeAd();

Ad.getAd(mediaId, spotNo);

#else

UnityAdStirNativeAd Ad = new UnityAdStirNativeAd();

Ad.getAd(mediaId, spotNo);

#endif

}

上記のように、NETFX_CORE条件内でUWP側のアプリケーションの挙動、もう一方でUnity側の挙動を記述します。参照 : HoloLens用 Managed Pluginの作り方

③ 広告を配置する

広告のオブジェクトの読み込みが終わり、空間のデータも蓄積されたとき、広告を壁に埋め込みます。

GameObject surface = hitInfo.transform.gameObject;
SurfacePlane plane = surface.GetComponent&lt;SurfacePlane>();
// 視線が壁のオブジェクトと一致したとき
if (plane != null && plane.PlaneType == PlaneTypes.Wall)
{
    Vector3 position = Camera.main.transform.position + Camera.main.transform.forward * 1.25f;
    Quaternion rotation = Quaternion.identity;
    position = surface.transform.position + (plane.PlaneThickness * plane.SurfaceNormal);

    // 広告オブジェクトを壁表示する位置を取得する
    position = AdjustPositionWithSpatialMap(position, plane.SurfaceNormal);
    rotation = Quaternion.LookRotation(surface.transform.forward, Vector3.up);

　 // 広告オブジェクトの生成
    GameObject spaceObject = Instantiate(spaceObjectPrefabs, position, rotation) as GameObject;
    spaceObject.transform.parent = gameObject.transform;
}

GameObject surface = hitInfo.transform.gameObject;

SurfacePlane plane = surface.GetComponent<SurfacePlane>();

// 視線が壁のオブジェクトと一致したとき

if (plane != null && plane.PlaneType == PlaneTypes.Wall)

{

Vector3 position = Camera.main.transform.position + Camera.main.transform.forward * 1.25f;

Quaternion rotation = Quaternion.identity;

position = surface.transform.position + (plane.PlaneThickness * plane.SurfaceNormal);

// 広告オブジェクトを壁表示する位置を取得する

position = AdjustPositionWithSpatialMap(position, plane.SurfaceNormal);

rotation = Quaternion.LookRotation(surface.transform.forward, Vector3.up);

　 // 広告オブジェクトの生成

GameObject spaceObject = Instantiate(spaceObjectPrefabs, position, rotation) as GameObject;

spaceObject.transform.parent = gameObject.transform;

}

ここでは、広告オブジェクトを配置する壁オブジェクトを決定します。視線方向を取得して、壁のオブジェクトと視線がぶつかる位置に広告を配置します。壁のオブジェクトが広告を配置することができない幅の壁のオブジェクトが、選択された場合は、壁オブジェクトの再配置を行います。

④LPを表示する

OnEable()
{
    InteractionManager.OnTapped += ReceivedTapped;
}

private void ReceivedTapped(GameObject go, InteractionManager.InteractionEventArgs e)
{
    if (CanvasObjCollider(go))
    {
        if (launchedDialog)
            return;
        var buttons = SimpleDialog.ButtonTypeEnum.Yes | SimpleDialog.ButtonTypeEnum.No;
        var title = "alert";
        var message = "Open the link in Edge.Is it OK ?";
        SimpleDialog simpleDialog = SimpleDialog.Open(dialogObj, buttons, title, message);
        simpleDialog.OnClosed += OnClosed;
        StartCoroutine(LaunchDialog(simpleDialog));
        launchedDialog = true;
    }
}

OnEable()

{

InteractionManager.OnTapped += ReceivedTapped;

}

private void ReceivedTapped(GameObject go, InteractionManager.InteractionEventArgs e)

{

if (CanvasObjCollider(go))

{

if (launchedDialog)

return;

var buttons = SimpleDialog.ButtonTypeEnum.Yes | SimpleDialog.ButtonTypeEnum.No;

var title = "alert";

var message = "Open the link in Edge.Is it OK ?";

SimpleDialog simpleDialog = SimpleDialog.Open(dialogObj, buttons, title, message);

simpleDialog.OnClosed += OnClosed;

StartCoroutine(LaunchDialog(simpleDialog));

launchedDialog = true;

}

Sceneに配置した。IntractionManagerにTapしたイベントを登録します。IntaractionManagerは人がGestureなどをしたときにイベントをフックしてくれていて、このコンポーネントを使うことでGestureなどの処理を容易にしてくれています。 Tapをされたオブジェクトが、実際の広告のオブジェクトだった時に、ダイアログを表示して遷移するかを問い合わせます。遷移するが選ばれたとき、アプリケーションをサスペンド状態にして、Edgeを呼び出します。

参照 : MixedRealityToolkit-UnityのInputManagerを紐解く～Gesture編～

以上が一通りの動作となります。

まとめ

壁の認識、ジェスチャーからのインプットなどを利用したHoloLens用adstirMRSDKを作ってきました。現状はまだ、動作環境が制限されるなど様々な課題が多く実用的なSDKとはなっていませんが、ARKit/ARCoreなどのモバイルのARフレームワークの発展とともにコンテンツが成熟してきたとき、また、広告を利用していただける機会が生まれたとき、より良いユーザー体験ができるSDKを提供していきたいと思っています。

ユナイテッドでは、エンジニアを絶賛募集しています！このように最新技術にトライできる環境で一緒に働きませんか？

乱数の生成について

こんにちは。データ解析部の西岡(敏)と申します。広告配信の効果分析とアルゴリズム実装を担当しております。

今回は、システムで利用する乱数の生成について書きたいと思います。

乱数について

過去、コンピュータの無い時代には、サイコロを振ってその数字を記録したり、コインを投げてその裏表を調べたりすることで乱数を得ることができました。ただし、こういった乱数は大量のデータを必要とする実験では不便であるため、現在はそのような実験を行う際にはコンピュータを利用して擬似的に乱数を生成しています。

こういった形でコンピュータを用いて生成する乱数のことを、前述のサイコロの乱数などとは区別し、算術乱数や擬似乱数と呼びます。

プログラムで使われる乱数

前述のように、疑似乱数はコンピュータから生成を行うため、乱数が本来持っているべき予測不可能性を持ちません。ですが、出来る限りランダムであることが担保され、非常に長い周期を持てば、シミュレーションなどに利用するのに差し支えはないということになります。

疑似乱数は主にコンピュータ上で生成・使用され、そのために生成される乱数の質や周期の長さと同等に、コンピュータ上で生成される速度、プログラミングの簡素さ、使用する記憶領域の大きさが乱数生成器の評価のポイントとなっています。

以下では、システムで広く使われている代表的なアルゴリズムである、線形合同法と、メルセンヌ・ツイスタの説明を行っていきます。

線形合同法

概要

1948年頃にレーマーが提案した方法で、今日でも広く使われている手法になります。 BASICやPascal, UNIXなどでこちらの方法が乱数生成に使われています。

基本的な算出方法は、適当なX(0), 乗数a, 加数b, 法(modulus)Mを指定し、

X(n) = (a * X(n-1) + b) % M として、計算を行うものになります。

このとき、 bがMと互いに素 a-1がMを割り切るすべての素数の倍数 Mが4の倍数であれば

a-1も4の倍数

という条件が、X(n)の周期が最大値Mとなる必要十分条件になります。

pythonによるプログラミング

実際にpythonで線形合同法を用いて乱数の生成を行ってみます。

# 線形合同法

# Xの初期値, 定数a, b, Mを決めます。
# 定数a, b, Mに関してはNumerical Recipesという書籍で推奨されている値を入力しました。

list_x = [2497]
const_a = 1664525
const_b = 1013904223
const_m = 2 ** 32

for i in range(10):
    list_x.append((const_a * list_x[i] + const_b) % const_m)

print(list_x)

# 線形合同法

# Xの初期値, 定数a, b, Mを決めます。

# 定数a, b, Mに関してはNumerical Recipesという書籍で推奨されている値を入力しました。

list_x = [2497]

const_a = 1664525

const_b = 1013904223

const_m = 2 ** 32

for i in range(10):

list_x.append((const_a * list_x[i] + const_b) % const_m)

print(list_x)

Out[]:

[2497,
 875255852,
 3289380251,
 2302488126,
 419758213,
 2363619360,
 2816822527,
 2762503762,
 1176740233,
 3304832340,
 131795107]

Out[]:

[2497,

875255852,

3289380251,

2302488126,

419758213,

2363619360,

2816822527,

2762503762,

1176740233,

3304832340,

131795107]

線形合同法の問題点

上記のように、簡単に生成ができる一方で、この手法には周期の限界が2^32にあること、下位ビットには周期的な偏りがあること(上記の例でも、下一桁が奇数→偶数の繰り返しになっている) などの問題点があります。そのため、下位ビット部分を切り捨てて実装を行っている場合もあるそうです。

メルセンヌ・ツイスタ

概要

続いては、今日様々なシステムで利用されているメルセンヌ・ツイスタに関して説明を行います。メルセンヌ・ツイスタはM系列を松本眞氏と西村拓士氏が発展させたものであり、疑似乱数生成器として現在最も高い評価を得ています。

メルセンヌ・ツイスタの説明

メルセンヌ・ツイスタの説明に当たって、まずM系列の説明を行います。

M系列

M系列とは、Xn = Xn-p + Xn-qで表現できる1ビットの数列になります。

例えば、 X0 = 1, X1 = 0, X2 = 0, p = 3, q = 1 としたとき、

X3 = X0 + X2 = 1 + 0 = 1
X4 = X1 + X3 = 0 + 1 = 1
X5 = X2 + X4 = 0 + 1 = 1
X6 = X3 + X5 = 1 + 1 = 0

という計算によって算出できます。

TGFSR

M系列を発展させたものに、 Twisted General Feedback Shift Register(TGFSR)があり、こちらは

Xn+p = Xn+q + Xn *A (p>q>0)

を用いて表現することができます。ここで、Aはω✕ωの行列、Xはω次元の行ベクトルになります。また、このときωは32, 64など、生成する乱数のビット数となります。

メルセンヌ・ツイスタ

さらにTGFSRを

Xn+p = Xn+q + Xn+1 *B + Xn *C (p>q>0)

として改良したものがメルセンヌ・ツイスタになります。

特に、こちらのこちらの最大周期が2^19937 - 1のときにMT19937と呼ばれ、通常のメルセンヌ・ツイスタはこのMT19937を指します。

具体的な計算方法としては、 init_genrand関数とgenrand_int関数を準備し、 init_genrand関数を実行した後に、生成する乱数分genrand_int関数を実行するという流れになります。それぞれの関数は以下の形式となります。

init_genrand 関数
- 0 番地に種を代入する。
- i 番地の値を得るために、i-1 番地の値を使用する。シフト、EXOR、乗算、加算の順に演算を行う。
- 623 番地まで計算を繰り返し、値を代入する。
genrand_int 関数 a. (624n+1)回目の場合(n=0,1,2,…)
- 0 番地から 227 番地まで AND、OR、EXOR 演算を行う。
- 228 番地から 623 番地まで、別の数値を使って AND、OR、EXOR 演算を行う。
- 再び別の数値を使って、AND、OR、EXOR 演算を行い、y に代入する。 b. 毎回共通の部分
- y にシフト、AND、EXOR 演算を行い、生成された乱数を返す。

pythonによるプログラミング

実際にpythonで線形合同法を用いて乱数の生成を行ってみます。

# MT19937の記述

w = 32
n = 624
m = 397
r = 31
u = 11
s = 7
t = 15
l = 18
a = 0x9908B0DF
b = 0x9D2C5680
c = 0xEFC60000

# print(a, b, c)

i = 0
x = []

mask_u = 0
mask_l = 0

for i in range(w):
    mask_u = (mask_u &lt;&lt; 1) + (i &lt; w - r)
    mask_l = (mask_l &lt;&lt; 1) + (i >= w - r)

# print(mask_u)
# print(mask_l)

def init_genrand(seed):
    x = [seed & 0xFFFFFFFF]
    for j in range(1,n):
        x.append((1812433253 * (x[j-1] ^ (x[j-1] >> 30)) + j) & 0xFFFFFFFF)
    # print(x)
    return x


x = init_genrand(20150919)

# 乱数の生成
def genrand_int(x, i):
    # step1
    z = (x[i] & mask_u | x[(i+1) % n] & mask_l)
    # step2
    if (z & 1) == 0:
        d = 0
    else:
        d = a
    x[i] = x[(i+m) % n] ^ (z >> 1) ^ d
    # step3
    y = x[i]

    y = y ^ (y >> u)
    y = y ^ ((y &lt;&lt; s) & b)
    y = y ^ ((y &lt;&lt; t) & c)
    y = y ^ (y >> l)

    # 乱数を返す

    return y

i = 0
for j in range(1, 11):
    print(genrand_int(x, i))
    # カウンタをインクリメント
    i = (i + 1)%n

# MT19937の記述

w = 32

n = 624

m = 397

r = 31

u = 11

s = 7

t = 15

l = 18

a = 0x9908B0DF

b = 0x9D2C5680

c = 0xEFC60000

# print(a, b, c)

i = 0

x = []

mask_u = 0

mask_l = 0

for i in range(w):

mask_u = (mask_u << 1) + (i < w - r)

mask_l = (mask_l << 1) + (i >= w - r)

# print(mask_u)

# print(mask_l)

def init_genrand(seed):

x = [seed & 0xFFFFFFFF]

for j in range(1,n):

x.append((1812433253 * (x[j-1] ^ (x[j-1] >> 30)) + j) & 0xFFFFFFFF)

# print(x)

return x

x = init_genrand(20150919)

# 乱数の生成

def genrand_int(x, i):

# step1

z = (x[i] & mask_u | x[(i+1) % n] & mask_l)

# step2

if (z & 1) == 0:

d = 0

else:

d = a

x[i] = x[(i+m) % n] ^ (z >> 1) ^ d

# step3

y = x[i]

y = y ^ (y >> u)

y = y ^ ((y << s) & b)

y = y ^ ((y << t) & c)

y = y ^ (y >> l)

# 乱数を返す

return y

i = 0

for j in range(1, 11):

print(genrand_int(x, i))

# カウンタをインクリメント

i = (i + 1)%n

Out[]:

909111202

1775643220

3451004665

2605250323

792343277

400246621

276307393

1249125328

2330296989

749880315

以上のような形で、乱数の生成を行うことができました。このメルセンヌ・ツイスタに関しては、線形合同法であったような、周期が短かったり、値によっては下一桁が偶数, 奇数と交互に出現するといった問題点が発生せず、非常に優れた擬似乱数生成器として、様々なプログラミング言語に幅広く利用されています。

まとめ

このように乱数と言っても複数の種類があり、それぞれの乱数には特徴があります。現在は、高い精度の検証をする場合などは、線形合同法で作成されたよりもメルセンヌ・ツイスタの方が良いと言われています。もし何かしらの実験を行っているときなどに、こういった乱数に着目すると、意外な解決策が見いだせるかもしれません。

参考

http://www.math.sci.hiroshima-u.ac.jp/~m-mat/TEACH/ichimura-sho-koen.pdf http://www.soi.wide.ad.jp/class/20010000/slides/03/sfc2002.pdf http://random.ism.ac.jp/info01/random_number_generation/node6.html http://www.utp.or.jp/book/b300859.html http://www001.upp.so-net.ne.jp/isaku/rand.html https://narusejun.com/archives/5/