Blogs

「Go言語でつくるインタプリタ」を読み終えました。

この本は、私のようなGo言語初心者のでも理解しやすい構成で、プログラミング言語の構築に必要な字句解析や構文解析の基本的な概念を学びながら、動くインタプリタを完成させることができます。 プログラミング経験があるものの「字句解析」や「構文木」という用語はなんとなく聞いたことがある程度だった私にとって、実際に手を動かしてプログラムを動作させることができ、学びと共に達成感もひとしおの一冊でした。

多くのハンズオン本では、「インタプリタの一部だけ」を実装することに留まりがちで、最終的に完成するのが部分的なコードであったり、動作しない例も多くあります。 しかし、この本は最初から最後まで実際に動くインタプリタを作り上げることができ、その後にいろいろ試してみたりもできます。 実際私はREPLだけでなく、ファイルを読み込んで実行させるように改造して楽しみました。 このような「動くものを作り上げる」体験は、プログラミングやシステム理解にとって重要なポイントだと思います。

具体的には、以下のようなセクションごとの段階的なステップを経て、コード全体が構築されていきます。

字句解析: 「トークン」や「字句解析器（Lexer）」という基本的な概念から始まり、最初のコードの実装を通してデータがどのように分解されるのかが明確になります。

構文解析: 「構文解析器（Parser）」を使ってAST（抽象構文木）を生成することで、プログラムの文法構造がどのように理解されるかを深く掘り下げます。

評価器: 最後に、構築したASTを評価し、インタプリタとして動作する一連のフローを構築していきます。

Go言語初心者でも学びやすい理由の一つに、各章でのコードと概念の説明が一貫して丁寧に行われている点が挙げられます。 例えば、字句解析の章では一つひとつのトークンの分解と構造を詳しく説明しており、構文木（AST）の組み立ての理解が進むような工夫がされています。 構文解析や評価に進むに従い、少しずつ難易度が上がりますが、各章で順を追って理解していけるため、初心者にとっても安心です。

また、テスト駆動開発であることも理解が非常に進みました。 私は普段コードをあまり書いていないので、テスト駆動開発の恩恵に預かれていなかったのですが、先にどういうアウトプットが想定されるかが提示されるのは理解の助けになりました。

時間はかかりましたがなんとか読破できましたので、個人的なおすすめの読み方を提示しておきます。

コードを写経する

手を動かしてコードを書くことは、多くの人にとってこの本を読む上では必須だと思います。手稲にコードが記載されていますが、それでもぼんやり写経してると飛ばしたり書くところを間違えたりして、それを解消するプロセスがまた理解につながったりします。 電子書籍だとコードをコピー&ペーストできてしまうこともありますが、ここは実際に自分の手で入力して試行錯誤するのがおすすめです。

理解のズレに注意

インタプリタの実装は、エラーが起こりやすい分野でもあります。 例えば、字句解析のセクションではトークンの扱いを間違えやすく、構文解析ではASTの構造の間違いが結果に影響します。 実行時にエラーが発生した場合、なぜそのエラーが起きたかを丁寧に追っていくことが大切です。 その点、ユニットテストも提示してくれているところは非常にありがたいです。ユニットテストを頼りに自分の理解を補強することができます。

各章ごとに振り返る

一気に読み進めず、各章で実装した機能がどのように働くのかを自分なりに理解する時間を取ることをお勧めします。各章が次の章の基礎になるため、途中で理解が浅い部分があると後の章でつまづきやすくなります。

総評

「Go言語でつくるインタプリタ」は、Go初心者から中級者に向けた実践的な内容で、言語構築の基礎を学びたい方に非常におすすめの一冊です。 特に、プログラムが動作し、抽象的な概念を実際の動くコードとして実現できる達成感を味わうことができます。 字句解析や構文木のような「聞いたことはあるけれど、実はよく知らない」という分野に一歩踏み込めて、さらにGo言語の勉強にもなる本書は、自分のレベル感にはちょうど良い一冊でした。 ちなみにジェネリクスは無い時代の本なので、その辺を改良するのも良い勉強になるかもですね。

また、本書には続編があります。 https://compilerbook.com

今度はインタープリターではなくコンパイルして仮想マシン上で動作させるそうです。 こちらはまだ日本語訳されていませんが、時間を見つけて取り組んでみたいと思います。

2年間履修した帝京大学の科目履修生コースを終えました。今年も引き続き科目履修生として登録するか、正科生への切り替えを検討していましたが、家庭の事情やら仕事の忙しさもあり、今年はこれまで学んできたことの復習を隙間時間に独学でやっていこうという事にしました。

一区切りということで2年間の科目履修生を振り返ってみたいと思います。

履修科目

2020年度履修科目

• 論理数学
• 数理統計学
• 離散数学
• オートマトンと計算理論
• コンピュータネットワーク
• データ構造とアルゴリズム

平日の仕事終わりや土日の多くの時間は授業を受けたりレポート書いたりしてました。仕事は比較的早めに上がれたのと、家族の協力あってのことでした。

おかげさまで全部単位は取れましたが、統計学が痛恨のC評価で、一番勉強した科目だったのですがセンスのなさが露呈しました。

年度の前半でほぼ試験も終わったので、詰め込めばもっと受けられたと思います。とはいえ1年に10〜15科目くらいが限界かなと思います。ちなみに暇だった年度の後半は勉強熱が高まっていたので、AWSの資格試験をバンバン受けてました。

2021年度履修科目

• オペレーションズリサーチ
• コンピュータシミュレーション

2年目は2教科履修しました。2年目に教科数を減らしたのは仕事が忙しくなりそうと予想したからで、案の定忙しくなったので2教科にしておいてよかったです。この年も年度の前半でほぼ授業は終わりました。後半は仕事に追われてました。

授業の感想

レポートも細かいところまでしっかりみてもらえ、赤ペンでいっぱい書き込まれたフィードバックをもらった時は、間違えた気恥ずかしさはありつつもやってよかったなと思えました。 メールで何度か質問させてもらいましたが、こちらも丁寧に回答いただきました。

テストはコロナ禍と言うこともあり全てオンラインで受けました。なので本当にリアルでの大学生活は皆無でした。通信大学とはいえテストくらいは大学生気分を味わえると思っていたのですが、まぁ仕方がないですね。

お金に見合った価値があったかと言われると、自分の場合はあったと言い切れます。 ただ、大学は通えば良いというものではないことは、義務教育・高校・大学と進んだ若かりし頃を思えば想像ついていましたし、多くの方も実感としてあるのではないでしょうか。そこでどのように振る舞うかが重要です。単位はもちろん取れたら嬉しいですが、単位を取るためだけを目的にしていたら、あまり価値はないものになるでしょう。

仕事に好影響を与えるかと言われると微妙です。SREとして従事している今の仕事に直接影響はなかったかなと思います。そもそも自分の場合は仕事に生かすというよりかは、今従事しているコンピュータというお仕事の根幹をなす学問領域に興味があったというのが動機でした。

ただ、オンプレの経験がある人がパブリッククラウド上でのインフラ構築時に勘所に優れていたり、低レベルなプログラミング言語を操っていた人が高水準のプログラミング言語でも優れた能力を発揮されているケースを何度かみてきました。コンピュータサイエンスのような低レベルな知識が間接的にコンピュータに携わる上で好循環を与えるのかなとも考えています。

これから

冒頭にも書きましたが、今年は独学でコンピュータサイエンスを学んだりしながら授業で習ったことを今一度復習したいなと考えています。本当にたくさんの事を知らないんだなと痛感することができました。仕事やプライベートに余裕が出たら改めて門をたたきたいと思います。

ここ数年HHKBをずっと使っていて、非常に気持ちよく利用しています。巷の噂に違わず入力するのが気持ちよくなるキーボードです。

ヘビーユーザーの証の一つである尊師スタイルの存在は以前から知っていて憧れがありました。ただ自分はHHKBを机に置いたノーマルなスタイルで利用することが多かったため、踏み込んでいませんでした。

そんな自分が尊師スタイルに入門したきっかけや感じたメリット・デメリットを記しておこうと思います。これからHHKBを購入するか迷っている人、尊師スタイルに憧れている人の参考になればと思います。

きっかけは在宅続きの運動不足からスタンディングデスクでPCを利用したくなったところから始まりました。私の場合は予算の関係上Moft Zを利用したスタンディングデスクを導入しました。1万円以内でスタンディングデスクが導入できるので非常にお得です。

Moft Zを利用する場合、MacとHHKBを両方置くスペースは流石に無く、キーボードはどうしてもMacのものを利用することになります。キーボードの打鍵感もさることながらバッククォートの位置がMacのキーボードとHHKBで異なることにストレスを感じるようになりました。最近はドキュメントを書く時やSlackでコミュニケーションするときもインラインコード入力時にバッククォートを使うことが多いです。

スタンディングデスク利用時のキーボード入力にストレスがあると、健康のためというだけではだんだんやらなくなるのは目に見えていました。せっかく買ったMoft Zも無駄になり、座りっぱなしの日常では不健康まっしぐらです。

そこで思いついたのが尊師スタイルです。尊師スタイルであればスタンデイングデスクの時もHHKBをMacの上に置けるのでストレスから解放されます。私のような理由以外でもただ単にカッコ良さそうとか、デスクのレイアウトの関係上尊師スタイルが向いている方もいるかと思います。自分なりに尊師スタイルのよかった点と悪かった点をまとめてみました。

メリット

一番の目的としていた入力インターフェースがHHKBに統一される点は、想像以上に快適でした。効率化は計り知れないと感じます。おかげで当初の目的であるスタンディングデスクでの作業も快適です。

デメリット

TouchIDが使えない

Macの場合Touch IDがキーボードブリッジに覆われてしまいます。PCログイン時やサイトのパスワード入力時に重宝していたTouch IDが使えなくなるのは思ったより痛手でした。Touch ID付きHHKBが出たら即ポチってしまいます。が、今は大人しくパスワードを入力するしかないですね。Apple Watchを持っている方ならMacへのログインはApple Watchが補助してくれます。

持ち運びが不便

これでどこでも尊師スタイルでHHKBと一緒と思いきや、HHKB、キーボードブリッジ、パームレストの3点セットを持ち歩くのはやはり面倒でした。コーディングならいざ知らず、ちょっとした作業や買い物をする際はTouchIDが使えないことも相まって面倒さが増します。 いい感じのキャリングケースがあれば緩和されるかもしれません。

尊師スタイルに必要なもの

尊師スタイルに必要なものも挙げておきます。

• キーボードブリッジ
• セパレート式パームレスト

キーボードブリッジが無くても置き方によってはいけそうな気もしましたが、注意深くキーボードを置くのは面倒なのでお金で解決しました。また、MacのTouchPadを使うのであればセパレート式のパームレストも必要です。

まとめ

ちょっぴり憧れだった尊師スタイルに入門できただけでも良い経験でした。この調子でHHKB沼へどんどんはまっていこうと思います

1年ほど積読にしてた「入門監視」を読み終えました。

この本は監視ツール個別の話はほとんどなく、監視をどう設計すべきかが記されています。 私自身、監視を色々見たり設定したりしてきましたが、考え方の答え合わせができて楽しかったです。

まずはアンチパターンをチェックする

この本を買う前に1章だけでもみてもらって、自分がアンチパターンにハマっていないかを是非チェックして欲しいと思います。

流行りのツールにOSメトリクスでとりあえずアラート仕掛けたりしていた自分には胸をえぐられる章でした。 心当たりのある人は目次を見ただけでもグッとくると思います。

 1.1　アンチパターン1：ツール依存
    1.1.1　監視とは複雑な問題をひとくくりにしたもの
    1.1.2　カーゴ・カルトなツールを避ける
    1.1.3　自分でツールを作らなければならない時もある
    1.1.4　「一目で分かる」は迷信
1.2　アンチパターン2：役割としての監視
1.3　アンチパターン3：チェックボックス監視
    1.3.1　「動いている」とはどういう意味か。「動いている」かどうかを監視しよう
    1.3.2　アラートに関しては、OSのメトリクスはあまり意味がない
    1.3.3　メトリクスをもっと高頻度で取得しよう
1.4　アンチパターン4：監視を支えにする
1.5　アンチパターン5：手動設定

自分の経験を答え合わせをする

自分はこれまでWebシステムを中心に監視の設定を見たり行ったりしてきました。 世の中には色々なシステムがあり、色々な事情もあるとは思いますが、これまで経験してきた監視の経験を本書で答え合わせをしてみました。

所属組織の標準的な監視パッケージに乗っかったこともあれば、自分で1から設計したこともあります。

どこかで見聞きしたものなのかも忘れましたが、自分で監視を設計するときは下記のポリシーを大切にしてきました。

• 利用者より先に気づく事
• ユーザーの実際の不便を監視する事

「動いている」を監視する

サービスが元気に動き続けることがベストではありますが落ちることも必ずあるのがシステムの世界です。そんな時にサービスの不調にいち早く気づくことが監視の一番の役目であると考えています。

なのでURL監視を中心とした監視を組みつつ、可能であればビジネス的に一番熱いユースケースやユーザーの不便に最もつながりそうなユースケース、 つまりユーザーの信頼を失うことにつながりそうなユースケースを実際にトレースするような監視を作ってきました。

また、この監視が動いていればとりあえずシステムとしては健全であるというラインが担保できていると、システム運用にかかるストレスが非常に低減されます。

例えばポイント系のシステムでポイント付与までのユーザーの行動をトレースするようなものを、ユーザー登録が命のようなサービスであれば実際にユーザー登録を行うようなものを外部から実際に自動的に行うような監視を行ったりしました。

なので本書でいう「動いている」を監視できていたと思います。

お手製の監視システムでOSメトリクスを監視

一方で、OSのメトリクスに関してはついついCPU 80%とかロードアベレージ10とかを閾値とした監視を手癖のように設定しがちではありました。振り返ってみてもCPUやロードアベレージなどのメトリクスによるアラートにより実際に行動を起こす事は少なく、 URL監視などの外形監視の補助として使うことが多いです。

そして監視ツールの選定ですが、監視Saasが成熟してきている今日、未だに本格導入を試みたことがありません。AWSだとCloud Watchをベースにした監視などを行ったことはありますが、当時はカスタマイズ性にいまいち満足がいかなかったという記憶があります。 そういった古い記憶から無意識に監視Saasそのものを敬遠しがちになっていたのは否めません。

直近ではPrometheusを導入して、それはそれで満足はしていますが、同じことを監視Saasでやれなかったとはいえません。監視Saasも時が経つにつれ洗練されているはずで、コスト対効果を鑑みてもっと積極的に導入を考えても良かったと反省する部分もあります。

この辺りは本書でいうアンチパターンにハマっていたと思います。

まとめ

自分の経験を本書と照らし合わせて答え合わせのできる楽しい読書でした。

監視は作り手の自由度が高い分野であると思います。下手をするとアラートだらけになってしまう所を本書の指針に沿って設計する事で、みんなが幸せになれる監視を設計できるかと思います

Pythonプログラミングイントロダクション(20章) まで読みました。 20章はベイズ統計についです。

スパムフィルターでおなじみベイズフィルターですが、ベイズ統計を理解する上で重要な条件付き確率に関して問題が出ていました。

本書ではベイズ統計に関する例を取り上げつつ、Pythonで実装を提示してくれています。 そこは本書をみていただくとして、条件付き確率の問題について取り上げたいと思います。

条件付き確率

本書掲載の問題は下記になります。

あなたは森の中を歩き回り,美味しそうなキノコの群生を見つけた.
カゴいっぱいにキノコを採り,家に帰って旦那様のために,キノコ料理の準備を始める.
しかし,調理を始める前に,彼はあなたに,野生のキノコが毒キノコかどうかを本で調べるようにお願いした.
本によると,おしなべて野生キノコの80%は毒キノコとのことだ.
あなたは手元のキノコを本の中の写真と見比べ,95%の割合でそれが安心して食べられるキノコと判断した.
あなたはどの程度安心して,キノコ料理を旦那様に作ってあげられるだろうか.

採ってきたキノコが毒キノコの確率は80%になります。なので毒キノコで無い確率は20%です。これをP(B)とします。

写真をみて95%の確率で毒キノコで無いと判断していますが、実際採ってきたキノコが毒キノコの場合とそうで無い場合で、写真から判定される確率は異なると思います。 つまり写真と見比べて毒キノコで無いと判断した際の確率P(A)は毒キノコを採ってきたかそうで無いかという確率P(B)と独立でないと考えます。

仮に独立だとしても P(A) = P(A|B)です。

仮に、とってきたキノコが毒キノコで無いという条件のもと、写真と見比べた結果毒キノコで無い確率を95%とした場合、 採ってきたキノコが毒キノコでなく、写真で判断しても毒キノコで無い確率はいくらになるでしょうか。

写真判定シロ かつ 採ってきたキノコは毒でないを求めたいわけです。式で表すと

P(A)とかP(B)を使うと

つまり19%の確率で毒キノコではないということになります。

あなたは食べますか？

僕は食べません。

まとめ

条件付き確立の問題について取り上げてみました。 ここまでの章で学んできた頻度統計学とは一線を画するベイズ統計について興味を持たせてくれる章でした。

Pythonプログラミングイントロダクション(19章) まで読みました。 19章では統計的仮説検定についてです。

自転車競技の記録を伸ばす薬であるPED-XとPED-Yのどちらが優れているか無作為試験を行った結果から、得られた平均タイムについて検討していきます。 本当に平均タイムに有意な差があるのか、偶然発生した差なのかを検定をする方法が説明されています。

フィッシャーの統計的仮説検定

本書で述べられているフィッシャーの統計的仮説検定の手順は下記になります。この手順により観測結果が偶然起こったのか確率を評価していきます。

1. 帰無仮説(null hypothesis)と対立仮説(alternative hypothesis)を立てる.
2. 評価する標本についての統計的仮定を理解する.
3. 適切な検定統計量(test statistic)を計算する.
4. 帰無仮説における検定統計量の確率を得る.
5. その確率が、帰無仮説が偽であると推測できるほどに十分小さいかどうか,すなわち,帰無仮説を棄却(reject)するかどうかを決定する.

1. 帰無仮説(null hypothesis)と対立仮説(alternative hypothesis)を立てる.

PED-XとPED-Yの平均値を例にとった場合、下記のようになります。

帰無仮説：PED-XとPED-Yの平均の差が0

対立仮説：PED-XとPED-Yの平均の差が0でない

2. 評価する標本についての統計的仮定を理解する.

PED-XとPED-Y使用者のフィニッシュタイムは正規分布となり、標本は無限母集団から抽出した物でした。

3. 適切な検定統計量(test statistic)を計算する.

t値を計算します。t値とは二つの平均の誤差が標準誤差で見て0からどれくらい差があるかを表した値です。 本書では-2.13165598142となりました。

標準誤差は17章で出てきたものです。

\( \sigma_{a}^{2} : 母集団の分散 \\\ \sigma_{a} : 母集団の標準偏差 \\\ n : 標本サイズ \)

\[ SE = \frac{\sigma}{\sqrt{n}} \]

4. 帰無仮説における検定統計量の確率を得る.

p値を計算します。p値とは帰無仮説が成立する前提で、統計量が観測された値以上の極端な値が得られる確率のことです。 つまり、PED-XとPED-Yがの平均の差が0であったときに差が現れるような極端な値が得られる確率のことになります。 本書では0.0343720799815、つまり3.4%程になりました。

5. その確率が、帰無仮説が偽であると推測できるほどに十分小さいかどうか,すなわち,帰無仮説を棄却(reject)するかどうかを決定する.

本書ではp値が3.4%となり、PED-Xの方が優れていそうと言う結論になりかけますが、この後p値が信用ならないと言う話につながっていきます。

p値を何度も計算してみる

本書を読み進めていくと、この例はプログラムにより作り出されたデータであると言うたねあかしがされます。 平均値119.5、標準偏差5.0の分布と平均値120.0、標準偏差が4.0の分布からランダムにサンプリングした値が使われていました。 random.gauss関数でデータを作っています。

実際に何度か下記のようなプログラムを流してみると、p値の値は乱高下します。

 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

import random
from scipy import stats


def get_t_value():
    treatmenDist = (119.5, 5.0)
    controlDist = (120, 4.0)
    sampleSize = 1000

    treatmentTimes, controlTimes = [], []

    for s in range(sampleSize):
        treatmentTimes.append(random.gauss(treatmenDist[0], treatmenDist[1]))
        controlTimes.append(random.gauss(controlDist[0], controlDist[1]))

    controleMeans = sum(controlTimes) / len(controlTimes)
    treatMeans = sum(treatmentTimes) / len(treatmentTimes)

    twoSampleTest = stats.ttest_ind(treatmentTimes, controlTimes, equal_var=False)
    print('Treatmeant mean - control mean = ',
          treatMeans - controleMeans, ' minutes')
    print('The t-statistic from two-sample test is ', twoSampleTest[0])
    print('The p-value from two-sample test is ', twoSampleTest[1])


for i in range(10):
    twoSampleTest = get_t_value()

結果は下記のようになります。p値が1%台の時もあれば60%台の時もあります。

Treatmeant mean - control mean =  -0.6843325129559474  minutes
The t-statistic from two-sample test is  -1.0716961219246917
The p-value from two-sample test is  0.285226164284344
Treatmeant mean - control mean =  -0.39645994742086543  minutes
The t-statistic from two-sample test is  -0.6751950331252392
The p-value from two-sample test is  0.5004085432884583
Treatmeant mean - control mean =  0.25269386654427706  minutes
The t-statistic from two-sample test is  0.3950053957541872
The p-value from two-sample test is  0.6932658871131978
Treatmeant mean - control mean =  0.50100925426743  minutes
The t-statistic from two-sample test is  0.7660580410192988
The p-value from two-sample test is  0.4445624847162438
Treatmeant mean - control mean =  -0.9383152222142996  minutes
The t-statistic from two-sample test is  -1.3925614858916715
The p-value from two-sample test is  0.16555904568337207
Treatmeant mean - control mean =  -0.4504161747373985  minutes
The t-statistic from two-sample test is  -0.7134525526996994
The p-value from two-sample test is  0.4764536661236095
Treatmeant mean - control mean =  -0.517887936147531  minutes
The t-statistic from two-sample test is  -0.818800036206608
The p-value from two-sample test is  0.41390883875821693
Treatmeant mean - control mean =  -0.7535810558372305  minutes
The t-statistic from two-sample test is  -1.2462075073793941
The p-value from two-sample test is  0.21417712992121332
Treatmeant mean - control mean =  -1.5033637291589912  minutes
The t-statistic from two-sample test is  -2.347673847276862
The p-value from two-sample test is  0.019992408614120653
Treatmeant mean - control mean =  -0.3150896341041829  minutes
The t-statistic from two-sample test is  -0.5046312441817462
The p-value from two-sample test is  0.6144229896604907

と言うことで、200程度のサンプリング数だとp値がこれだけ上下することが示されました。

まとめ

仮説検定についての流れを追ってみました。取り上げた流れは2標本の両側検定と言うものです。 本書では他にも片側検定と1標本検定についての話題も取り扱っています。

本書ではp値が低いのは帰無仮説が実際に低いかもしれないし、母集団の代表的な例でない場合もあると言われています。 そして、このp値の問題を解決する考え方として20章のベイズ統計に続いていくわけです。

Pythonプログラミングイントロダクション(17章) まで読みました。 17章では中心極限定理について説明されています。

この章の目標の一つは「一つの標本から標本平均の信頼度を出す」ことにあります。 そのための道具として中心極限定理を利用します。

中心極限定理

中心極限定理についての引用は下記になります。

- 同じ母集団から抽出された標本サイズが十分に大きいと、標本の平均(標本平均)はおおよそ正規分布に従う
- この正規分布の平均は母集団の平均に非常に近い。
- 標本平均の分散は母集団の分散を標本サイズで割ったものに非常に近い

中心極限定理により標本平均は正規分布に従います。

また、標本平均の分散は母集団の分散を標本サイズで割ったものに非常に近いという事は、標本平均の標準偏差は母集団の標準偏差を標本サイズの平方根で割ったものに非常に近いということになります。 これを式にすると、標本平均の標準偏差を求める式は下記のようになります。

母集団の標準偏差を標本サイズの平方根で割ったものは標準誤差とも言われます。 そしてこのとき必要な母集団の標準偏差は、ある一定数の標本サイズがあれば標本の標準偏差を用いて良いとされているそうです。

まとめると、下記の手順で中心極限定理が正しいかを見ていきます。

1. 標本の値から標準偏差を求め、母集団の標準偏差とする。
2. 標本平均の標準偏差を求める。
3. この標準偏差から95%信頼区間を導きます。
4. 実際の母集団の平均値と標本平均の誤差を求める

アニメのレーテイングデータで検証する

kaggleのサイトでアニメのレーティングに関するデータを利用させていただきました。 本書と同様の手順でこのデータから標本サイズ200のデータを元に信頼区間を導き、 実際の母集団から求めた平均との絶対誤差を元に信頼区間の値が妥当を見ていきます。

 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

import random
import pylab
import csv

def variance(X):
    """分散"""
    mean = sum(X)/len(X)
    tot = 0.0
    for x in X:
        tot += (x - mean)**2
    return tot/len(X)


def stdDev(X):
    """標準偏差(分散の平方根)"""
    return variance(X)**0.5

def read_anime_score(path):
    """データを読み込みスコアのリストを返す"""
    scores = []
    with open(path, newline='') as f:
        reader = csv.reader(f)
        # skip header
        next(reader)
        for row in reader:
            scores.append(float(row[15]))

    return scores

data = read_anime_score('./anime_cleaned.csv')

# 母平均
popMean = sum(data) / len(data)
# 標本サイズ
sampleSize = 200
numBad = 0
for t in range(10000):
    sample = random.sample(data, sampleSize)
    # 標本平均
    sampleMean = sum(sample) / sampleSize
    # 標準誤差
    se = stdDev(sample) / sampleSize**0.5

    # 母平均と標本平均の誤差から信頼区間外のデータ数を求める
    if abs(popMean - sampleMean) > 1.96 * se:
        numBad += 1
print('Fraction outside 95% confidence interval =', numBad / 10000)

: Fraction outside 95% confidence interval = 0.05

概ね5%程度の結果が求められます。 中心極限定理で定義されている通り、標本サイズが十分大きいと標本平均が正規分布に従い、母集団の平均に非常に近いと言えそうです。

まとめ

本章は話の流れがわかりにくかったので最初は何を言っているのか正直わかりませんでしたが、 まとめてみると最終的に標本平均から信頼区間を求めたいのだとわかりました。

母集団の標準偏差が標本の標準偏差で代用できる点は気持ちの悪いところですが、本書でも下記のように説明されています。

多くの統計学者が言うには、母集団の分布がおおよそ正規分布であれば標本サイズが30から40で十分に大きい。
標本サイズが小さい場合にはt分布と呼ばれるものを使って区間計算したほうが良い。

t分布はこの後出てくるようなので一旦読み進めていくとします。

Pythonプログラミングイントロダクション(16章) まで読みました。 16章ではモンテカルロ・シミュレーションについての説明になります。

モンテカルロ・シミュレーションって大そうな名前がついていますが、この章までにやってきたランダムネスを利用したシミュレーションのことをモンテカルロ・シミュレーションといいます。 この章でもサイコロゲームを題材にランダムネスを利用したモンテカルロシミュレーションの実装がいくつか例示されています。

ランダムネスで円周率を求めるプログラミングなんかは痺れました。

ジャンケンのシミュレーション

ここではジャンケンのモンテカルロ・シミュレーションをPythonで実装して理解を深めていきたいと思います。

 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

import pylab
import random


def janken():
    JANKEN_PATTERN = ['G', 'C', 'P']
    while True:
        a_t = random.choice(JANKEN_PATTERN)
        b_t = random.choice(JANKEN_PATTERN)
        if a_t != b_t:
            break
    return janken_judge(a_t, b_t)


def janken_kai():
    JANKEN_WIN_PATTERN = ['G', 'G', 'C', 'P']
    JANKEN_LOSE_PATTERN = ['G', 'C', 'C', 'P']
    while True:
        a_t = random.choice(JANKEN_WIN_PATTERN)
        b_t = random.choice(JANKEN_LOSE_PATTERN)
        if a_t != b_t:
            break
    return janken_judge(a_t, b_t)


def janken_judge(a_t, b_t):
    """ a_tとb_tは違う値の想定。
    aが勝っていたらTrue, 負けていたらFalseを返す """

    win_pattern = [('G', 'C'), ('C', 'P'), ('P', 'G')]
    return (a_t, b_t) in win_pattern


def show_plot(trials, target):
    a_rate = []
    a_win = 0
    for k in range(1, trials + 1):
        if target():
            a_win += 1

        a_rate.append(a_win / k * 100)
    pylab.plot(a_rate)
    pylab.ylim(0, 100)
    pylab.xlim(1, trials)
    pylab.xlabel('times')
    pylab.ylabel('rate')
    pylab.title('Junken')
    pylab.show()


def get_stdDev(trials):
    a_rate = []
    for i in range(0, 100):
        a_win = 0
        for k in range(0, trials):
            if janken():
                a_win += 1
        a_rate.append(a_win / k * 100)
    return stdDev(a_rate)


def variance(X):
    """分散"""
    mean = sum(X)/len(X)
    tot = 0.0
    for x in X:
        tot += (x - mean)**2
    return tot/len(X)


def stdDev(X):
    """標準偏差(分散の平方根)"""
    return variance(X)**0.5


show_plot(1000, janken)
show_plot(1000, janken_kai)
print('10 trials stdDev: ', get_stdDev(10))
print('100 trials stdDev: ', get_stdDev(100))
print('1000 trials stdDev: ', get_stdDev(1000))

janken関数では、a,bそれぞれのグー(G)、チョキ(C)、パー(P)の中からランダムに選択した上で勝ち負けを判定しています。 あいこの場合はもう一度ランダムに選択しなおします。

完全にランダムにグー、チョキ、パーを選ぶ場合は試行回数が増えるについれて勝率は5割に近づいていきます。 15章で学んだ大数の法則の通りですね。

人間がジャンケンをする場合は完全にランダムとはいかず、グー、チョキ、パーのどれかに選択が偏ることがあるでしょう。 ランダムではありますが、aさんがグーを出しがち、bさんがチョキを出しがちな場合はaさんの勝率が上がっていきます。そんなjanken_kaiを作ってみました。 この場合はaさんの勝率が5割より上回っていますね。

続いて標準偏差を出してみました。 試行回数を増やすにつれ、標準偏差が減っていきます。

: 10 trials stdDev:  17.91612832955234
: 100 trials stdDev:  5.337158102067871
: 1000 trials stdDev:  1.6880153857382443

まとめ

これまで学んだ統計の知識を織り交ぜつつ、ジャンケンを題材にモンテカルロ・シミュレーションをやってみました。 ここ何章かは思ったより統計学の内容が多いことに若干戸惑いもありましたが、内容も楽しいですしトレンド的にも良い題材なのかと思います。 Pythonを学ぶついでに統計学も学べるお得な一冊だと思います。

Pythonプログラミングイントロダクション(15章)まで読みました。 15章では確率統計を題材に、正規分布や関連する話題として信頼区間について述べられています。

様々な確率統計の題材について、Pythonでの実装とpylabを用いたグラフ表示付きで説明されているので、写経をしていくことで確率統計の知識が深まります。 その中でも確率分布についてみていこうと思います。

様々な確率分布

この章では様々な確率分布を例にPythonでの実装が例示されています。 この章で学ぶ確率分布の一部を列挙してみました。

2項分布の練習問題

上記にあげた確率分布の中でも、2項分布については練習問題がついていたのでみていこうと思います。 練習問題は下記になります。

k回の公平なサイコロ投げでちょうど3が2回出る確率を計算する関数を実装せよ。
この関数を用いて、kを2から100まで変化させたときの確率をプロットせよ。

まず、サイコロを投げて3が出る確率は 1/6 になります。 k回サイコロを投げて2回しか3がでず、残りは3以外の数字が出るということになります。 k回サイコロを投げて3が2回出る組み合わせを考える必要があります。

k = 2 のときは

k = 3 のときは、1,2回目が3の時、1,3回目が3の時、2,3回目が3の時の3通りの確率を足し合わせます。

このような2項分布と呼ばれる確率分布の式は以下で表されます。

Pythonでの実装は下記のようになります。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

import math
import pylab

def binomical_distribution(k, p, t):
    c = math.factorial(k) / (math.factorial(t) * math.factorial(k - t))
    return c * p**t * (1 - p)**(k - t)


def show_plot(p, t):
    x, y = [], []
    for k in range(t, 101):
	x.append(k)
	y.append(binomical_distribution(k, p, t))

    pylab.plot(x)
    pylab.plot(y)
    pylab.ylim(0, 0.5)
    pylab.xlim(2, 100)
    pylab.xlabel('k')
    pylab.ylabel('p')
    pylab.title('Binomical Distribution')
    pylab.show()

show_plot(1 / 6, 2)

まとめ

2項分布の指練習について、実際にPythonで実装してみました。 この章で学んだ信頼区間についてはこの後でも深掘りされていくようなので楽しみです。

「統計学が最強の学問である(数学編)」を読み終わりました。機械学習を学ぶための数学の知識を叩き込んで、機械学習の入り口まで連れて行くぞという著者の意志を強く感じる本でした。

記号論理学から始まって、連立方程式、線形代数、微分積分まで、幅広い知識の中から機械学習に必要なエッセンスを詰め込みながら、機械学習・ディープラーニングの入り口まで連れて行ってくれる一冊になっています。

内容

この本の流れを簡単にまとめると下記のようになります。

1. 簡単な確率、記号論理学、集合から最終的には迷惑メールフィルタでも用いられていたナイーブベイズ
2. 1次、2次の関数を学び、ちょうど良い数を求めるための平方完成
3. 二項定理から二項分布、対数、ネイピア数について
4. たくさんの数を効率的に表記するためにΣ、ベクトル、行列
5. ちょうど良いところを探すための微分、確率密度関数を理解するための積分
6. 線形代数と微積分を同時に使うことで機械学習、ディープラーニングの入り口へ！

やはりクライマックスはこれまでの知識を積み上げて機械学習の入り口が見えてくる6章です。特に行列、偏微分のありがたみがわかった最小二乗法の解法はしびれました。

行列、偏微分を用いた最小二乗法とnumpy

2章では訪問回数と契約件数のデータを元に比例関係の式を最小二乗法で2章で学んだ最小二乗法は平方完成により答えをめました。 6章ではこれまで学んできた知識、つまり微分、行列を使った解法を解説してくれています。

下記が2章の例示されていたデータになります。

上記の表で表されるモデルが訪問回数をx, 契約件数をyとした場合に直線で表されるとします。 この場合にモデルを表す式は下記のようになります。

εが誤差でこの誤差を最小にするa, bを求めます。 データは3つありますので、それぞれx, yに代入すると3つの式が出てくるわけですが、これを行列で表現すると3つの式ではなく一つの行列式で表すことができます。

\[ y = X\beta + \varepsilon \]

X, y, β,εはそれぞれ下記のような行列になります。

εの二乗和が最小になるようなa, とbを求めるのが最小二乗法です。 行列の場合、誤差の二乗の和を求める方法としてはεの転置行列とεの内積で表現することができます。

\[ \varepsilon^{ \mathrm { T } }\varepsilon \]

この二乗和が最小になる時のa,bを求めるためにβで偏微分してやります。

\[ \frac{ \partial \varepsilon^{ \mathrm { T } }\varepsilon } { \partial \beta } \]

ここからは、この章にたどり着くまでに得た行列や偏微分の知識を駆使して進めて行くのですが、 詳しい解法は是非本書で確認してみてください。結果として下記の式が得られます。

\[ \beta = (X^{ \mathrm { T } }X)^{-1}X^{ \mathrm { T } }y \]

\[ (X^{ \mathrm { T } }X)\beta = X^{ \mathrm { T } }y \]

この式に今回の例を代入すると、下記のようになります。

ここまでくるとPythonで機械学習をやる上で定番のnumpyを使い、逆行列を用いて傾きと切片を求めることができます。 この本で学んだ数学の知識が、自分がこれまで学んできた機械学習の定番、numpyと出会った瞬間でした。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

import numpy

x = numpy.array([[1, 1], [1, 2], [1, 3]])
y = numpy.array([10, 40, 82])

left = numpy.dot(x.T , x)
right = numpy.dot(x.T,  y)
inv_left = numpy.linalg.inv(left)
beta = numpy.dot(inv_left, right)
print(beta)

結果は下記のようになります。

1

[-28.  36.]

データ(今回の場合x, y)が増えていっても、この式に当てはめれば容易に最小二乗法を用いたモデルを導くことができます。

まとめ

自分自身はこの本を読む前に機械学習の本を読み始めてしまっていましたが、これから学習する方には是非ともお勧めしたい一冊です。 この本で得た知識を携えて、なんとなく通り過ぎているこれまでの知識をもう一度ブラッシュアップしたいと思います。