スタートページ

数学ライブラリ(mathlib.js)

私用を目的にしたJavascriptのライブラリですので、保証はしません。誤差を考慮していないので厳密解、実務的な大規模な問題には使えません。
 (トレース版)は、作成時のデバッグに使ったものですが、処理経過がわかり、初心者向けの教育用に使えると思います。
 mathlib.jsの本体はこちらにあります。コピーしてお使いになって結構です。
 外部jsとして使うときは、次のように指定してください。
   <script src="http://www.kogures.com/hitoshi/javascript/mathlib.js"></script>



分野・機能一般形備考
●代数系
方程式の根(はさみうち法) x = xFx(fx, xmin, xmax[, eps])
方程式の根(ニュートン法) x = xFxNewton(fx, x0[, eps])
 (トレース版) x = xFxNewton(表示場所, fx, x0[, eps])
3次以下の方程式の複素根 ary = xFx3(a0, a1,..., a3) 一つの実根をニュートン法で得て2次方程式にする
 (トレース版) ary = xFx3(表示場所, a0, a1,..., a3)
数式 f(x) を最小にするx x = xminFx(fx, xmin, xmax [, eps])
数式 f(x) を最大にするx x = xmaxFx(fx, xmin, xmax [, eps])
微分(はさみうち法) a = dFx(fx, x [, d, eps])
定積分(シンプソン則) x = sFx(fx, xmin, xmax[, n2, eps]) n2はキザミ数
●多点系
補間(配列入力) polArray(n, m, xa, ya) n:次数、m:補間するxの個数、xa, yaは配列
補間( (x,y)入力) y = polPoints(x, x0,y0, x1,y1, x2,y2, ..., xn,yn) x:補間点、(xi,yi):曲線が通る点
最小二乗法(高次式:y=a0+a1x +...+anxn LeastSquare(n, p, f, xa, ya, aa) n:次数、p:観測点個数、f:内外挿点個数、aa:出力配列,係数ai
 (トレース版) LeastSquareTrace(表示場所, n, p, f, xa, ya, aa) トレース表示場所:HTMLで<div id="表示場所"></div>で指定
●配列・行列系
配列のソート(行) sortMatrix(a, keycol [, order [, fixrow]])
連立一次方程式 GaussJordan(n, a, b, x) 配列xに解が入る
連立一次方程式(特殊入力形式) GaussJordanArray(a) 全データを一つの配列 a で与える。a は維持されない。解は a[i][0] に入る
 (トレース版) GaussJordanArrayTrace(表示場所, a) 経過がHTMLの%lt;div id="〇〇"></div>に表示される
逆行列 invMatrix(n, a [, x]) xに逆行列が入る。省略するとaに入る
 (トレース版) invMatrixTrace(表示場所, n, a)
線形計画法(LP) ary = LP(m, n, a) m:制約式数、n:変数個数、a:行列配列(a[0][j]:目的関数、a[i][0]:定数項
 (トレース版) ary = LP(表示場所, m, n, a)
輸送問題 ary = TLP(m, n, a) m:供給地数、n:需要地数, a:輸送費用(a[i][0]:供給量,a[0][j]:需要量)
 (トレース版) ary = TLP(表示場所, m, n, a)
●確率分布系
基本統計量:一変数 ary = basicStat(x1, x2, ..., xn) aryは連想配列:ary["合計"], ary["標準偏差"] など
基本統計量:一変数(配列入力) ary = basicStatArray(x[, n]) x[0]~x[n-1]が対象のときnを指定、省略時はx.length
基本統計量:2変数 ary = basicStat2(x1,y1, x2,y2, ..., xn,yn) aryは連想配列:ary["共分散"], ary["相関係数"] など
基本統計量:2変数(配列入力) ary = basicStatArray(x, y[, n])
正規分布:確率密度 p = pNormalDist(x [, μ, σ])
正規分布:累積確率関数(標準形) s = sNormalDistStd(x) 正規分布(0,1)の範囲-∞~xの累積値
正規分布:累積確率関数(一般形) s = sNormalDist(x [, μ, σ])
正規分布:区間内累積確率 s = sNormalDistMinmax(xmin, xmax [, μ, σ])
正規分布:累積確率密度の逆関数(標準形) x = invNormalDistStd(p) 正規分布(0,1)の累積確率がpとなる点x
正規分布:累積確率密度の逆関数(一般形) x = invNormalDist(p [, μ, σ])
二項分布・二項係数(組合せの数)nCx m = nCx(n, x) nCx = n!/x!(n-x)!
二項分布・確率密度 p = pBinDist(x, n, p) 発生確率pの事象がn回の試行でx回発生する確率
二項分布・累積確率 s = sBinDist(x, n, p)
負の二項分布・確率密度 p = pBinNegDist(x, n, p) 発生確率pの事象がはじめてx回発生したときの試行回数nの確率分布
負の二項分布・累積確率 s = sBinNegDist(x, n, p) n回以内である確率
幾何分布・確率密度 p = pGeoDist(n, p) 発生確率pの事象がはじめて発生したときの試行回数nの確率分布
幾何分布・累積確率 s = sGeoDist(n, p)
超幾何分布・確率密度 p = pHyperGeoDist(a, b, α, β) Aがa個、Bがb個から、α+β個とりだしたときAがα個である確率
超幾何分布・累積確率 s = sHyperGeoDist(a, b, α, β)
指数分布・確率密度 p = pExpDist(t, λ) 平均λの指数分布で、t=t のときの確率密度
指数分布・累積確率 s = sExpDist(t, λ) 区間 t=0~t の累積確率
指数分布・区間内累積確率 s = sExpDistMinmax(tmin, tmax, λ) 区間 t=tmin~tmax の累積確率
指数分布・累積確率の逆関数 t = invExpDist(p, λ) sExpDistの逆関数。累積確率がpとなるtの値
ポアソン分布・確率密度 p = pPoissonDist(x, λ) 平均λの発生事象がx回発生する確率
ポアソン分布・累積確率 s = sPoissonDist(x, λ)
ポアソン分布・累積確率の逆関数 x = invPoissonDist(p, λ)
三角分布・確率密度 p = pTriDist(x, xmin, xmax, xmod) 最小値xmin, 最大値xmin, 最頻値xmin
三角分布・累積確率 s = sTriDist(x, xmin, xmax, xmod)
三角分布・累積確率の逆関数 x = invTriDist(p, xmin, xmax, xmod)
●乱数系
一様分布・連続 x = xRandomUni(xmin, xmax)
一様分布・整数 i = iRandomUni(imin, imax)
正規分布・連続 x = xRandomNormal(μ, σ) 正規分布累積確率逆関数 invNormaiDist を利用
正規分布・連続 x = xRandomNormal12(μ, σ) 12個のMath,Random()の和が正規分布になることを利用
正規分布・整数 i = iRandomNormal(μ, σ)
正規分布・連続・区間指定 x = xRandomNormalMinmax(xmin, xmax, μ, σ)
正規分布・整数・区間指定 i = iRandomNormalMinmax(imin, imax, μ, σ)
指数分布・連続 x = xRandomExp(λ) 単位時間内の発生回数が平均λの事象が、x単位時間に発生する確率の乱数
指数分布・整数 i = iRandomExp(λ)
指数分布・連続・区間指定 x = xRandomExpMinmax(xmin, xmax, λ)
指数分布・整数・区間指定 i = iRandomExpMinmax(xmin, xmax, λ)
三角分布・連続 x = xRandomTri(xmin, xmax, xmod)
三角分布・整数 i = iRandomTri(xmin, xmax, xmod)

代数方程式 f(x) = 0 の根: x = xFx(fx, xmin, xmax[, eps])

fx= xmin= xmax= eps= 解=
はさみうち法による f(x) = 0 の解(参照:「方程式の根」)を求めます。
xmin~xmaxの間に唯一の実根があることを前提にしています( f(xmin) と f(xmax) が同符号ではいけません。
二分法xmin~xmaxで範囲を狭めていき、相対誤差 |f(x)/F| < eps (F=|f(xmin)| + |f(xmmax)|)になるまで繰り返します。epsを省略すると、0.0001にしています。

方程式の根(ニュートン法) x = xFxNewton(fx, x0{, eps])
  トレース版 x = xFxNewtonTrace(表示場所, fx, x0{, eps])

ニュートン法により、方程式の実根を求める。
   xn = xn-1 - f(xn-1)/f'(xn-1)   (f'=df/dx)
を繰り返して根に近づく。参照:「方程式の根」
 x0は初期値。その値により異なる解に達したり、解が得られないことがある。

トレース版では、表示場所を次のようにする。
  メイン:  xFxNewtonTrace("〇〇", fx, x0, eps);
  HTML: <div id="〇〇">ここに表示される</div>

fx= x0= eps= 解=

3次方程式の複素数根 a0 + a1x + a2x2 + a3x3 = 0
    a3を省略すると2次方程式、さらにa2を省略すると1次方程式になる。

2次方程式は根の公式による。3次方程式は、特殊なケースを除き、ニュートン法により一つの実根を求め、式の変形により2次方程式にして求める。
  xfx3(2, 1)    x + 2 = 0     1次方程式
  xfx3(1, 1, 1)   x2 + x + = 0    2次方程式 根の公式
  xfx3(0, 1,-1, 1) x(x2 -x1 +1) = 0  ( )内は2次方程式
  xfx3(2, 0, 0, 1) x3 = -2      3乗根の公式を使う
  xfx3(2, 0,-1, 1) x3 - x2 + 2 = 0   一般的な3次方程式

戻り値は2次元配列である。
   var ary = xfx3(~);
とすると、
   ary[i][0] 解iの実数部、ary[i][0] 虚数部
が入る。

a0 + a1x + a2x2 + a3x3 = 0

x0= + i
x1= + i
x2= + i

数式 f(x) を最小にする x = xminFx(fx, xmin, xmax [, eps])

fx= xmin= xmax= eps= 解=

数式 f(x) が最小となるxの値を戻します。
 xmin~xmaxを等分してdxとし、for (x=xmin, x<=xmax; x=x+dx) でf(x) を最小にする x=xpot を求め、次に xopt-dx~xopt+dx で同様な操作を繰り返すことにより、最終的な xopt を求めます。epsを省略すると0.0001になります。

数式 f(x) を最大にする x = xmaxFx(fx, xmin, xmax [, eps])

fx= xmin= xmax= eps= 解=

数式 f(x) が最大となるxの値を戻します。
 xmin~xmaxを等分してdxとし、for (x=xmin, x<=xmax; x=x+dx) でf(x) を最大にする x=xpot を求め、次に xopt-dx~xopt+dx で同様な操作を繰り返すことにより、最終的な xopt を求めます。epsを省略すると0.0001になります。

代数式 f(x) の微分係数: x = dFx(fx, x [, d, eps])

fx= x= d= eps= 解=

fxのxにおける微係数aを求めます。xの近傍の2点 x-d*xとx-d*xでの勾配をa0とし、dを1/2にしたときの勾配aとの差が、相対誤差 |(a-a0)/a| < eps になるまで繰り返します。
dを省略すると0.01、epsを省略すると0.001にしています。

代数式 f(x) の定積分(シンプソン則): x = sFx(fx, xmin, xmax [, n2, eps])

fx= xmin= xmax= n2(キザミ数)= eps= 解=

シンプソンの公式(参照:「定積分」)により計算しています。
最初のn2で計算した値をx0、次にn2を2倍(すなわちキザミ幅を1/2)にして計算した値をxとして、相対誤差 |(x-x0)/x| n2を省略すると20、epsを省略すると0.0001にしています。

ラグランジュの補間法 polArray(n, m, xa, ya)

n+1組の点xa[0]~xa[n]、ya[0]~ya[n] を通るn次式をラグランジュの補間法により求め、m 個の、xa[n+1]~xa[n+1} における ya[n+1]~ya[n+m] を戻します。

4(=n+1)点(0,1), (2,5), (3,8), (5,10) を通る3(=n)次式により、
2(=m)点、ya[4] = y(4), ya[5] = y(6) を求める
    var n = 3;
    var m = 2;
    var xa = [0, 2, 3,  5,    4, 6];
    var ya = [1, 5, 8, 10];
    polArray(n, m, xa, ya);
解: ya[4] = y(4) = 10.1 ya[5] = y(6) = 6.6
n= m=
xa[0]= xa[1]= xa[2]= xa[3]= xa[4]= xa[5]= xa[6]= xa[7]= xa[8]= xa[9]=
ya[0]= ya[1]= ya[2]= ya[3]= ya[4]= ya[5]= ya[6]=

ラグランジュの補間法 y = polPoints(x, x0,y0, x1,y1, x2,y2, ..., xn,yn)

n+1個の点、(x0,y0),(x1,y1),(x2,y2)~(xn,yn) を通るn次式をラグランジュの補間法により求め、xにおけるn次式の値を戻す。
  y = polPoints(2, 1,2, 3,6) (1,2)と(3,6) を通る1次式の x=2 の値が y に入る
  y = polPoints(2, 0,0, 2,0, 3,3) 3点 (0,0).(2,0),(3,3)を通る2次式の x=-1 値が y に入る

 
polPoints()

最小二乗法(高次式)LeastSquare(n, p, f, xa, ya, aa)  (トレース版)LeastSquareTrace(表示場所, n, p, f, xa, ya, aa)

p組の(x,y)から、最小二乗法により、n次の近似式 y=a0+a1x+a2x2+...+anxn の係数 ai を求め。この近似式により、f個の x に対応する y を計算します。x,y,a は、それぞれ一次元配列 xa, ya, aa でメイン側で定義します。
参照:「回帰分析・最小二乗法」「多項式、高次式での最小二乗法」

近似方程式次数(n)= 実績組数(p)= 予測期間(f)=
xa=
ya=

配列のソート(行を変える) sortMatrix(a, keycol [, order [, fixrow]])

配列 a[i][j] を、j=keycol の列の大小によりソートします。
order:省略時あるいは正の値なら昇順、負の値なら降順
fixrow:省略時あるいは0ならすべての行を対象。先頭のn行を固定する(ソートの対象にしない)ときはnの値
(fixrowを与えるときは、orderは省略できない)

入力データ
a[0] = ["C", 30];
a[1] = ["B", 20];
a[2] = ["E", 50];
a[3] = ["A", 10];
a[4] = ["D", 40];
keycol = 1;

order= fixrow=

連立一次方程式(ガウス=ジョルダンの掃出し法)GaussJordan(n, a, b, x)

出力配列 x に解が入ります。入力配列 a と b は、関数呼出後も値が変わりません。
 アルゴリズムは「連立方程式」を参照。誤差の考慮はしていないので、nが大きな問題は解けません。

    function GaussJordanTest1() {
        var n = 3;
        var a = new Array();
        a[0] = [0, 3, 2];       //  0x0 + 3x1 + 2x2 = 12
        a[1] = [2, 4, 0];       //  2x0 + 4x1 + 0x2 = 10
        a[2] = [3, 1, 5];       //  3x0 + 1x1 + 5x2 = 20
        var b = [12, 10, 20];
        var x = new Array();    //  出力配列(解が入る配列)もメイン側で定義しておく
        GaussJordan(n, a, b, x);
        alert("x[0]=" + x[0] + ", x[1]=" + x[1] + ", x[2]=" + x[2]);  // 解の表示
    }
変数(式)の個数n=
x0 x1 x2 x3 x4
0
1
2
3
4

連立一次方程式(特殊入力形式) GaussJordanArray(a)
 (トレース版) GaussJordanArray(表示場所, a)

GaussJordan() と同様のガウス-ジョルダン法を用いた連立一次方程式の求解ですが、入力データの形式が異なります。
 一つの二次元配列 a にすべてのデータをまとめます。

   var a = new Array();
   a[0] = [ 3,  0, 0, 0];  // n=a[0][0]:変数・式の個数
   a[1] = [12,  0, 3, 2];  // 12 = 0*x1 + 3*x2 + 2*x3
   a[2] = [10,  2, 4, 0];  // 10 = 2*x1 + 4*x2 + 0*x3
   a[3] = [20,  3, 1, 5];  // 20 = 3*x1 + 1*x2 + 5*x3
   GaussJordanArray(a);    // 引数は a だけ
解は、a[i][0] に入ります。
   x1 = a[1][0] = 1
   x2 = a[2][0] = 2
   x3 = a[3][0] = 3

関数実行後 a は壊されます(単位行列になっています)。
 アルゴリズムは「連立方程式」を参照

n=
定数項12345

逆行列 invMatrix(n, a [, x])
 トレース版 invMatrix(n, a)

行・列のサイズn の配列 a[i][j] の逆行列を、ガウス=ジョルダンの掃出法で求めます。メイン側で定義した二次元配列 x を与えれば、x に逆行列が入り、a は壊されません。xを省略すると a が逆行列になります。

  n = 3;            invMatrix(n, a, x);のとき         invMatrix(n, a);のとき
    a[0] = [ 1, 2. 1];     a[0} = [ 1, 2. 1] x[0] = [ -0.4, 0.6, 0.2]   a[0] = [ -0.4, 0.6, 0.2]
  a[1] = [ 2, 1. 0];  →  a[1] = [ 2, 1. 0] x[1] = [  0.8,-0.2,-0.4]   a[1] = [  0.8,-0.2,-0.4]
  a[2] = [ 1, 1, 2];     a[2] = [ 1, 1, 2] x[2] = [ -0.2,-0.2, 0.6]   a[2] = [ -0.2,-0.2, 0.6]
  x = new Array();が必要か?    必要                     不要
n=

線形計画法: ary = LP(m, n, a);  トレース版 any = LPTrace("表示場所", m, n, a)

トレース版を使うときは、
script: var ary = LPTrace("〇〇", m, n, a);
  HTML:  <div id="〇〇"></div>
とします。

シンプレックス法により線形計画法を解きます。それに加えてシャドウプライス(レジューズドコスト)や価格・定数項のレンジなどの感度分析もしています。参照:「シンプレックス法」「シャドウ・プライスとレンジ」

            変数
       0 1 2 ...  n  
      ┌─┬─────────┐   m,n,aは左図を参照してください。
     0│ │  目的関数   │
      ├─┼─────────┤   a[0][0] には "max"(最大化) "min"(最小化)が入ります。
     1│ │         │
   制 2│定│         │   定数項はすべて正数で、定数項≧~ の形式だけを扱います。
   約  │数≧         │
   式 :│項│         │   誤差の累積などは考慮していません。
     m│ │         │   小規模な問題(mもnも10以下)を想定しています。
      └─┴─────────┘

関数内部では、次のようなaを拡張した作業用配列を用いるので、aの値は保持されます。
                                 nx+1          n=nx+m
             0      1     2      ↓ 3     4    5 ↓
          +------+-------------+-------------------+-----
        0 |      |             |                   |  :←目的関数
          +------+-------------+-------------------+  : 
        1 |      |             |                   |  :
          | 定   |  基幹部     |   スラック部      | m+1
        2 | 数   |             |                   |  :
          | 項   |             |                   |  :
    m = 3 |      |             |                   |  :
          +------+-------------+-------------------+-----
                 |----- nx ----|-------- m --------|
          |--------------  nx + m + 1 -------------|

解は2次元配列になります。
メイン側で、 ary = LP(m, n, a); とすると、
ary の内容は次のようになります。
   ary[3][1] = 72, ary[2][5] = 3

  ↓j k→   0   1    2    3      4    5
   0   変数名 目的関数 解の値 レジューズド レンジ レンジ
           定数項      コスト   下限  上限
   1   X1    3    16    0      2    8
   2   X2    2     8    0      0.75   3
   3   λ1    72     0    0.333    60   96
   4   λ2    48    0    0.833    36   52.364
   5   λ3    48    8    0      NaN   NaN

変数名は、基幹変数(与えた変数)はXj、制約式(スラック変数)はλiになります。
基底に入らないX、基底に入ったλのレンジは計算していません。

変数の個数: 制約式の個数:
目的関数

制約条件







輸送問題: ary = TLP(m, n, a);  トレース版 any = TLPTrace("表示場所", m, n, a)

供給地iの供給量a[i][0] (i=1~m), 需要地jの需要量a[0][j] (j=1~n) と、i→jの単位輸送費用a[i][j]を与えて、総輸送費用を最小化するように、i→jの最適輸送量を計算します。
 北西隅の方法により初期解を得て、U-V法により解の改善を行い、最適解を求めます。
 すべての入力値は正の整数で、Σa[i][0] = Σa[0][j] でなければなりません。

入力データ
  供給地数(m) = 3, 需要地数(n) = 4
  a[i][j]  j→ 0   1    2    3    4=n 
             需要地1 需要地2 需要地3 需要地4
  ↓i           12    20    16    18 ←需要量
  1 供給地1   18   7    3    2    10
  2 供給地2   22   9    3    6    8
  3 供給地3   26   8    7    8    6
  ↑m       ↑供給量
関数呼び出し
  var ary = TLP(m, n, a);
        トレース版を使うときは
          script: var ary = TLPTrace("〇〇", m, n, a);
        HTML:  <div id="〇〇"></div>
戻り値(aryの内容)は、ary[i][j]の2次元配列になります。
 i→   0   1   2
 ↓  供給地 需要地 輸送量 
  0   1   1   2
  1   1   3   16
  2   2   1   2
  3   2   2   20
  4   3   1   8
  5   3   4   18
  ↑
  ary[0].length = 6
総費用 = Σ(a[ary[0][j]][ary[1][j]]*ary[2][j]) = 296
供給地の個数: 需要地の個数:
需要地1需要地2需要地3需要地4需要地5需要地6需要地7
供給量↓需要量→
供給地1
供給地2
供給地3
供給地4
供給地4

基本統計量・一変数
  リスト入力 ary = basicStat(x1, x2, ..., xn)
  配列入力  ary = basicStatArray(x[, n])

基本統計量を得る。
リスト入力の場合 var ary = basicStat(40, 20, 50, 10, 30);
配列入力の場合  var x = [40, 20, 50, 10, 30];
         var ary = basicStatArray(x);   配列xの全要素が対象
         var ary = basicStatArray(x, 4);  先頭の4要素が対象(30 は対象外)
ary は統計量名称が要素とした連想配列となる。
ary["標本数"} =   5    n
ary["合計"} =    150   Σx
ary["平方和"} =   5500   Σx2
ary["平均"} =    30    μ = Σx/n
ary["偏差平方和"} = 1000   Σ(x-μ)2 = Σx2 - (Σx)2/n<
ary["不偏分散"} =  250   σ2 = 偏差平方和/(n-1)
ary["標準偏差"} =  15.81  σ = √不偏分散
ary["変動係数"} =  0.527  標準偏差/μn
ary["最小値"} =   10
ary["最大値"} =   50
ary["範囲"} =    40    最大値 - 最小値
ary["中央値"} =   30    メディアン
参照:「平均、分散、標準偏差」

基本統計量・2変数
  リスト入力 ary = basicStat(x1,y1, x2,y2, ..., xn,xn)
  配列入力  ary = basicStatArray(x[, n])

基本統計量を得る。
リスト入力の場合 var ary = basicStat2(26,54, 25,62, 23,51, 27,58, 28,63,
                    25,65, 22,59, 26,63, 25,65, 23,60);
配列入力の場合  var x = [26, 25, 23, 27, 28, 25, 22, 26, 25, 23];
         var y = [54, 62, 51, 58, 63, 65, 59, 63, 65, 60];
         var ary = basicStatArray(x);   配列xの全要素(10組)が対象
         var ary = basicStatArray(x, 8);  先頭の8要素が対象((25,65) (23,50) は対象外)
ary は統計量名称が要素とした連想配列となる。
  any["標本数"] =    10    n
any["平均X"] =    25    μx = Σ(x)/n
any["平均Y"] =    60    μy = Σ(y)/n
any["平方和X"] =   6282   Σ(x2)
any["平方和Y"] =   36194  Σ(y2)
any["積和"] =     15023  Σ(x・y)
any["偏差平方和X"] = 32    Sxx = Σ(x-μx)2
any["偏差平方和Y"] = 194   Syy = Σ(y-μy)2
any["偏差積和"] =   23    Sxy = Σ(x-μx)(y-μy)
any["共分散"] =    2.556  Sxy/(n-1)
any["回帰係数X"] =  0.718  Bx = Sxy/Sxx  y = a + bx のb
any["回帰係数Y"] =  0.119  By = Sxy/Syy  x = a + by のb
any["寄与率"] =    0.085  r2 = Bx・By
any["相関係数"] =   0.291  r = √(r2)
X、Yは大文字全角

正規分布:確率密度関数 p = pNormalDist(x, [μ, σ])

x= μ= σ= 解=

平均μ、標準偏差σの正規分布において、値xの確率密度関数の値pを戻します。
  p = 1/√(2πσ2) * e(-(x-μ)2/2σ2)
μ、σ2を省略すると、μ=0、σ2=1の標準正規分布になります。

正規分布:累積確率関数 標準形(0,1)s = sNormalDist(x)  ⊆ 一般形(μ、σ) s = sNormalDist(x, [μ, σ])

x= μ= σ= 解=

平均μ、分散σ2の正規分布において、-∞~xの確率密度関数の値sを戻します。μ、σ2を省略すると、μ=0、σ2=1の標準正規分布 sNormalDistStd() になります。
上の一般形の数値での意味:男性成人の平均身長は170cmで標準偏差が10cmであるとする。160cm以下の割合は? 解:15.9%

sNormalDist()は、確率密度関数 pNormalDist() を数値積分 sFx() で計算してもよいのですが、直接計算する近似式が知られています。ここではHastingsの式を使いました。
(引用先:http://qiita.com/gigamori/items/e17e6f9faffb78822c56  by Yanai Takamichi)

正規分布:区間内の累積確率 s = sNormalDistMinmax(xmin, xmax [, μ, σ])

xmin= xmax= μ= σ= 解=

正規分布(μ、σ)での xmin~xmax の累積確率 s を戻します。
 例えば、「身長が平均(μ)170cm、標準偏差(σ)10cmの正規分布に従うとき、身長が160cm~180cmの割合は? → 答 = 68.3%

 -∞~xの累積確率を求める関数 sNormalDist() を用いて、s = sNormalDist(xmax) - sNormalDist(xmin) をしているだけです。
 μ、σを省略すると(0、1)の標準形正規分布になります。

正規分布:累積確率分布の逆関数
標準形: x = invNormalDist(p)  ⊆ 一般形: x = invNormalDist(p, [μ, σ])

p= μ= σ= 解=

正規分布の累積確率関数 sNormalDist() の逆関数です。累積確率がpになる点xを戻します。
 例えば、「身長が平均(μ)170cm、標準偏差(σ)10cmの正規分布に従うとき、確率(p)90%の人はxcm以下である」かを求めるとx=183cmだとなります。

標準形 invNormalDistStd() は、μ=0、σ=1の場合です。 sNormalDist() でμとσを省略した場合も同じです。
 p=0.01 → x=-2.32, p=0.05 → x=-1.65, p=0.50 → x=0, p=0.55 → x=1.65, p=0.99 → x=2.32
となります(p<0.5 のときは x が負になります)。
 invNormalDistStd() での逆関数の近似式は、「https://potatodigger.wordpress.com/2013/06/30/正規分布の確率密度関数pdf、累積分布関数cdf、逆累/  by potatodigger」を参考にさせていただきました。

指数分布

関数一般形式機能計算式
p = pExpDist(t, λ)t=t における確率密度p = λ*e(-λ*t)
s = sExpDist(t, λ)区間 t=0~t の累積確率s = 1-e(-λ*t)
s = sExpDistMinmax(tmin, tmax, λ)区間 t=tmin~tmax の累積確率s = e(-λ*tmin)-e(-λ*tmax)
t = invExpDist(λ, t)累積確率が p となる t の値t = -(1/λ)log(1-p)

参照:「ポアソン分布と指数分布」

λ=
t= p=
t= s=
tmin= tmax= s=
p= t=

ポアソン分布

確率密度 p = pPoissonDist(x,λ), 累積確率 s = sPoissonDist(x,λ), 累積確率の逆関数 x = pPoissonDist(p,λ)

一定(微小)時間内にある事象が発生する回数が平均λであるとき、発生回数がxとなる確率pは、
    p(x) = (λx/x!)*e ← pPoissonDist(x,λ)
のポアソン分布に従う(参照:「ポアソン分布と指数分布」
 累積確率は、単純に s(x) = s(x-1) + p(x) で求めた。
 逆関数は、x=0, 1, 2, ... として累積確率 s(x) を計算し、s(x) ≦ p である最大のxを戻している。すなわち「確率pでx回以下だといえる」こととした。

λは正の実数、pは0≦p≦1、xは0以上の整数でなければならない。

λ= x= p= s=
λ= p= x=
   

二項分布 二項係数:m = nCx(n, x) 確率密度:p = pBinDist(x, n, p) 累積密度:s = invBinDist(x, n, p)

発生確率がpの事象が、n回の試行でx回発生する確率。
  確率密度: x回発生する確率   P(x) = nCx・px(1-p)n-x = pBinDist(x, n, p)
  累積密度: 0~x回発生する確率 S(≦x) = ΣP(x) = sBinDist(x, n, p)
    nCx = n!/x!(n-x)! = nCx(n, x) 二項係数(n個からx個を取り出す場合の数(場合の数))
参照:「順列・組合せ」「二項分布」「2項分布の計算プログラム」

ここの例では、xを与えず、x=0~nの表を掲げます。

二項分布の計算は計算量が大きく、大きい数と小さい数の計算が必要になるので誤差が生じます。計算を容易にするために、
  ・np > 5 かつ n(1-p) > 5 ならば、平均μ=np、標準偏差σ=√np(1-p) の正規分布で近似できる。
   確率密度:区間 x-0.5~x+0.5 の正規分布累積確率 sNormalDistMinmax(x-0.5, x+0.5, μ, σ)
   累積密度:≦x+0.5 の累積確率 sNormalDist(x+0.5, μ, σ)
  ・n>50, np≦5 ならば、平均λ=np のポアソン分布で近似できる。
などの方法がありますが、ここでは用いていません。

p= n=

負の二項分布:確率密度 p = pBinNegDist(x,n,p)、累積密度 s = sBinNegDist(x,n,p)

負の二項分布とは、ある事象が起こる確率がpのとき、n回試行したときに、その事象がはじめてx回になる確率P(n)の確率分布です。
   P(n) = n-1x-1・px(1-p)n-x
で計算できます(参照:「二項分布」)。
x=1 のときは、幾何分布 pGeoDist と一致します

x= n= p=

幾何分布: 確率密度 p = pGeoDist(n, p)、 累積確率 s = sGeoDist(n, p)

発生確率pの事象がはじめて発生したときの試行回数nの確率分布(負の二項分布で、x=1とおくと、幾何分布に一致)
  確率密度 p = pGeoDist(n, p) = (1-p)n-1*p
  累積確率 s = sGeoDist(n, p) = 1-(1-p)n
参照:「二項分布」

n= p=

超幾何分布
  確率密度 p = pHyperGeoDist(a, b, α, β)
  累積確率 s = sHyperGeoDist(a, b, α, β)

確率密度:袋の中にAがa個、Bがb個入っている。α+β個を取り出したとき、Aがα個、Bがβ個である確率
   p = aαbβa+bα+β
で計算できる(参照:「二項分布」)。
 確率密度:αが0~α個である確率。「α以上である確率」は sHyperGeoDist(a, b, α-1, β+1) になる。

  A  B
袋の中の個数: a = b =
取り出した個数: α= β=

三角分布

3点 (xmin,0), (xmod, h), (xmax,0) を結ぶ三角形となる分布です。
ここでは、いろいろな点の値を計算して表にしています。

xを整数値にしたいときがあります。連続量の場合は、xmin および xmaxでは発生確率が0になりますが、整数の場合はその点でも0ではない値であるのが通常です。そのときは、やや乱暴ですが、xmin-0.5~xmax+0.5とします。そして累積確率では、x+0.5 の点を累積確率とし、逆関数のときは得られた値を四捨五入ば、それらしい値が得られます。

最小値:xmin= 最大値:xmax= 最頻値:xmod=

一様分布乱数:連続 x = xRandomUni(xmin, xmax), 整数 i = iRandomUni(imin, imax)

下の表では、連続のときは、得たxを四捨五入して表示しているので、最小値と最大値の割合は他の半分になっています。
整数のときは、区間を xmin~xmin+1として、乱数を発生させ、得たxを切り捨てているので、同じ値になっています。

最小値= 最大値= 発生数=
  

正規分布乱数
  出力    区間指定なし          区間指定あり
  連続量 x = xRandomNormal(μ, σ)  x = xRandomNormalMinmax(xmin, xmax, μ, σ)
  整数値 i = iRandomNormal(μ, σ)  i = iRandomNormalMinmax(imin, imax, μ, σ)
  連続量 x = xRandomNormal12(μ, σ)

0~1の値を持つ乱数から、正規分布の累積確率の逆関数を用いて、値xを得ます。
 整数値のときは、単に四捨五入しているだけです。区間指定ありの場合は、区間外の値になったときは。区間内になるまで繰返し乱数を発生させます。

xRandomNormal12はxRandomNormalと全く同じ機能ですが、考え方が異なります。
中心極限定理により、多数の一様分布に従う標本の和は正規分布になります。Math.random()は0~1の一様分布です。それを12個を合計した値は、μ=6、σ=1の正規分布に近づきます。xRandomNormal(μ, σ)が内部でinvNormalDistを用いているのに対して、xRandomNormal12は非常に簡単なアルゴリズムで実現しています。精度は、xRandomNormalのほうが良いようです。

μ= σ= 発生数=
最小値= 最大値=
解=
  
  
  
  

指数分布
        区間指定なし      区間指定あり
   連続量 x = xRandomExp(λ)   x = xRandomExpMinmax(xmin, xmax, λ)
   整数値 i = iRandomExp(λ)   i = iRandomExpMinmax(imin, imax, λ)

単位時間中にある事象が発生する平均回数をλとするとき、その事象の発生間隔がx単位時間である確率は指数分布に従う。
 指数分布累積確率の逆関数は、x = -1/λ*log(1-p) であるから、指数分布乱数は、x = -1/λ*log(1-Math.random()) で得られる。xの平均値は 1/λ になる。
 なお、事象発生は稀であるとの前提から、λの値が小さくなるよう単位時間を設定することが望ましい。
参照:「ポアソン分布と指数分布」

ここでの出力について

回数=
λ=
最小値= 最大値=
平均= 最頻値=
  
  
  

三角分布: 連続量 x = xRandomTri(xmin, xmax, xmod)
      整数値 i = iRandomTri(xmin, xmax, xmod)

xmodは最頻値、平均は (xmin+xmax+xmod)/3 になる。
 三角分布の累積確率の逆関数を利用してもよいが、二つの一様分布乱数、u, v (u>vとする)として、
   (1-c)*v + c*u    c = (xmod-xmin)/(xmax-xmin)
とすることにより、xmin=0, xmax=1, xmod=c の三角分布乱数が得られる。このほうが簡単なので、ここではこれを採用した。

回数= xmin= xmax= xmod=
平均= 最頻値=