無機野

ここはお前の日記帳

matplotlibで二項分布を表示する(だけ)

先の記事で書いたコードだけどnumpy.random.binomial()で実際にサンプリングした数字からなるリストを作ってそっからヒストグラムつくるみたいなよく考えたら(よく考えなくても)頭悪いことしてたの今更気づいて、自分で単純なベルヌーイ試行だっつってんだからそのまま二項分布のグラフ書けばいいじゃねえかと思って書いた
range(x)が0から始まるリストを返すのにだけ気をつけて書くとこんな感じ

import matplotlib.pyplot as plt
import scipy.misc as scm #組合せ計算用

p = 1.5/100 #レジェンド
pack = 50
time = pack*8

def unko():
    aaa = []
    for j in range(time-1):
        elm = scm.comb(time,j) * (p**(j)) * ((1-p)**(time-j))
        aaa.append(elm)
    print(aaa,sum(aaa)) #表示される数は*100すると%表示になる。そっちのが見やすいか?
    plt.plot(aaa)

unko()

[0.014428036055035228, 0.043833398877860821, 0.088556816007048775, 0.13384665972638971, (以下略)
1.0
f:id:moppii:20161130115543p:plain

オッ

import matplotlib.pyplot as plt
import scipy.misc as scm #組合せ計算用

p = 1.5/100 #レジェンド
pack = 50
time = pack*8

def unko():
    aaa = []
    for j in range(time-1):
        elm = scm.comb(time,j) * (p**(j)) * ((1-p)**(time-j))
        aaa.append(elm)
    print(aaa,sum(aaa)) #表示される数は*100すると%表示になる。そっちのが見やすいか?
    plt.plot(aaa)
    plt.xlim([0,20]) #x軸を0から20までに絞る

unko()

<f:id:moppii:20161130115706p:plain

ちなみにこの50パックという数字は9800円課金してそのままガチャ回した時の確率である
この分布の平均と標準偏差はそれぞれ

m = time*p
s = (time*p*(1-p))**(1/2)
print(m-s,m+s)

3.5689508437713564 8.431049156228644

なので諭吉放りこむとおよそ3回に2回は3~8枚引ける計算になる
これやっぱリセマラすべきなんじゃ…

pythonで解くソローモデル 初級編

経済学を勉強していくとソローsolowモデルというものに出会う
例えば学部で学ぶケインジアンのISLMなんかは静学分析といって当期だけの均衡を見るが、現実問題として投資すれば投資するだけ来期のGDPとか消費とかが伸びることは想像に難くない
雑に言うとそーゆー時間経過の影響も見ましょうよっていうのがソローモデルである
ところでいきなりいろいろ未知数とかパラメータもりこんで複雑な連立方程式を解こうとしても挫折するのでシンプルなモデルからやってみましょうよ、ついでに自力で計算するのも面倒なのでpythonで説いてみましょうねというのがこのエントリの趣旨である
ちなみに自分は大して勉強してない(なんせ留年してる)ので話半分に聞いておくと吉だ


変数は世に(というかモデルに)たくさんあるが、内生変数をYとKに絞れば(つまりいろんなパラメータとか人口Lとか技術水準Aとかを常に既知とすると)2つの方程式で古典派の成長理論を表すことができる
とりあえずこれを最もシンプルなソローモデルとでも呼んでみよう

{ \displaystyle
Y_t = {{K_t}^{\alpha}} * ({A_t}*{L_t})^{1-\alpha}
}
{ \displaystyle
K_{t+1} = (1-\delta)*K_t + s*Y_t
}

前者は生産関数、後者は資本蓄積式である
繰り返すがここで未知数は各期のYとKのみ、便宜的に{ \displaystyle
A_t
}{ \displaystyle
L_t
}と表記しているがこれらはすべての期で一定、加えてα(資本分配率)やδ(資本減耗率)やs(貯蓄率)も一定とするとこれは0期のKの値だけ与えてやれば解ける
というわけで解いてみよう

# -*- coding: utf-8 -*-
    
import numpy as np
import matplotlib.pyplot as plt
 
t = 100 #何期までやるか
Y = np.empty(t)
K = np.empty(t)

alpha = 0.6 #コブダグラス型の労働分配率あるいは資本分配率。数字は適当
delta = 0.1 #資本減耗率。こちらも適当
s = 0.2 #貯蓄率。適当

K[0] = 1 #単位は便宜的に万円ということで
A = 10 #技術水準あるいは全要素生産性、今回は外生値。単位は知らん
L = 30 #人口。単位は便宜的に万人で。


class Solow:
    """最もシンプルなソローモデル"""
     
    def __init__(self,Y,K): 
        self.Y, self.K = Y,K
        
    def calc(self):
        for i in range(t-1):
            self.Y[i] = (self.K[i]**alpha) * ((A*L)**(1-alpha)) #コブダグラス型生産関数
            self.K[i+1] = (1-delta)*self.K[i] + s*self.Y[i] #資本蓄積
            
            
obj = Solow(Y,K)
ans = obj.calc()
print(Y)

plt.plot(Y,"blue")
plt.plot(K,"red")

まあこんな感じであろうか
numpyのempty関数でYとKのt期分の空きリストを作っておいて、0期のKと各パラメータを外生的に与えてやり、あとは適当に方程式を与えてやれば解ける、みたいな
実行結果

runfile('C:/Users/a/Desktop/temp1.py', wdir='C:/Users/a/Desktop')
[ 9.79148362 18.38708522 29.40164841 42.40192563 56.99399415
   (中略)
815.21230173 816.53530369 817.80614335 819.02684748 820.19936649
821.32557714 822.4072851 823.44622751 824.4440754 0. ]
f:id:moppii:20161125115540p:plain

グラフの赤線は資本Kの額、青線は所得Yの額である
数字は適当だが100期ちょい過ぎで均衡値に達するであろうことが見て取れる
こーゆーのなんて言うんだっけ、定差方程式?違ったらごめん


これがソローモデルの最低限な部分かな 普通はこれのAとKを

{ \displaystyle
A_{t+1}=(1+g^{A})*A_t
}
{ \displaystyle
L_{t+1}=(1+g^{L})*L_t
}

こんな感じに内生化する。({ \displaystyle
g^{A}
}{ \displaystyle
g^{L}
}はそれぞれ技術進歩率、人口成長率である、念のため)
お察しの通りこのモデルには家計や企業の意思決定が介在していない(あと政府と海外部門も)
普通ソローモデルと呼ばれるものはここに企業の利潤最大化を盛り込んで資本蓄積させたものをいう(ハズ)
そうして組んだモデルはほっといてもそのうち均衡点に達するというのが教科書的な新古典派経済学の主張である(たぶん)


ついでに家計の効用最大化を盛り込むとラムゼーモデルになり、さらに失業と技術ショックをとりいれてRBC、硬直性や市場の不完全性などケインズ的な要素を取り入れてDSGEというようにいろんな現代的なモデルの骨子となっているので勉強しましょうという話

気が向いたら続きます
おわり

食い延ばし

食い延ばしという技術がある
愚形複合ターツから鳴いて待ちを良くするアレである
ex. 34568から25チーでカン4から47待ちに変える、みたいな
ある程度麻雀に慣れてくると自分がなんの牌で食い延ばすべきかは分かってくるものである
重要なのは他家の食い延ばしを見破ることである
これについてはググってもそんなにヒットしないようなので食い延ばした時の捨て牌から待ちを読む公式的なものを備忘録的に記しておく
理想は公式として覚えずに相手の手牌を想像して読むことだがいちいちそんなことに頭のリソースを割くほど余裕のないことが多いのでやはり機械的に覚えてしまうことが楽ではある
ただ他の捨て牌読みにも言えることだが一点読みは基本的に不可能なので他に安牌が無い時の参考までにとどめておいたほうがいいだろう


というわけで以下は愚形含みの5枚形(?)で食い延ばした時の公式である

① 鳴いた牌のスジを切ったらその裏スジが待ち
三筒四筒五筒六筒八筒チー二筒五筒八筒待ち四筒七筒
25で鳴いてスジの8が出たらその裏スジの47があたりになっているパティーン
これは比較的覚えやすいか
あえて名付けるならカンチャン整理型その1ってところか
以下全パティーン
元の形 喰い取り牌 切る牌 待ち(裏スジ)
12346  3      6   25
23457  14     7   36
34568  25     8   47
45679  36     9   58
13456  47     1   25
24567  58     2   36
35678  69     3   47
46789  7      4   58
ちなみにこの形では5は切られない

②の1 鳴いた牌の±4を切ったらそのスジが待ち
三筒四筒五筒六筒八筒チー四筒八筒待ち二筒五筒
±4と書くとわかりづらいけど要は喰いとった牌と切った牌の間に順子ができる場合ってこと
元の形 喰い取り牌 切る牌 待ち(裏スジ)
12346  2     6   3
23457  3     7   14
34568  4     8   25
45679  5     9   36
13456  5     1   47
24567  6     2   58
35678  7     3   69
46789  8     4   7
スジのうち外側か内側かについては、そもそも食い延ばしが上家が捨てた牌を手元に2枚温存する手段であるので喰いとった牌で両面を構成するほうのスジと考えるとわかりやすいか

②の2 鳴いた牌の±4を切ったらその跨ぎスジが待ち

④ 5切ったら鳴いたのと別のスジ待ち
(147チーして打5なら待ち369、369チーして打5なら147と大まかにおぼえればいい。。。のか?)

疲れた

シャドバでリセマラすべきか

最近はまってるシャドウバースというカードゲーム
他のスマホゲーよろしくリセマラを行うことでよりレア度の高いカードを引くことができる
では実際リセマラはするべきであるか?
回答をするならば”嫌になるのでリセマラすべきではない”
そもそも手慰みでゲームしてんのにいきなり単純労働を強いるのはあまりにマゾヒスティックではないか
だいたいシャドバはいらないカード砕いて自分の欲しいカードを生成するシステムがあるのでそこまで右手力は問われないし毎日ログインしてればログボとミッションでわりと引ける
それでもリセマラをしたい向きのためにどの程度の枚数レジェンドを引けばリセマラ成功といえるのかちょっと計算してみた

import matplotlib.pyplot as plt
import random
import collections as col

a = []
p = 1.5 #虹確率
pack = 36

def risemara(x):
    time = pack * 8
    for j in range(x):
        leg = 0
        for i in range(time):
            count = random.uniform(0,100)
            if count < p:
                leg += 1
        a.append(leg)

risemara(100000) #リセマラ回数

plt.hist(a,bins=max(a))
print(pack,"パックを",x,"回リセマラした場合に引くレジェンド(確率",p,")の枚数")
aaa = col.Counter(a)
for key ,value in sorted(aaa.items()):
    print(key,"枚",value,"回、その確率は",(value/x)*100,"%")

f:id:moppii:20161104221508p:plain
36 パックを 100000 回リセマラした場合に引くレジェンド(確率 1.5 )の枚数
0 枚 1348 回、その確率は 1.348 %
1 枚 5728 回、その確率は 5.728 %
2 枚 12214 回、その確率は 12.214 %
3 枚 17905 回、その確率は 17.904999999999998 %
4 枚 19260 回、その確率は 19.259999999999998 %
5 枚 17001 回、その確率は 17.000999999999998 %
6 枚 12158 回、その確率は 12.158 %
7 枚 7327 回、その確率は 7.327 %
8 枚 3913 回、その確率は 3.913 %
9 枚 1912 回、その確率は 1.9120000000000001 %
10 枚 785 回、その確率は 0.7849999999999999 %
11 枚 299 回、その確率は 0.299 %
12 枚 102 回、その確率は 0.10200000000000001 %
13 枚 32 回、その確率は 0.032 %
14 枚 12 回、その確率は 0.012 %
15 枚 4 回、その確率は 0.004 %

現在インストール直後のボーナスで36パック引けるようだが、それで7~9枚引ければ成功といえるのではないか

ところでこれって典型的なベルヌーイ試行である
というわけでnumpyのbinomial関数を使ってみる

import matplotlib.pyplot as plt
import numpy as np
import collections as col

#numpy.random.binomialを使った場合
p = 1.5/100 #レジェンド
pack = 36
time = pack*8
x = 10000000 #試行回数
leg = np.random.binomial(time,p,x) #リストとなる
legdict = col.Counter(leg)

plt.hist(leg,bins=max(leg)) #なぜかデフォでbins=10

print(pack,"パックを",x,"回リセマラした場合に引くレジェンド(確率",p,")の枚数")
for key ,value in sorted(legdict.items()):
    print(key,"枚",value,"回、その確率は",(value/x)*100,"%")

f:id:moppii:20161104221810p:plain
36 パックを 10000000 回リセマラした場合に引くレジェンド(確率 0.015 )の枚数
0 枚 128278 回、その確率は 1.28278 %
1 枚 564045 回、その確率は 5.64045 %
2 枚 1234020 回、その確率は 12.3402 %
3 枚 1792282 回、その確率は 17.92282 %
4 枚 1943714 回、その確率は 19.43714 %
5 枚 1680494 回、その確率は 16.80494 %
6 枚 1207038 回、その確率は 12.07038 %
7 枚 739904 回、その確率は 7.399039999999999 %
8 枚 396402 回、その確率は 3.96402 %
9 枚 187852 回、その確率は 1.8785199999999997 %
10 枚 79711 回、その確率は 0.79711 %
11 枚 30445 回、その確率は 0.30445 %
12 枚 10814 回、その確率は 0.10814 %
13 枚 3543 回、その確率は 0.035429999999999996 %
14 枚 1074 回、その確率は 0.01074 %
15 枚 287 回、その確率は 0.00287 %
16 枚 73 回、その確率は 0.0007300000000000001 %
17 枚 17 回、その確率は 0.00017 %
18 枚 5 回、その確率は 4.9999999999999996e-05 %
19 枚 1 回、その確率は 9.999999999999999e-06 %
20 枚 1 回、その確率は 9.999999999999999e-06 %

最初かいたコードとおよそ似たような結果になったので満足
ところで標準ライブラリのrandomを使うと100万回試行で軽食が取れる程度には時間食うのに、ふたつめのコードで使ったnumpyのやつだと1000万回回してもものの数秒で終わるのはどういう理由からなのだろう

kiwami high low攻略

kiwami high lowというフリゲがあって、まあ普通のハイロウゲームなんだけど、ちゃんと遊んだことなかったのでやってみるとコレが面白く、熱中してしまった、30分くらい
毎ゲームごとにプレイヤーがミスをすると、女の子(きわみちゃん?)相手にロシアンルーレットが始まり、運良く空発だと女の子が恐怖で失禁するのだが、この演出がなかなか良くて、酒のみながらプレイしてたのもあってとても興奮する そして運が悪いと普通に死ぬ
以下は攻略というほどでもないが、PCに眠らせておくのもなんかアレなので、ここに置いておく
 
---
そもそものルール、ハイローに成功するたびに所持金が2倍になっていき最終的に1億円を達成していれば勝利
ドロー又はパスの場合所持金は不変
また20億円を達成したら別endとなる
数字を外すとロシアンルーレットが始まり確率で死ぬ
ハイローを終えても1億に満たないとbad end
ハイローの代わりにセイムを成功させると所持金が10倍になる
 
35ターンあるのでそれまでに1億円を達成しないといけない
成功するたび金額は2倍になるので、逆に言えば1億達成するために必要な成功回数は対数で計算できる
log2(1億円)≒26.57 :つまりセイムなしで1億円達成するためには27回正解する必要があり、ミスとパスを計8回まで許容できる
ミスの確率を考えれば現実的に3パスまでか
セイムなしで20億達成するなら? log2(20億)≒30.89なので31回 は正解しないといけない つまり4ミスまで、厳しい?
5ミスで確定射殺
 
セイムはいつ出してもいい、なら初手セイムを成功させるリセマラしたほうが効率いい
一発目にセイムになる確率は3/35≒0.0857つまり8.6%
2連セイムは無理...(連続セイムは数字が続くのでセイム正解セイムを狙うとして、確かに金額達成という意味では楽になるが確率的にはかなり厳しい、リセマラするにしても時間の無駄
単純計算で2連セイムが0.86^2=0.73% はい無理)
3正解より1セイムのほうが価値あり、つまり初手セイムを決めればわりと休める
 
リセマラして初手セイムを成功させるとして、それぞれの成功条件は
log2(1億円/10)≒23.25
log2(20億円/10)≒27.57
つまり1億達成したいなら初手セイム成功のリセマラさえすればパスもしくはハイローの失敗を計10回まで、
20億達成したいなら、初手セイム成功のリセマラさえすればパスもしくはハイローの失敗を計6回まで繰り返せる
20億達成したいなら、5を全部パスなんて悠長なことはやってられない
 
結論→初回パスリセマラしろ、成功したらミスとパスの数を数えつつ後は運否天賦

f:id:moppii:20161014212953p:plain

f:id:moppii:20161014212957p:plain