SL7にrootを入れる

Scientific Linux7 にrootを入れてみる。

UbuntuやMacには以前入れたので以下を参照。

Ubuntuにrootを入れる - 物理の研究の備忘録

Mac High Sierraにrootを入れる - 物理の研究の備忘録

基本Macでやったのと同様にgitでやってみる。

最初の環境の設定以外は同じなので上のリンクも参照。

SL7での環境の設定

$sudo yum update

$sudo yum upgrade

で最新にする。

CERNの公式ページを見ながら必要なパッケージを入れる

Build Prerequisites | ROOT a Data analysis Framework

$sudo yum install git cmake gcc-c++ gcc binutils \
libX11-devel libXpm-devel libXft-devel libXext-devel

$sudo yum install gcc-gfortran openssl-devel pcre-devel \
mesa-libGL-devel mesa-libGLU-devel glew-devel ftgl-devel mysql-devel \
fftw-devel cfitsio-devel graphviz-devel \
avahi-compat-libdns_sd-devel libldap-dev python-devel \
libxml2-devel gsl-static

git

あとは同様にCERN公式のダウンロードページ

https://root.cern.ch/downloading-root

をみながら適当なバージョン選んで

 $git clone http://github.com/root-project/root.git

 $cd root
 $git checkout -b v6-13-02 v6-13-02

とかやった。

./configure

とうってみると

ROOT is built with CMake, see https://root.cern/building-root
Please run
mkdir obj; cd obj; cmake ..; make -j 4

とメッセージが出た。

言われた通りにobjという名前のディレクトリを作って

$cmake ..

$make -j 4

でうまくいった。

あとは.bashrcに

source ~/root/obj/bin/thisroot.sh

の１行を加れば毎回

$root

でrootが起動するようになる。

cmakeが古いと怒られる場合

ここで入れたcmakeがgitで./configureするときに最新じゃないよと後々怒られることがあるので、その場合以下のサイトをみながら最新のcmakeを入れる。

[Linux][cmake] 最新のcmakeをインストールする方法

↑最後のビルドの部分はsudoを入れないとこけたのでsudo必要っぽい。

ビルドが終わった後にためしに

$cmake

と打ってみるとpathが通ってないよと怒られた。

install_manifest.txtなるものを見てみると、/usr/local/以下にpathを通すと使えるようになった。具体的には~/.bashrcに

export PATH=$PATH:/usr/local

の一文を加えて

$source ~/.bashrc

とした。

TVをデュアルディスプレイでPCに繋いだ際に画面の端がはみ出る問題の対処法

HDMIケーブル等でテレビをPCに繋げば２画面で快適に作業することができる。ただ、テレビ側の画面の端がはみ出てしまうことがある。

原因はTV側の設定

テレビにはデフォルトでオーバースキャンという機能がonになっていることが多い。

これはテレビを受信する際に端の部分を見えないように少し拡大して表示する機能で、アナログテレビ時代の名残らしい。（昔は画面端の部分にノイズが乗ったりすることがあったらしく、拡大表示するのが普通だった）

テレビの設定でこのオーバースキャンをoffにすれば直る。

ちなみに、普通にテレビを見る場合でも、特に理由がなければオーバースキャンをoffにしておくといいかもしれない。

意外に知られてないみたいだねー

おわり

Mac High Sierraにrootを入れる

ここでいうrootはcernの開発してる解析ライブラリをさします。

最近SierraからHigh SierraにOSをアップグレードしたら今まで動いていたrootが（案の定）動かなくなってしまった。

多分pathとかの設定がいろいろ変わっちゃってるんだと思うんだけど、error message見てもよくわからなかったので、最初からインストールし直すことにした。

rootを入れる方法は大きく分けて3つあって

①自分でsouce持ってくる

②バイナリ形式

③git

の3つ。

今回はgitで。

インストールする前に以下のlink見ながらOSごとに必要なものを先にインストールしておく。

自分の場合はmax OS X

https://root.cern.ch/build-prerequisites

1.xcodeとか

https://root.cern.ch/build-prerequisites

xcode -select —install

でなぜかフリーズしてしまったので以下のlinkから手動で探す。

command line toolsで検索して自分のmacに入れる。

https://developer.apple.com/download/more/

2.Xquartz

次に以下のlinkからxquartzを入れる。こっちはすんなりいった。

https://www.xquartz.org/

3.その他

$brew install gcc

$brew install cmake

でいろいろ入れる。

ここで入れ忘れたのがあると後でmakeしたときにこける。

自分はcmakeは入れずに進んじゃって後から怒られた。

次にcern公式ページからrootをダウンロードする。

最新版よりいくつか古い奴を適当に選んで入れる。

https://root.cern.ch/content/release-61204

ターミナルで

$git clone http://github.com/root-project/root.git

と打つ。

$cd root

で移動した後

$git checkout -b v6-12-04 v6-12-04

$./configure

$make

と順番に打つ。ちょっと時間がかかるので気長に待つ。

これでインストールは完了。最後にpathを通す。

~/.bashrcに

export ROOTSYS=/Users/usrname/root

source $ROOTSYS/bin/thisroot.sh

と加える。usrnameのとこは自分のuser name

今回はhome directory 以下のrootに入れたのでそこにpathを通している。

ターミナルを開いたときに毎回thisroot.shをsourceするようにする。

これでターミナルを開きなおしてrootって打ったときにrootが起動すれば成功。

ロゴが毎回でるのはうざいので

alias root='root -l'

とかも加えとくと良い。

Linux環境でPython、機械学習を使って手書き文字を識別する

基本的には以下のサイトの通りに（言われるがままに）やっみた。

【機械学習】Python3 + scikit-learn で識別率99%の手書き数字の分類器を作った - 株式会社クイックのWebサービス開発blog

自分の環境はUbuntu 16.04 LTS

環境構築

1.apt-get使ってぱーっと入れる。

sudo apt-get install python3 python-dev python3-pip
sudo apt-get install virtualenv

sudo apt-get install python-pip

2.
mkdir image_classify
cd image_classify
virtualenv -p python3 env
source env/bin/activate

3.

教師データは↓のものを使うよ

MNIST handwritten digit database, Yann LeCun, Corinna Cortes and Chris Burges

# sudoで実行しない
pip install numpy
pip install scipy

# 先にライブラリを入れておく
sudo apt-get install gfortran liblapack-dev
pip install scikit-learn
4.

# データを入れておく場所
mkdir data
cd data
wget http://yann.lecun.com/exdb/mnist/train-images-idx3-ubyte.gz
wget http://yann.lecun.com/exdb/mnist/train-labels-idx1-ubyte.gz
wget http://yann.lecun.com/exdb/mnist/t10k-images-idx3-ubyte.gz
wget http://yann.lecun.com/exdb/mnist/t10k-labels-idx1-ubyte.gz
5.

gzip -dc train-images-idx3-ubyte.gz | od -An -v -tu1 -j16 -w784 | sed 's/^ *//' | tr -s ' ' > train-images.txt
gzip -dc train-labels-idx1-ubyte.gz | od -An -v -tu1 -j8 -w1 | tr -d ' ' > train-labels.txt
gzip -dc train-images-idx3-ubyte.gz | od -An -v -tu1 -j16 -w784 | sed 's/^ *//' | tr -s ' ' > test-images.txt
gzip -dc train-labels-idx1-ubyte.gz | od -An -v -tu1 -j8 -w1 | tr -d ' ' > test-labels.txt

6.
適当にimage.pyみたいなの作って以下のコードをかく

#!/usr/bin/env python

import numpy as np
from sklearn.ensemble import RandomForestClassifier

# learninig data
data_training = np.loadtxt('data/train-images.txt', delimiter=' ')
label_training = np.loadtxt('data/train-labels.txt')

# test data
data_test = np.loadtxt('data/test-images.txt', delimiter=' ')
label_test = np.loadtxt('data/test-labels.txt', delimiter=' ')

# learn
#estimator = LinearSVC(C=1.0)
estimator = RandomForestClassifier()
estimator.fit(data_training, label_training)

# test
print(estimator.score(data_test, label_test))

7.

ターミナルで

python image.py

とかやると

0.9990166666666667

と出てきた。とりあえずやってみた感が強いので、明日あたりにちゃんと理解して編集、更新します。

物理屋さんのiPhoneアプリ開発日記１~X codeを使ってみる~

はじめに

クソゲーを作りたい。

これは俺の長年の悲願である。

もともと単純作業だったりすごくくだらないクソゲーが好きなので自分でもiPhoneアプリを作ってみたいと思い、Swiftの勉強を始めることにした。

もはやブログのタイトルとまったく関係なくなってきているが気にしない。

前回の記事

c++使いのSwift備忘録 - 物理の研究の備忘録

でSwiftの超基本的な文法について学んだので、今回はX codeを使ってiPhoneアプリを

作るのに必要なiOS独自のクラスに触れてみたいと思う。

実際の開発に比べるとまだまだ初歩の初歩という感じだが、物理の研究の合間にiPhoneアプリの

開発の勉強をするのはとても楽しいし、いい息抜きになるのでよしとする。

準備

MacでX codeを起動してCreate a new Xcode projectをクリックする。

今回はiPhoneで動くiOSアプリを作るのが目的なので、iOSにチェックを入れて、Single View Appを選んでNextをクリック。

プロジェクト情報を入力していく。正直なんでもいいのだけど、Project nameを

MyAppとかにしてそれ以外を自分の名前とかにした。下の3つのチェックは外しておく。デバイスをiPhoneにして、言語はもちろんSwiftにする。（Objective-Cも選べる）

保存場所はどこでもいいけど、とりあえずデスクトップにしておく。

作業画面になった。左側にAppDelegate.swiftとかViewController.swiftとかMain.storyboardとかいろいろ出てくるが、とりあえずMain.storyboardをクリックするとiPhoneの画面っぽい枠がでてくる。なんかそれっぽいぞ！

このストーリーボードってとこで、画面のどこに何を配置してどんな動きをして...みたいなのを作ってくっぽい。

右下のところがらいろいろ部品が選べるので、とりあえずLabelってやつを選んでストーリーボードにぴーっとドラッグしてみる。画面にLabelが貼り付けられた。

↓こんなん。

 右上にある丸が二つ重なってるみたいなやる（Assistant Editor）をクリックすると、に画面編成になって右側にViewController.swiftがでてくる。このViewController.swiftってとこに書いてあるプログラミングが、画面に何をどう表示するかをコントロールしてるっぽい。貼り付けたラベルをコントロールするために、Ctrlを押しながらラベルをドラッグしてViewController.swiftのclass ViewController: ...とか書いてあるとこの下にドラッグする。名前は適当にlabelとかにする。すると@IBOutlet weak var weak label: UILabel!という１行が加わった。

 これでlabelというクラス変数（インスタンス）を通してストーリーボード上のラベルを制御できるようになった。

ラベルの文字を変える。

label.text = "わーい"

の行をviewDiaLoal()関数の中に加える。

その状態で左上の再生ボタン的な奴をクリックする。

シミュレーションの画面に「わーい」と出てくればおっけー。

class ViewControllerってなんや

そもそもここに書かれているプログラミングはなんなのか。

クラスについて知っている前提で話をする。（クラスについては別記事にするかも）

class ViewController: UIViewController{

}

と書いてあるのは、これはViewControllerというクラスについての説明で、UIVewControllerというクラスを継承しているということだ。

このクラスにはVewDidLoad()という関数とdidReceiveMemoryWarning()という二つの関数（メソッド）を持っている。didReceiveMemoryWarning()というのはメモリがやばくなった時に読み込まれるメソッドなので今は気にしない。VewDidLoad()は画面（ビュー）の準備ができたときに読み込まれるメソットで、とりあえずここに変更したいことを加えれば良いことになる。overrideと書いてあるのは継承元のクラスであるUIViewControllerに同じ名前のメソッドがあるけどこっち優先的に読み込んでね！という意味である。ただ、実際にはUIViewControllerクラス内のVewDidLoad()を先に読んで欲しいので、最初にsuper.VewDidLoad()と書いてある。こう書くことで、先に継承元のUIViewControllerのメソッドを読んで、それから変更部分を加えたViewController()クラスのメソッドを読んでくれるということだ。

@IBOutletってなんや

次にドラッグした時にできた@IBOutlet weak var weak label: UILabel!について。

@IBOutletっていうのはストーリーボード上からラベルの情報に参照するのに必要らしいけど、うまく説明できないし正確に理解できているかもわからないので今はなんとなくで理解しておく。

var a:Int = 0

とか打つとInt型の変数aが定義されるのと同じように、UILabel型のインスタンスlabelが定義されて、このlabelというインスタンスを通してラベルの情報にアクセスできる。

なのでlabel.text="..."と書くことで、ラベルのテキスト情報にアクセスしてラベルの文字を変更することができる。UILabel!の"!"は非オプショナル型といってnil(つまりデータのない空の状態)を入れちゃいけません、という意味。逆に"?"をつけるとオプショナル型となりnilを入れることが可能となる。weakと書いてあるのはクラスの参照とかメモリがどう解放されるかとかそのへんに関わってくるっぽいけど、とりあえずは気にしない。

正直細かいプログラムの定義は使いながらなんとなくニュアンス掴んでいく部分もあると思うし、現時点ではこれくらいの理解でいいかな。

まだラベルの文章変えただけだし、次はもう少しアプリっぽい何かに挑戦したい。

おわり。

 

c++使いのSwift備忘録

今日から趣味でSwiftの勉強を始めたので普段俺が使ってるc++との違いなどを備忘録的にまとめてみる。

Swiftとは2014年にAppleによって発表された割と新しいプログラミング言語のことだ。Mac関連やiOSアプリを開発しようと思ったらSwiftとお友達になるのが１番の近道だと思われる。Swiftが出てくる前まではObjective-Cという言語がiOSの開発で使われていた。

実は何年か前に一度Objective-Cを勉強しようと思ったことがあったんだけど、いろいろと癖が強すぎて早々に挫折。最近になってSwiftについて少し調べてみたらObjective-Cに比べて随分とっつきやすい気がしたので再度挑戦してみる。

現在の開発環境はMac OS Sierra 10.12.5, X code8.3.2

ターミナルでswiftとコマンドを打ち込んでインタープリタとして使うこともできるけど、せっかくなのでX code開いて作業することにする。

四則演算とか

コメントアウトは//または/* */なのでc++と同じ。四則演算も+-*/%等まったく同じ

に使えるので楽。この時点で妙に親近感を覚えるね。まぁこのへんはどの言語でも

似たようなものなのかもしれないけど。

実行結果の出力とか

print(“hoge”)

みたいな感じで出力できる。c言語でいうprintf、c++でいうところのcout

Pythonはprintだからpythonと同じって考えていいのかな？

ちなみに文字列は

print(“aaa”+”bbb”)

みたいな感じで足し算で簡単に繋げられるので便利。c++のstringっぽい感じで使えるぞ。

てかSwiftって文末に;いらないのね。なんかpythonみたい。

型はおなじみのInt, Float, Double, Bool, Stringなどが普通に使える。

キャストも(Double)3みたいな感じだからc++とほぼ同じ。これはSwiftに限った話じゃないけど、Intと打った場合にどれくらいの容量が確保されるかはマシンスペックに依存する。なので、厳密に管理した場合はInt64みたいな感じで明示的にbit数指定することもできる。

変数や定数の宣言

c++だと

int a= 1;

と書いていた部分、Swiftだと

var a : Int =1

てな感じになる。:Intの部分で型を指定してるけど

var a = 1

みたいな感じで省略すると勝手に型推論してくれる。ただし型推論ばっかりしてるとコンパイルに時間かかっちゃうみたいなので、できるだけ明示的に型指定したほうがいいのかな？

定数の場合はvarの代わりにletを使う。なので

let b : Float = 10.0

はc++でいう

const float b=10.0

と同じ。ちなみにc++だと10.0と書く代わりに10.と打ってもよかったけど、Swiftだとエラー出た。ちゃんと0まで書かないといけないのね。

配列

 c++だと、可変長の時はvectorを使って、普遍なときは普通の配列を使うって感じで雑多で統一感がない感じだったけど、Swiftの場合可変、不変はvarとletを使い分ければよくて、わりとすっきりしている。

・不変の場合

c++だと

const int array[3]={1,2,3};

Swiftだと

let array : [Int] =  [1,2,3]

・可変の場合

var array : [Int] = [Int]()

array.append(1)

array.append(2)

...

みたいな感じで要素追加する。

print(array[0])

とやるとちゃんと1が出力された。なんか[Int]()って書き方慣れないな。初期化的なことをしてるのかな？

関数

Int型の引数を２倍にして返す関数を考えてみる。c++だと 

int test(int a){

   int b = 2*a ;

   return b;

}

//in main function

cout<<test(5)<<endl;

出力結果：10

Swiftだと

func test(a : Int)->Int{

    var b = 2*a; //or var b:Int =2*a

    return b

}

//in main fuction

print(test(a:5))

出力結果：10

こんな感じ。慣れるまで時間かかりそうだな。２変数だと

func test(a : Int,b : Int)->Int{...}みたいな感じで適当に書いてみたら回った。

for文

1から10までの和をfor文で回すとする。足すたびに出力もする。

c++だと

int sum=0;

for(int i=0;i<10;i++){

  sum+=i+1;

  cout<<"sum="<<sum<<endl;

} 

Swiftだと

var sum=0

for i in 0...9 {

    sum+=i+1

    print("sum=\(sum)")

}

出力結果

sum=1

sum=3

sum=6

sum=10

sum=15

sum=21

sum=28

sum=36

sum=45

sum=55

繰り返し回数が0...9で10回足されてることに注意。

\(sum)って書いてあるのはInt型の変数を文字列にぶち込みたい時にこう書けるらしい。

print("sum="+String(sum))

としても同じ。

cと同じ文法でfor文が回せる...?

少し調べてみると、このfor文の使い方はRubyに近いらしく、c言語ユーザーに馴染み深い記法

for(var i=0;i<10;i++)...

でも回ると書いてあった。俺的にはこっちの方が慣れ親しんでるので「やったー」と思いながら打ち込んでみたら

error: C-style for statement has been removed in Swift 3

なるエラーが出た。えぇ。。この書き方はSwift3で廃止になったみたいだ。ひどいよ...

まだ比較的新しい言語だから色々変わってる最中なのかな。ちなみに++,--もSwift3から廃止らしく、「i++じゃなくてi+=1のような書き方をしてね！」ということのようだ。

if文

cとほぼ同じ。cだと

if(a<10){

 ...

}

と書くとしたら、Swiftだと

if a<10{

 ...

}

と書けばいい。()がなくなっただけだね。elseも同様。これはすぐ覚えられそう。

X codeでプログラムを打ってみる

x codeにはPlaygroundってやつがあって、プログラムを打つとその場で出力結果を教えてくれるモードがある。今までのプログラムを打つとこんな感じ。

今回は別にSwiftの文法に限った話で、特別iOS開発に関する事柄はでてこなかったので、次はiOS独自のクラスを使いながらiOS開発、オブエジェクト志向について学んでいこう。

おわり

負の二項分布（Negative Binomial Distribution）について学ぶ

負の二項分布とは

負の二項分布を使う機会があったのだけど、定義がいっぱいあってよくわからないので調べたことをまとめてみる。

日本語版wikipediaだと一つの定義しかないが、英語版wikipediaの項目（Negative binomial distribution）が割と充実していたので英語版wikipediaを見ながら学ぶ。

NBDと呼ばれることもあるね。

二項分布については使われる機会も多く、少し調べてみると充実した解説が山ほど出てくるが、負の二項分布についてはあまり充実した情報はないようだ。英語版wikiによると全部で5つの定義があるみたい。式書くのめんどくさいので以下スクショ。

↓リンク

Negative binomial distribution - Wikipedia

ここでは試行回数：n、成功した回数：k、失敗した回数：r、成功する確率：pと定義する。

当たり前だが k+r=n の関係が成り立つ。

５つはそれぞれ

１：r回失敗する前に何回成功したか

２：r回失敗するのに必要な試行回数

３：k回成功する前に失敗した数

４；k回成功するのに必要な試行回数

５：n回の試行で何回成功したか（二項分布）

となっている。

一番わかりやすいのは５番目の定義で、これは二項分布の定義そのもの。

n回の試行に対して何回成功したかという確率分布である。

コインを100回投げたとき何回表が出たか、みたいなやつ。

成功した回数の分布を見ているのでf(k;,n,p)とkの関数になっていることがわかる。

4つ目の定義については確か日本語版wikipediaで説明されていた気がする。

このように、二項分布の場合は定義が一つだが、負の二項分布の場合はそれぞれ定義が微妙に違うあるので、もし「負の二項分布を使って〜」みたいな説明を見たときはどの定義なのかちゃんと注意する必要がある。

これら５つの式は割と簡単に導くことができる。試しに１つ目の定義について考えてみる。

n-1 (k+r-1)回の試行でk回成功して、n回目で失敗すればいいので、

 \begin{align}f(k:r,p) &=\begin{pmatrix}{r+k{-}1}\\{k}\end{pmatrix}p^{k}(1-p)^{r-1}(1-p)\\ &=\begin{pmatrix}{r+k{-}1}\\{k}\end{pmatrix}p^{k}(1-p)^{r}\end{align}

 てな感じ。他も同様に導ける。そんなに難しくないね。

せっかくなのでもう少し踏み込んでみよう。

実数に拡張する

 この定義式だとrもkも整数しかとれない。成功回数が何回か、みたいな議論をしてるので当然だけども。これを正の実数に拡張してみよう。そうすれば連続的な関数として取り扱うことができるようになるはずだ。全部の定義でやるのはかったるいので、１番目の定義に絞ることとする。

まずは階乗の部分をガンマ関数に

ガンマ関数についてはぐぐればいくらでも出てくるので詳しくは書かないけど、これを使えば階乗の概念を連続的な量に拡張できる。

\begin{align}\Gamma{(n+1)=n!}\end{align}

であることを考えれば

\begin{align}\begin{pmatrix}{r+k{-}1}\\{k}\end{pmatrix}&=\frac{(r+k{-}1)!}{k!(r-1)!}\\&=\frac{\Gamma{(r+k)}}{\Gamma{(k+1)}\Gamma{(r)}}\end{align}

となるね。

分布の平均値をパラメータとして使う

さて、これで一通り準備は揃ったような気もするけど、最後に分布の平均値（期待値）を計算してそいつをパラメータにしよう。

なんでそんなことするのか。

ポアソン分布にしろ、ガウス分布にしろ、平均値を表すパラメータが式の中に入ってる。なので成功確率pをそのまま使って表すよりも平均値を用いたほうが分布の形状を考えるとき扱いやすいよね。まぁ理由なんてなんでもいいんだけど、とにかく平均値を考える。

と言っても考え方は二項分布のときと同じで、試行回数にpをかければ成功回数になるし、(1-p)をかければ失敗回数になるので以下の式が成り立つはずだ。

\begin{align}\frac{\langle{k}\rangle}{p}=\frac{r}{1-p}=k+r\end{align}

\begin{align}{\langle{k}\rangle}=\frac{pr}{1-p}\equiv{m}\end{align}

 平均値を表すパラメータをmと定義した。

あれ、kには期待値記号が付くのにrに付かないのはなぜ？と思う人もいるかもしれないが、現在の定義だとkについての関数になっていてrは固定値なので期待値記号は付かない。もちろん別の定義のNBDの定義だと変数も変わる。

これらを用いれば最終的には

\begin{align}f(k;r,m)&=\frac{\Gamma{(r+k)}}{\Gamma{(k+1)}\Gamma{(r)}}\left(\frac{m}{r+m}\right)^{k}\left(\frac{r}{r+m}\right)^{r}  \end{align} 

あるいは 

\begin{align}f(k;r,m)&=\frac{\Gamma{(r+k)}}{\Gamma{(k+1)}\Gamma{(r)}}\left(\frac{m}{r}\right)^{k}\left(\frac{m}{r}+1\right)^{-(r+k)}  \end{align} 

とかくことができる。

上の方の式は英語版wikipediaにも載っているので参照にされたし。

そもそもなぜ”負”の二項分布なのか？

これでなんとなく負の二項分布についての定義はわかったと思うが、なぜこのような名前がついているのだろうか。実は負の二項分布には負の二項定理なるものを使えるのだ。

通常の二項定理は

\begin{align}(1+p)^{n}=\begin{pmatrix}{n}\\{k}\end{pmatrix}\sum_{k=0}^n{p^{k}}\end{align}

と書くことができる。通常nは自然数なでなければならないが、これを負に拡張して

\begin{align}(1+p)^{-n}=\begin{pmatrix}{-n}\\{k}\end{pmatrix}\sum_{k=0}^n{p^{k}}\end{align}

なるものを考える。ここでやっかいになるのが \begin{align}\begin{pmatrix}{-n}\\{k}\end{pmatrix}\end{align}であるが、これを変形すると

\begin{align}\begin{pmatrix}{-n}\\{k}\end{pmatrix}&=\frac{(-n)(-n-1)...(-n-k+1)}{k!}\\&=\frac{(-1)^k(n+k-1)(n+k-2)...(n)}{k!}\\&=(-1)^{k}\begin{pmatrix}{n+k-1}\\{n-1}\end{pmatrix}\end{align}

 と変形することができる。これを負の二項分布の式に適用すると、

\begin{align}\begin{pmatrix}{r+k{-}1}\\{k}\end{pmatrix}p^{k}(1-p)^{r}=\begin{pmatrix}{-r}\\{k}\end{pmatrix}p^{k}(1-p)^{r}\end{align}

と、すっきりする形に変形することができた。この変形によってどれだけすっきりしたかイメージできるだろうか。例えば、通常の二項分布における\begin{align}\begin{pmatrix}{n}\\{k}\end{pmatrix}\end{align}の部分は、nが固定の値でkが変数なので、例えば和を取る場合など二項定理の適用がしやすい。一方で負の二項分布の最初の定理の場合だと、\begin{align}\begin{pmatrix}{k+r-1}\\{r-1}\end{pmatrix}\end{align}の部分において、rが固定、kが変数になっているので、二項定理と比べて定数、変数の場所が逆になってしまっている。変形後は\begin{align}\begin{pmatrix}{-r}\\{k}\end{pmatrix}\end{align}となっているので\begin{align}\begin{pmatrix}{定数}\\{変数}\end{pmatrix}\end{align}の形に揃えることができた。これはこのあと出てくる母関数の計算においても用いる。

母関数を計算して平均と分散を計算してみる

さて、負の二項分布の平均と分散を計算してみよう。平均値についてはさっき求めたけど、母関数を計算すればn回微分でn次のモーメントorキュムラントが求まるので便利。１次のモーメント、キュムラントは平均値にあたり、2次のモーメント、キュムラントが分散にあたる。

最後に二項分布と負の二項分布、何が違うのか比べてみる。

二項分布の場合

二項分布の式は上にあるように

 \begin{align}f(k;n,p)=\begin{pmatrix}{n}\\{k}\end{pmatrix}p^{k}(1-p)^{k-1}\end{align}

である。

モーメント母関数M（θ）の定義は

\begin{align}\sum_{k}e^{\theta{k}}f\end{align}

一応確認しておくと、fは確立密度関数である。fの部分に二項分布の式をぶち込むと

\begin{align}M(\theta)=\sum_{k}\begin{pmatrix}{n}\\{k}\end{pmatrix}(p{e^\theta})^{k}(1-p)^{n-k}=(pe^{\theta}+1-p)^n\end{align} 

等式の変形で二項定理を用いた。

これをn回微分すればn次のモーメントが求まる。モーメントから計算してもいいんだけど、今回はキュムラントからアプローチしてみる。キュムラント母関数の定義はK(θ)=lnM(θ)なので、モーメント母関数のlogをとるだけ。モーメント同様キュムラント母関数をn回微分すると、n次キュムラントが出てくる。ｎ次キュムラントをCnとすると、

\begin{align}C_{n}=\frac{\partial}{\partial\theta}K(\theta)|_{\theta=0}\end{align} 

となる。

平均と分散を調べたいのであれば、１階微分と２階微分を調べればいい。M(θ)のlogをとると肩のnが落っこちてくるだけなので

\begin{align}K(\theta)=n(pe^\theta+1-p)\end{align} 

となる。これをθで微分していく。細かい計算は省略するけど、計算すると

\begin{align}\mu=np\end{align} 

\begin{align}\sigma^2=np(1-p)\end{align} 

ここでμ、σ^2はそれぞれ平均、 分散を表し、C1, C2にあたる。また、

\begin{align}\sigma^2=\mu(1-p)\end{align}

の関係が成り立つ。1-pをεと定義しなおせば

\begin{align}\sigma^2=\mu\epsilon\end{align}

というすっきりした関係が導ける。ここでpは確率を表すので、ε<1となる。つまり、

二項分布においては平均は分散より大きいということがわかる。

負の二項分布の場合

負の二項分布でも同じことをやってみよう。式は変形後の

 \begin{align}f(k;r,p)=\begin{pmatrix}{-r}\\{k}\end{pmatrix}p^{k}(1-p)^{r}\end{align}

を用いる。

モーメント母関数は

\begin{align}M(\theta)=\sum_{k}\begin{pmatrix}{-r}\\{k}\end{pmatrix}(p{e^\theta})^{k}(1-p)^{r}=(\frac{1-p}{1-pe^{\theta}})^r\end{align} 

等式の変形で負の二項定理を用いた。

同様にM(θ)のlogをとると肩のｒが落っこちてくるだけなので

\begin{align}K(\theta)=r(\frac{1-p}{1-pe^{\theta}})\end{align} 

となる。同様の計算すると

\begin{align}\mu=\frac{pr}{1-p}\end{align} 

\begin{align}\sigma^2=\frac{pr}{(1-p)^2}\end{align} 

ここで

\begin{align}\sigma^2=\frac{\mu}{(1-p)}\end{align}

の関係が成り立つ。二項分布の場合は平均に1-pを掛けると分散になるのに対して、負の二項分布の場合は割り算になっているという非常に綺麗な関係を導くことができた。さらに1/1-pをεと定義しなおせば

\begin{align}\sigma^2=\mu\epsilon\end{align}

と、二項分布と全く同じ形で表すこともできる。

このとき、負の二項分布においてはε>1となる。つまり、

負の二項分布においては平均は分散より小さいということがわかる。

途中まで全然違う計算をしているようで、最後に綺麗に同じような形で表すことができるのは驚きだ。ちなみに、ポアソン分布においては平均と分散は同じなので、無理やりまとめるとすれば

二項分布：ε>1 (平均>分散)

ポアソン分布：ε=1 (平均=分散)

負の二項分布：ε<1 (平均<分散)

となるね。よくできてるなぁ。

ちなみに、モーメント母関数から計算しても、別の定義の負の二項分布から計算しても同様の結果を導くことはできる。さらに言えば、平均と分散を計算しただけなら別に母関数を使わなくてもできる。ただ、俺はキュムラント母関数が便利で応用がきくので好きだ。

追記：３次以上の場合

上では平均と分散について、εを定義すれば二項分布と負の二項分布がまったく同じ関係式で表すことができるを示した。

実はこれ、３次以降のすべての次数のキュムラントにおいても一致する。

文明の利器（mathematica）を使って、キュムラント母関数を微分して6次のキュムラントまで計算してみると、どちらもまったく同じ結果になった。すっきり。

\begin{align}&C_{1}=\mu\\&C_{2}=\mu\epsilon\\&C_{3}=\mu\epsilon(2\epsilon-1)\\&C_{4}=\mu\epsilon(6\epsilon^{2}-6\epsilon+1)\\&C_{5}=\mu\epsilon(2\epsilon-1)(12\epsilon^{2}-12\epsilon+1)\\&C_{6}=\mu\epsilon(120\epsilon^{4}-240\epsilon^{3}+150\epsilon^{2}-30\epsilon+1)\end{align} 

おわり。