大学初年級の数学計算までならば、大多数を配布しているデフォルトの Python sf ワンライナー計算式で書ける。それを具体例で示す。たまたま面白そうな未定係数法に関連する問題が 大学数学のラグランジェ未定数乗数法が分かりません」で質問されていた。これを題材に使って Pythn sf one-liner 数式の実用性を例示する。

以下は私自身の問題を解く過程での思考過程を、誤った小道に入り込んで戻る過程も含めて示す。人間の思考手順とワンライナー計算が互いに影響しあう様子を示せるからだ。同時にこれが Python sf の short tour にもなっている。


# 半径 1 の球および それと接する平面がと原点が作る様子を図示する
dct={}; for idx,(φ,θ) in enmitr(klsp(0,2pi), klsp(0,pi)):dct[idx]=[sin(θ) cos(φ),sin(θ) sin(φ), cos(θ)];tempDct1:= dct

dct={}; vx,vy,vz = sqrt(3) kzrs([3,3])^0; prm=klsp(0,1);for idx,(u,v) in enmitr(prm, prm):dct[idx]=u vx + v vy + (1-u-v) vz;tempDct2:=dct

=:tempDct1, tempDct2; renderFaces(tempDct1, blMeshOnly=True); renderFaces(tempDct2, blMeshOnly=True, meshColor=red);drawAxis(2)



球または三角形の板を表示するだけならば、一つの Python sf on-liner で表示させられる。でも球と三角の二つを一つの one-liner で表示させるのはできるけれど、 one-liner とは言いにくい長さになってくる。なので下の三つの one-liner で一つの図を表示させるようにした。One-liner として扱うのはディスプレーの一行の範囲で表示できる範疇までだと思う。

  • 球のデータを作る one-liner
  • 三角データを作る one-liner
  • 球と三角と x,y,z 軸を表示する one-liner

この oneliner で klsp(0,pi) は下のようなベクトルデータを作っている。ようは [0,pi] の範囲を 50 等分するベクトルだ。0 と pi の両端を含む長さが 51 のベクトルデータなのwで、図示するときに端が欠けることがなくて便利だ。


Python sf expression
klsp(0,pi)
===============================
[ 0.          0.06411414  0.12822827  0.19234241  0.25645654  0.32057068
  0.38468481  0.44879895  0.51291309  0.57702722  0.64114136  0.70525549
  0.76936963  0.83348377  0.8975979   0.96171204  1.02582617  1.08994031
  1.15405444  1.21816858  1.28228272  1.34639685  1.41051099  1.47462512
  1.53873926  1.60285339  1.66696753  1.73108167  1.7951958   1.85930994
  1.92342407  1.98753821  2.05165235  2.11576648  2.17988062  2.24399475
  2.30810889  2.37222302  2.43633716  2.5004513   2.56456543  2.62867957
  2.6927937   2.75690784  2.82102197  2.88513611  2.94925025  3.01336438
  3.07747852  3.14159265]
---- ClTensor ----

二つの klsp(..) ベクトルを引数とする enmitr(..) 関数は Python sf が用意している下のような複数ループを一つで済ますための generator 関数だ。行列のための整数インデックスも組で返しており、今回のような三次元に分布する二次元データを作るのに便利だ。Python sf では変数にギリシャ文字が使えるので、リスト内包表記が分かっている理数系の Python 使いならば「for idx,(φ,θ) in enmitr(klsp(0,2pi), klsp(0,pi)):dct[idx]=[sin(θ) cos(φ),sin(θ) sin(φ), cos(θ)]」が球面データを生成していることが読み取れるはずだ。数学では、このような単純な多重ループが頻発する。多重の for ループは Python one-liner で書けない。enumitr(..) 関数のような小さな仕掛けを使うことで多重ループでも One-liner が使え、また可読性を保てる。


Python sf expression
sc.source(enmitr)
    """ 列挙インデックス付きの多次元繰り返しイタレータ
        multiple dimentional iterator with enumerate index
    e.g.
    s=set(['a','b']);list(enmitr(s,s))
    ===============================
    [((0, 0), ('a', 'a')), ((0, 1), ('a', 'b')), ((1, 0), ('b', 'a')), ((1, 1), ('b', 'b'))]

    s=set([1,2]);list(enmitr(s,s,s))
    ===============================
    [((0, 0, 0), (1, 1, 1)), ((0, 0, 1), (1, 1, 2)), ((0, 1, 0), (1, 2, 1)),
     ((0, 1, 1), (1, 2, 2)), ((1, 0, 0), (2, 1, 1)), ((1, 0, 1), (2, 1, 2)),
     ((1, 1, 0), (2, 2, 1)), ((1, 1, 1), (2, 2, 2))]
    '"""
    head, tail = args[0], args[1:]

    if 'index' in kwarg.keys():
        index = kwarg['index']
    else:
        index = ()

    if type(head) in [int, long, float]:
        head = range(int(head))

    if tail:
        if len(tail) == 1 and hasattr(tail[0],'next'):
            # to avoid multiple use of one iterator
            tailAt = (tuple(tail[0]), )
        else:
            tailAt = tail

        for idxI, i in enumerate(head):
            for idxJ, j in enmitr(index=index, *tailAt):
                if len(tail) == 1:
                    yield index + (idxI,idxJ),  (i, j)
                else:
                    yield index+(idxI,)+idxJ, (i,)+j
    else:
        for idxI, i in enumerate(head):
            yield idxI, i

===============================
None
Python sf expression

上の Python sf 式で、 sc.source(..) 関数は、引数に与えた関数やクラスまたはモジュールなどの Python ソース・コード文字列を返す関数だ。Python sf では 90% 以上のソースを公開しており、分からない関数などがあったら、そのソース・コードを sc.source(..) を使うことで用意に覗ける。

なお Python sf ては sc ラベルに numpy モジュールを割り振っている。だから sc.source は numpy モジュールが備える関数を呼び出しているだけだ。この sc を経由することで numpy の全てを Python sf 式に記述できる。numpy の info(..) 関数を下のように使えば、下のようにソースを含まない enmitr の introspection 部分を取り出せる。


Python sf expression
sc.info(enmitr)
In file: D:\my\vc7\mtCm\pysf\basicFnctns.py

def enmitr( *args, **kwarg):
python -u sfPP.py "sc.info(enmitr)"
 enmitr(*args, **kwarg)

列挙インデックス付きの多次元繰り返しイタレータ
    multiple dimentional iterator with enumerate index
e.g.
s=set(['a','b']);list(enmitr(s,s))
===============================
[((0, 0), ('a', 'a')), ((0, 1), ('a', 'b')), ((1, 0), ('b', 'a')), ((1, 1), ('b', 'b'))]

s=set([1,2]);list(enmitr(s,s,s))
===============================
[((0, 0, 0), (1, 1, 1)), ((0, 0, 1), (1, 1, 2)), ((0, 1, 0), (1, 2, 1)),
 ((0, 1, 1), (1, 2, 2)), ((1, 0, 0), (2, 1, 1)), ((1, 0, 1), (2, 1, 2)),
 ((1, 1, 0), (2, 2, 1)), ((1, 1, 1), (2, 2, 2))]
'
===============================
None

sc.source や sc.info 関数が扱える対象は Python の全てだ。 Python sf 数式だけに限らない。Python sf が Python に対して upper-compatible だからだ。下のように wx モジュールの関数の関数のソースを覗くこともできる。


Python sf expression
import wx; sc.source(wx.Slider.SetTickFreq)
In file: D:\lng\Python26\lib\site-packages\wx\_controls.py

    def SetTickFreq(*args, **kwargs):
        """SetTickFreq(self, int n, int pos=1)"""
        return _controls_.Slider_SetTickFreq(*args, **kwargs)


import wx; help(wx.Slider.SetTickFreq)
SetTickFreq(*args, **kwargs) unbound wx._controls.Slider method
    SetTickFreq(self, int n, int pos=1)

===============================
None

下のような動作に自身を持てない Python コードを one-liner で簡単に試すこともできる。print 分を書かなくても、勝手に Python sf 式の最終式の値をコンソールに打ち出すので便利だ。実際に Python sf を使っていると、print debug では追いつかなくなったときに限って Python debugger は余り立ち上げるようになる。


Python sf expression
(1,2,3) == [1,2,3]
===============================
False

(1,2,3) == (1,2,1+2)
===============================
True

これらの使用例で分かるように Python sf の one-liner 式は Python sf 数式だけに限らず、Python コード一般に役立つ。

「tempDct1:= dct」のように Python には存在しない「:=」記号を使っている。これは Python sf one-liner 数式特有の記号であり、dct データでtempDct1 ファイル変数を作ることを意味する。実際にこの Python sf 式により tempDct1.pvl ファイルが作られている。 tempDct1, tmpDct2 ファイル変数は「=:」記号により Python sf 式内で読み戻され、Python オブジェクトとして扱える。すなわち Python sf 式の対象として扱える。「読み出されたファイル変数は「=:tempDct1, tempDct2; renderFaces(tempDct1…);…」のように表面描画関数 renderFaces(..) のデータとして与えられ、球や三角形のメッシュ図形が描かれる。


DOS command
dir tempDct1.pvl
 ドライブ D のボリューム ラベルがありません。
 ボリューム シリアル番号は 5C35-A5B7 です

 D:\my\vc7\mtCm のディレクトリ

2010/11/26  10:17           552,235 tempDct1.pvl
               1 個のファイル             552,235 バイト
               0 個のディレクトリ  80,429,912,064 バイトの空き領域

小沢一郎検察審査会のモンテカルロ・シミュレーション

小沢一郎の東京第五検察審査会11人の審査員の平均年齢が 30.9 歳/ 33.9歳 と公表され、それが起こりえない確立ではと論じられています。

この平均年齢の確率分布を、Python sf によるモンテカルロ・シミュレーションで計算します。

まず年齢分布のデータとして、これ使います。この左端の都民の年齢分布のデータを使います。

○ 東京都民の年齢分布を mtPopulation.pvl ファイル変数にする。mtPopulation.pyl を後の Python sf 式で何度も再利用する。後の Python sf 式は mtPopulation のみに文脈依存することになる。

//@@
putPv([
93600,
94675,
95104,
96624,
96689,
98882,
96400,
96520,
95391,
94189,
95569,
94794,
91678,
92579,
91973,
93776,
95875,
99257,
122509,
151551,
163834,
170901,
173758,
173652,
177597,
186120,
190836,
196995,
199645,
207634,
214539,
224749,
230934,
227457,
224010,
218822,
216925,
211284,
209625,
169360,
200967,
185644,
174963,
164762,
158810,
156817,
152725,
145950,
139638,
141526,
145809,
144130,
150631,
159957,
169527,
186817,
206177,
210925,
204355,
130395,
139546,
170477,
168667,
169831,
164901,
151532,
130880,
138083,
141234,
144215,
137808,
124627,
124319,
116547,
109099,
101779,
94921,
90497,
84778,
79382,
73248,
61693,
54876,
50681,
45240,
42453,
31918,
28931,
25745,
22723,
19735,
16572,
13416,
10640,
8134,
6282,
4491,
3190,
2285,
1358,
860,
564,
324,
206,
130,
58,
38,
19,
11,
3,
0,
1,
1,
0,
0
], 'mtPopulation')
//@@@
//copy c:\#####.### temp.py /y
//python sfPP.py -fs temp.py
//python -u __tempConverted.py

この年齢分布の Python sf でグラフ化します。

=:mtPopulation; plotGr(mtPopulation)

<== 18 を超えた所で急激に分布が変化している。これは現在でも 18 を超えた人間が全国から東京に移り住んでいることを意味する。

—————–

この年齢分布から、randint 関数を使って、無作為に 11 人を選んだときの平均年齢を計算せねばならない。どうするか。

○手がかりに、年齢の累積分布を描いてみる。
v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];plotGr(vv)

この年齢の累積分布の右端は総人口になっている。その総人口までの整数値で random 変数を発生させてやれば、無作為に人間を選択することのシミュレーションができる。東京都の全人口に整数番号を割り振って、その整数番号を random に発生させることになる。

○東京都の総人口
=:mtPopulation; sum(mtPopulation)
===============================
12415786

randint(..) 関数は上限を指定して random ベクトルを発生させられる。こんな具合だ。0–9 までの整数値を 11 個 random に発生させられる。

seed(0); randint(10, size=11)
===============================
[5 0 3 3 7 9 3 5 2 4 7]
—- ClTensor —-

seed(0) は random データの種を指定している。これを省略すると、種をコンピュータが勝手に選ぶことになる。結果としてランダム・データが計算の度に変わる。計算の再現性が

seed(0);v=:mtPopulation; randint(sum(v), size=11)
===============================
[ 8325804 1484405 2215104 5157699 10638075 8222403 7644169 10608339
5853461 6739698 374564]
—- ClTensor —-

欲しいのは、このような大きなランダム整数値の分布ではなく、年齢分布のデータだ。先の累積分布関数のベクトル vv と、この番号分布のデータから年齢分布の関数に直す必要がある。こんな関数 f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1) で求められる。次のような具合だ。

;;v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); f(vv[20])
===============================
20

再帰λ関数に Python のキーワード引数を使って初期値を与え、また三項演算子を使っている。この Python sf 式に可読性を感じるには、相当 Python に慣れてもらう必要がある。でも逆に、この関数を for ループ文と if .. return 構文で書いたら、ワンライナーでは可読性がなくなってしまう。再帰関数によるループ処理に慣れていない方はワンライナーを諦めて、素直に Python sf ブロック式で書くべきだろう。できたら、この機会に再帰関数による繰り返し処理にも慣れて欲しい。この場合で言えば、三項演算子と再帰関数の組み合わせのクールさ・可読性を味わって欲しい。

ここは、上の Python sf 式を何回か計算してみることで年齢を求める関数 f の働きを納得してもらうことを前提に、ワンライナーでの記述を進める。下の様な具合だ。

;;v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); f(vv[21] )
===============================
21

20 才以上の都民をランダムに 11 人選ぶのは次の Python sf 式で可能だ。

seed(0);v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))]; randint(vv[20], vv[-1],size=11)
===============================
[10477273 3635874 4366573 7309168 10373872 9795638 8004930 8891167
2526033 4984452 7082099]
—- ClTensor —-

vv[20] で 20 歳以上の都民の背番号の最初になり、vv[-1] で都民の背番号の最後の番号になっていることを利用して randint(vv[20], vv[-1],size=11) と呼び出している。

この randint(..) 関数の前に、都民の背番号から年齢を求める関数 f を入れてやることで、下のように 11 人の年齢分布を求められる。

○20歳以上からランダムに 11 人を選び出したときの年齢分布;;seed(0);v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); ([f(x) for x in randint(vv[20], vv[-1],size=11)])
===============================
[66, 28, 31, 46, 65, 62, 51, 56, 22, 34, 45]

○20歳以上からランダムに 11 人を選び出したときの平均年齢;;seed(0);v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); sum([f(x) for x in randint(vv[20], vv[-1],size=11)])/11
===============================
46.0
<== 意外と高い年齢になる。

この平均年齢を求める試行を N 回繰り返したときの平均年齢の長さ N のベクトルは、randint(.., size=[N,11]) が N x 11 のランダム値行列を返すことより、次のような Python sf 式で計算できる。

seed(0);v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); [sum([f(x) for x in vvv])/11 for vvv in randint(vv[20], vv[-1],size=[N,11])]

python -m sfPP "N=10;seed(0);v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); [sum([f(x) for x in vvv])/11 for vvv in randint(vv[20], vv[-1],size=[N,11])]"
===============================
[46.0, 47.727272727272727, 41.81818181818182, 50.454545454545453, 49.63636363636
3633, 41.636363636363633, 53.545454545454547, 48.272727272727273, 48.18181818181
818, 44.909090909090907]

N=10^4 に増やしてやれば、その年齢分布の histogram を計算してやれば、都民から 20 歳以上の人間を 11 人選んだときの平均年齢の分布のしかたのモンテカルロ・シミュレーションしたことになる。numpy に histogram 関数が備わっているので、下のように計算できる。

10^4 回試行したときの平均年齢の分布を numpy.hestogram(..) を使って求める
;;N=10^4;v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); seed(0);sc.histogram([ sum([f(x) for x in vvv])/11 for vvv in randint(vv[20], vv[-1],size=(N,11))], range(115))
===============================
(array([ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 1, 3, 7, 7, 24, 40, 76, 119, 144,
214, 296, 374, 481, 534, 571, 680, 753, 774, 738, 686, 649, 555,
504, 438, 321, 275, 209, 174, 127, 83, 49, 48, 15, 13, 8,
5, 4, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0]),

array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12,
13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25,
26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38,
39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51,
52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64,
65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77,
78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90,
91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103,
104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114]))
<== N280 netbook で 2分程度
<== range(115) は histgram を計算させるときの slot 分布を指定する tuple データです。0 — 115 才まで一切ごとに区切った slot を指定しています。

10^4 回試行したときの平均年齢の分布hestogram(..) をグラフ化する
;;N=10^4;v=:mtPopulation; vv=[sum(v[:(k+1)]) for k in range(len(v))];f=λ n,k=0:k if vv[k]<=n<vv[k+1] else f(n,k+1); seed(0);plotGr(sc.histogram([ sum([f(x) for x in vvv])/11 for vvv in randint(vv[20], vv[-1],size=(N,11))], range(115))[0])

——————-
–***************–
——————-
この結果からすれば、11 人の平均年齢が 34 才であったとしても、めったにおきることではないことが解る。

平均年齢が 31 才以下になる確率は
(1+3)/sum([1, 3, 7, 7, 24, 40, 76, 119, 144, 214, 296, 374, 481, 534, 571, 680, 753, 774, 738, 686, 649, 555,504, 438, 321, 275, 209, 174, 127, 83, 49, 48, 15, 13, 8,5, 4, 1])
===============================
0.0004

平均年齢が 34 才以下になる確率は
(1+3+7+7+24)/sum([1, 3, 7, 7, 24, 40, 76, 119, 144, 214, 296, 374, 481, 534, 571, 680, 753, 774, 738, 686, 649, 555,504, 438, 321, 275, 209, 174, 127, 83, 49, 48, 15, 13, 8,5, 4, 1])
===============================
0.0042

です。この数値からすると、無作為に選択された 11 人から、年齢の高い人たちに辞退者が多く発生したことが推測されます。

——————-
–***************–
——————-
ここまで計算してきた Python sf 式は mtPopulation ファイル変数のみに依存する式です。何ヶ月もして、式の意味を忘れても、mtPopulation ファイル変数さえ残っていれば計算できます。ワンライナーですから、そのまま Python sf に計算させるだけです。

ワンライナー式の変形はモンテカルロ・シミュレーションを行おうとして私が考えていった式のままです。この Python sf ワンライナー式の変形は、私の思考過程になっています。Python sf 式ワンライナーを作り上げるのはプログラムを作るのではなく、式を考え合わせていくことが解ってもらえますでしょうか。このような モンテカルロ・シミュレーションを行うためにデバッグが必要ないことを理解してもらえますでしょうか。

Hello world!

Welcome to WordPress.com. This is your first post. Edit or delete it and start blogging!

html link test for Sky Drives

test

今日 2001.12.25 から blog を始めます。とりあえずテストです。