5 統計編
この章では，エコノメにまつわるテクをざっと眺める．エコノメに縁がない人も，行 列・ベクトルの操作などは参考になるかもしれない.*3
言うまでもなく，統計分析をするなら他の統計用ソフトを使ったほうがよい． Mathe-maticaは外部言語としてRも動かすことができる(RLink)が，直接Rを使わない理由 がよくわからない．行列の計算は明らかにMatlabの方がフレンドリーだ．
version 10 からはDatasetを用いることで，関係データベースを操作しやすくなって いる．が私はよく知らないのでここでは触れない．
る．これならお手軽だ．
5.1.3 行列の条件抽出がしたい
パネルデータを扱う上でよくあるのが，「i= Japanの行だけ取り出したい」，「t= 2000 の行だけ取り出したい」，などの操作だ．この手の操作はStataが手軽だからStataを使 おう．Rだとdpylr，Pythonではpandasを使いこなせばこれまた早い．
Mathematicaの場合，特定の値を探すにはSelect またはPositionを使う．たとえ ばSelectを使って(虚数を含まない)ベクトルxの中から104 以上になっているものだ けを抽出するには
Select[x, # > 10^4 &]
とする．
Positionを使う練習として，適当な行列を作り，その中で何らかの条件を満たす行を
抽出してみよう．
mat = Prepend[
MapThread[Prepend, { RandomReal[{-1, 1}, {5, 3}],
{"Antras", "Bagwell", "Costinot", "Dixit", "Eaton"} }] (* add 1st row *) , {"name", "a", "b", "c"}] (* add 1st column *)
mat // MatrixForm
2行目後半では，-1以上1以下の一様分布から5×3行列を作れ，と命じている．あとは それに適当なラベルを付けている．以下のような行列が生成される．
name a b c
Antras 0.389616 0.141336 0.984982 Bagwell 0.782802 0.261349 −0.798015 Costinot 0.675463 0.401904 0.134947
Dixit 0.965734 −0.58665 0.545197 Eaton −0.44715 −0.502962 0.24322
この行列のうち，2列目(a)が1/2以上になっているものだけを取りだそう．
mat[[
Flatten@Position[mat[[All, 2]], x_ /; x > 1/2]
]]
% // MatrixForm
結果は以下の通りになる:
Bagwell 0.782802 0.261349 −0.798015 Costinot 0.675463 0.401904 0.134947
Dixit 0.965734 −0.58665 0.545197
探す条件が量的変数でなく質的変数の場合も同様にできる．たとえばdataという行列 の第3列目にCountryの情報が入っているとする．Countryが特定の値を取る行だけ抜 き出したい場合，その場所を探す関数を定義しておくと便利だ:
selectCountry[i_] := Flatten@Position[data[[All, 3]], i]
こうしておけば，selectCountry["Japan"]で，i=Japanとなるところだけを抜き出 すことができる．
同様に，第四列目にYearの情報が入っている場合，
selectYear[t_] := Flatten@Position[data[[All, 4]], t]
などとする．年月・国・産業など組み合わせていくことももちろん可能．
データからJapanだけを取り出したものを頻繁に使う場合，
japan[xit_] := xit[[selectCountry["Japan"]]]
などとする．japan@xitで，ベクトルxitからi =”Japan”だけを取り出したベクトル ができあがる．
反対に，Japanを取り除きたい場合は
xJapan[xit_] := Delete[xit, Partition[selectCountry["Japan"],1] ]
などという関数を定義しておけば，xJapan@xitで，Japan成分を除いたものができあ がる．
なお，データの頭部に系列名のラベルが入っていたりして行番号が 1 つずれるなら，
selectCountryなどから1を引けばよい．これに限らず行列番号が1つずれたりするこ とはしばしばあるが，MostやRestやLastをうまく使えば対応は難しくない．
5.1.4 簡単なベクトル・行列を生成したい
Mathematicaは行列・ベクトルの操作性が低いので，あらかじめ関数を定義しておく
とよい．
たとえばn次元1が並ぶベクトルを作りたい場合，
ones[n_] := Table[1, n]
などとすればよい．ones[5]とすれば，{1,1,1,1,1}ができる．コードの中にTableが 乱立し可読性が落ちるのを避けたいときに便利である．Tableはこのように何かを複製す るのにも有用だ．なお，単純に定数(やベクトル)を並べるだけならTableでなく
ones[n_] := ConstantArray[1, n]
のようにConstantArray関数を使うほうがよい．
TableとConstantArrayの微妙な違いは，Tableでは毎回評価されることにある．た とえば
Table[RandomReal[], 5]
ConstantArray[RandomReal[], 5]
を走らせてみればそれぞれ次のような出力を得る．
{0.538507,0.840555,0.0808903,0.632761,0.973275} {0.930344,0.930344,0.930344,0.930344,0.930344}
評価の仕組みについて詳しくは，Mathematicaのチュートリアル「評価の原理」を見よ．
連立方程式を解くとき，スカラーを複製したベクトルを用いる必要がないことも多い．
たとえば
Solve[{x+2y, -x+y} == 0, {x,y}]
と書けば，右辺は0ベクトル(0,0)と勝手に認識され，x+ 2y = 0と−x+y= 0という 2本の連立方程式として勝手に認識してくれる．右辺は0以外のスカラーでも同様．
ダミー変数を作るときも，Tableにはよくお世話になる．固定効果モデルを回したいと きに使う差分オペレーターもTableを工夫すればよい．
n次元の単位行列がほしい場合は，自前の関数を定義せずとも IdentityMatrix[n]
でよい．同様に，対角行列はDiagonalMatrixを使うとよい．たとえば ( 1−x1 0
0 1−x2 )
のような対角行列を作りたい場合は次のように書けばよい: DiagonalMatrix[1 - #] & [{x1, x2}]
すべての要素が同じ値となるような行列を生成する場合，ConstantArrayを用いる．
ConstantArray[0, {3, 4}] (* 3*4 matrix of 0 *)
定 数 で は な く よ り フ レ キ シ ブ ル に ，i, j 要 素 に 関 数 f(i, j) を 適 用 し た い 場 合 ， Table[f[i,j],{i,2},{j,3}]のようにしてもよいが，Arrayを用いることもできる．
Array[f, {2, 3}]
を実行すると次のような出力となる．
( f(1,1) f(1,2) f(1,3) f(2,1) f(2,2) f(2,3)
)
下三角行列・上三角行列はLowerTriangularize, UpperTriangularize, で作る．た とえば下半分がすべて1の4×4行列は次のように作る:
LowerTriangularize[Array[1 &, {4, 4}]] // MatrixForm
特定の要素以外はすべて0，という行列を作るにはSparseArrayを使う．
乱数で対称行列を作成したいとき，次のようにする．
mat = RandomReal[1, {3, 3}] (* generate a 3*3 matrix *)
matLT = Transpose@LowerTriangularize@mat (* lower triangle *)
mat = ReplacePart[mat, {i_, j_} /; j > i -> matLT[[i, j]] ] ; (* merge *) mat // MatrixForm (* symmetric matrix of random var *)
結果は次の通り:
0.657057 0.555953 0.170305 0.555953 0.716725 0.719343 0.170305 0.719343 0.0238045
このコードの3行目については第3.1.5, 5.1.9章も参照．j > iになるという条件，つま り上三角の部分だけ，ReplacePartで下三角の部分を代入している．乱数を含んだ正方 行列を生成している1行目はmat = Array[RandomReal[] &, {3, 3}] と書いても同 じ結果になる．
あるベクトル{x1, x2, x3, . . . , x10}を，{x1, x1, x1, x2, x2, x2, x3, x3, x3, . . . , x10, x10, x10}, という形で増殖させたいときがある．こういう場合は
n = 3 (* 増殖したい数 *);
Flatten@KroneckerProduct[x, Table[1, {n}]]
と 1 の並んだベクトルのクロネッカー積を計算するとすっきりする．Table 部分は 上で定義した ones[n] に代替してもよい．同じことだが，Table でたくさん並べて，
Flatten@Transposeをかけてもよい．
0か1を取るような何らかのindicatorベクトル・行列を作成したいときは上述の通り
Tableを工夫すればよいが，べき集合やそれっぽいものがほしければ Tuplesも知って
おくと便利だ．たとえば「長さnで，0か1を取るあらゆるベクトル」を生成したいと する．たとえばn= 4のとき，(0,0,0,0),(1,0,0,0),(0,1,0,0),· · · ,(1,1,1,1)のように 24 = 16個のベクトルがほしいとする．こういう場合，
n = 4 (* 長さ *);
Tuples[{0, 1}, n]
とすれば24×4行列としてほしいベクトルがまとめられて出力される．これを転置して MatrixFormで見ると次のようになる:
0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 0 0 1 1 0 0 1 1 0 0 1 1 0 0 1 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1
5.1.5 列和・行和を計算したい
Totalを使うと，列和が計算される．Transposeと組み合わせて転置すれば，行和も
計算することができる．
似た関数に級数和 Sum があるが，こちらは ∑
ixi に対応する．こちらは ∞ まで合 計させることもできる．数値的に近似解がほしければNSumを使う．同様に，級数積は Productを用いる．
Totalは列に単純和という操作を施す，ということだが，単純和以外の操作を施したい
ときもある．要はRでいうapply()を使いたい場合だ．こういうときはMapと組み合わ せる．たとえば行ごとに単純平均を求めたい場合，
Map[Mean, x]
とする．列ごとにしたい場合は転置すればよい．Mapは/@と基本的には同じ操作である．
なんらかの関数を@でつなげたときの出力とは行と列が入れ替わる．Map は類似の関数 (MapThreadやMapAllなど)がいくつかあるので，詳しくはヘルプを参照．
5.1.6 変数名に連番を振りたい
{x1, x2, x3, . . . , x10}, のように連番になる変数を作りベクトルに格納したいときがあ る．私は次のように生成している
nN = 10 (* 連番の長さ *);
vecx = ToExpression@Table["n" <> ToString[i], {i, nN}]
これでvecxに，{x1, x2, x3, . . . , x10}が格納される．行列の添え字を作りたい場合も同 様である．
関数としてx[1], x[2], x[3]. . .x[10], のような形のベクトルを生成したい場合は Arrayを使う:
Array[x, nN]
5.1.7 ロング形式・ショート形式を自在に行き来したい
地域i, 時間tについてある変数のパネルデータがあるとする．データをロング形式(た とえばIT 次元ベクトル)にもショート形式 (I ×T 行列)にもできると，何かと便利で ある．
ロング形式からショート形式にするときにはPartitionを使えばよい．逆はFlatten だ．詳しい使い方はそもそものデータの並び方による．
5.1.8 欠損値を消したい
version 10よりDeleteMissingを使うと欠損値を削除できるようになった．
5.1.9 リストの要素を書き換えたい
行列やベクトルの一部を別の値に置き換えたいときにはいくつか方法がある．シンプル なのは，単に値を割り当て直すことだ．たとえば行列xの(2,3)要素を0に置き換えたい 場合，
x[[2,3]] = 0 とする．
上のようにして割り当て直すとx自体が書き換えられるが，xはそのままにして別の新 しいリストを生成したい場合，ReplacePartを使う．たとえば行列x の(2,3)要素を0 に置き換えたものをyとする場合，
y = ReplacePart[x, {2,3} -> 0]
と書く．
書き換えたい位置を探すにはPositionを使うとよい．負の値，虚数，Indeterminate， など置き換えたいものがある場合はそれぞれどこにあるのかパターンマッチングで探す．
たとえばxの要素のうちどこに負の値や不定量があるかを知りたいなら
Position[x, _?(# < 0 &)] (* positions where is a negative value *) Position[x, Indeterminate] (* positions where is Indeterminate *)
などとする．複素数Complexesは当然ながら実数も含んでしまうので，虚数を探すとき は逆に実数の位置をPosition[x, _?(Element[#, Reals] &)]などで探す．虚部がゼ ロじゃない，という条件で以下のように探すこともできる:
Position[x, _?(Im[#] != 0 &)] (* positions where is an imaginary value *) ただしこれらのコードは∞や−∞などが含まれていると思い通りに行かないこともあ るので注意．
置き換えるのではなく，行を新たにどこかに挿入したい，という場合は Insert や Append, Prependを使う．列を挿入するなら，Transposeと組み合わせつつInsert を 使えば良い．
ついでにデータ・クリーニングにおける初歩的な注意しておくと，欠損値を恣意的に何 かに置き換えることは慎重であるべきだ．たとえば年収のデータが欠けているとき，0を 代入し年収ゼロとみなすのはおかしい．
