ネコでもわかる Stata 入門
東京大学大学院経済学研究科博士課程
別所俊一郎
1
はじめに
Stata は，広く用いられている統計ソフトパッケージのひとつです．Stata の web site1 の
表現を借りれば，
「Stata is used by medical researchers, biostatisticians, epidemiologists,
economists, sociologists, political scientists, geographers, psychologists, social scientists,
and other research professionals needing to analyze data. 」ということになります．
統計や計量経済学の専門家に言わせれば，最新の推定手法等が Stata のような canned
package に搭載されるまでには時間がかかるために，この手の統計ソフトはそれほど有用
ではないそうです．しかし，応用計量経済学を専門とするわれわれからみれば，標準的な
手法がほとんど組み込まれており，複雑な収束計算等をプログラミングすることなく分析
を進めることができる，非常に使用価値の高いソフトです．このように考えると，Stata を
一通り使うことができることは，研究を進めるうえで有効なスキルのひとつだといえます．
しかし，他にも統計ソフトがいろいろとあるなかで，ここでなぜ Stata を使うのか，疑
問に思われるかもしれません．計量経済学の分野では，TSP，Rats，SHAZAM，E-views，
LIMDEP 等がありますし，SPSS，SAS といったソフトもよく知られています．行列言語
の Ox や C，Matlab や GAUSS 等もしばしば用いられます．これらのソフトはそれぞれに
特色を持っています．Ox や C，あるいは Matlab や GAUSS は，行列計算や収束計算を自
分でプログラミングしなければなりませんが，推定手法の自由度はもっとも高いといえま
す（もちろん，計量経済学の教育効果もあります）．TSP や Rats，E-views は時系列解析
には非常に有効です．Stata は対照的に，個票（マイクロデータ）の処理に強みを持ってい
ます．
Stata が広く用いられ，世銀等が個票を Stata で提供している理由はこれだけではあり
ません．第１に，Stata はサンプル数に制限がありません．変数の数も 2000 程度まで用意
できます2 ．したがって，メモリの許す限り，いくらでも大きなデータを扱うことができま
す．第２に，Stata のメインプログラムがそれほど大きくないということが挙げられます．
メインプログラムが大きくない代わりに，さまざまな統計手法が ado ファイルというかた
ちで提供され，web 上で入手できます．ado ファイルはサブルーチンであり，メインプロ
グラムから呼び出されて利用されます．この仕組みにより，Stata は高い拡張性を持ってい
ます．多くの統計学者が ado ファイルを作成し，Stata corp. は日々そのデータベースを更
新しています．
前置きはこれくらいにしてさっそく Stata を使ってみましょう．統計ソフトは，計量経
済学や統計学の問題を解き，実際に使っていくことで慣れていくものでしょう．
1 http://www.stata.com/
2 Stata
7 special edition だと制限がありません．
1
インストールしたらこれはやっておこう
2
起動の方法は PC 等の環境によって異なるので，ここでは省略します．さて，インストー
ルして最初に起動したら，いくつかやっておくとよいことがあります．
2.1
アップデートしよう
もし PC がネットにつながっている場合には3 ，まずアップデートをしましょう．先に書
いたとおり，Stata 本体や ado ファイルは日々更新されています．これを取り込んでおか
ない手はありません．
メニューバーの「Help」から「Official Update」を選びましょう．開いたウィンドウの
下のほうに「Recommendation : compare these dates with what is available from」とあ
り，その下に「http://www.stata.com/」とありますので，これをクリックしましょう．す
ると，しばらく経つと「update query」が出てきます．ここで，
「currently installed」と
「last available」が一致していないと4 ，一番下に再び「Recommendation」があって，イン
ストールしなさい，というようなことが書いてあります．あとは，指示に従ってダウンロー
ドとインストールをしましょう．ダウンロードは，場合によってはかなりの量になります
から，ダイヤルアップ等でのダウンロードには気をつけましょう．本体（wstata.exe ）が
アップデートされている場合には，もとのファイルの名前を変えてファイルを移動して…
とちょっと手間がかかりますが，がまんしましょう．
このようなアップデートは，インストール直後だけではなく，数ヶ月に一度は行うとよい
でしょう．もちろん，アップデートされた ado ファイルをすぐに使う可能性は低いですが…．
ado ファイルは，Stata Corp. が公式に認定したものもありますが，それ以外にも Stata
ユーザが作成したものもあります．統計手法によっては，公式な ado ファイルではなく，こ
のようなユーザ作成の ado ファイルを使うと便利なことがあります．たとえば，Box-Cox
回帰などはユーザ作成の ado ファイルで提供されています．さて，このような ado ファイ
ルを使うときは，メニューバーの「Help」から「STB and User-written programs」
を選びましょう．そこを適当にたどっていけば，所望のファイルを入手できるでしょう．
ユーザ作成の ado ファイルなんて使うの？と思うかもしれません．たしかに，よく使用さ
れる ado ファイルは公式 ado ファイルとしてアップデートされています．help で「こんな
ことできたらいいな」と探してコマンド名まで分かったにもかかわらず実行するとエラー
が出る場合，そのコマンドはユーザ作成 ado ファイルで提供されている可能性があります．
そのようなときには，上述の方法に従って ado ファイルを入手する必要があります．
3 つながっていない場合は……どうしましょうか？
4 バージョンが変わると，アップデートがとまります．02 年 5 月現在，Stata 6 についてはアップデートが行
われていません．バージョンをアップするといいよ，という広告が出ます．
2
2.2
カスタマイズしよう
Stata を起動すると，４つのウィンドウが開きます．それぞれ，
「Stata Results」
「Review」
「Stata Command」という名前がついています．それぞれのウィンドウが
「Variables」
何を意味するかはあとで述べることにしますが，気になるひとは，これらのウィンドウの
大きさを変えたり，フォントを変えたりしましょう．フォントは，それぞれのウィンドウ
の左上のボックスをクリックし，一番下の「Fonts」から変更できます．フォントやウィン
ドウの位置の設定は，メニューバーの「Prefs」→「Save Windowing preferences」で
保存できます．
では次に，画面の説明をしておきましょう．
画面の説明
3
Stata を起動すると，なかに４つのウィンドウが開きます．ログ（記録）やグラフはまた
別のウィンドウになります．４つのウィンドウはそれぞれ以下のような役割を果たし，利
用されます．
3.1
Stata Results ウィンドウ
起動時に右上にある大きなウィンドウです．ここには処理の結果が表示されます．昔の
バージョンではスクロールバーがなく，出力結果が多い場合にはログをとっていないと以前
の処理が見えなかったということがありましたが，いまはスクロールバーがあるため，あ
る程度なら遡って出力結果を見ることができます．ここに直接入力することはできません．
3.2
Stata Command ウィンドウ
起動時に右下にある細長いウィンドウです．Stata は基本的には対話型のソフトですか
ら，ここにコマンドを入力します．カーソルが使えなかったりしますが，あんまり気にす
ることはありません．Stata で唯一，キーボードからのコマンド入力を受け付けてくれる
ウィンドウです．
3.3
Review ウィンドウ
起動時に左上にあるウィンドウです．それまでに Stata Command ウィンドウに入力し
た結果が順次表示されていきます．ここに表示されている以前の入力をクリックすると，
Stata Command ウィンドウにそれが入力されます．同じ，あるいは似たような入力をする
ときには便利です．
3
Variables ウィンドウ
3.4
起動時に左下にあるウィンドウです．現在メモリ上にあるデータに含まれる変数の一覧
が表示されます．左が変数名，右が変数につけられるラベルです．表示されている変数を
クリックすると，Stata Command ウィンドウにそれが入力されます．変数名の入力ミスを
防ぐことができます．
作業を始める前にやっておこう
4
起動してデータを読み込むまえ5 に，次のようなことをやっておくとストレスがたまりま
せん．
4.1
使用メモリを確保しよう
Stata は起動時にデフォルトでは 3MB のメモリ領域を確保しますが，とくに個票などの
大きなデータセットを使用するときには，これでは足りないことが往々にしてあります．自
分の使おうとするデータが大きいと思うときには，最初にメモリ領域を確保する命令を出
しておきましょう．
set mem 30m
これで，メモリを 30MB 確保できます．データセットと使用するメモリ領域の大きさの正
確な関係はよくわからないのですが，これくらいあればたいがいは大丈夫です．
4.2
ログを取ろう
画面の説明のところで述べたとおり，結果は Stata Results ウィンドウに表示されます
が，結果が長くなると遡って見ることができなくなります．これに対処するために，ログ
をとっておきましょう．
メニューアイコンの左から４番目，カーソルを当てると「Open Log」と表示されるアイ
コンをクリックしましょう．すると，ウィンドウが開くので適当な名前をつけてログファ
イルを開きましょう．もちろん，既存のログファイルを使うこともできます．そのときには
「Append（付け足していく）」か「Overwrite（上書きする）」かを選ぶことになります．
ログファイルのウィンドウが開くと，以後の結果はここに記録されていくことになりま
す．このウィンドウは Stata Results ウィンドウと異なり，いくらでも遡って結果をみるこ
とができます．ログファイルを見るときには，さきほどのアイコンの右隣のアイコンをク
リックしましょう．もちろん，ログファイルはただのテキストファイルですから，Stata を
起動させなくても，メモ帳や秀丸で見ることができます．
ログを取りたくないときには，開くときに使ったアイコンをクリックして「Close log
file」を選びましょう．
5 Stata
データセットと関連付けておくと，自動的に Stata が起動しますが…．
4
ログファイルウィンドウを選んだ状態で印刷アイコンをクリックすると，ログファイル
を印刷してくれます6 ．ほかの編集ソフトから印刷した場合と異なり，入力部分が太字にな
るなど，やや見やすく印刷を行うことができます．
4.3
画面を流れさせよう
さて，ログをとる準備ができたら，次のコマンドを入力するといいかもしれませんが，こ
れは好みによります．
set more off
Stata Results ウィンドウでは，遡ってみることのできる量に限度がありますが，その代わ
り，出力が長くなるときには，いったん出力を停めて
-- more -という表示が出ます．
「ここまでの結果はいいですか？よければ次を表示しますよ」という
ことですが，いちいちこれに反応するのがめんどうなばあいなどには，さきほどのコマン
ドを打っておくと，一気に出力してくれます．ログを取っていれば，結果はそちらで確認
できるので，do ファイル（後述）などを使う場合には有力です．
データを入力しよう
5
Stata は，Stata 形式のデータしか処理できません．なので，データを Stata 形式に変換
する必要があります．Stata 形式のデータは拡張子.dta となります．データの入力方法に
はいくつかありますが，ここではそのいくつかだけ紹介します．その前にまず，データの
形式について説明しておきましょう．
5.1
データの型
Stata は文字データも受け付けてくれるそうですが，僕が使ったことがないのと，使わな
くてもそんなに不便ではないことから，データは全て小数点表示の実数とします．桁数につ
いても制限があるかもしれませんが，気にする必要があるケースはほとんどないでしょう．
データはスプレッドシート１枚に納められます．それぞれの変数に対してサンプルが対
応します．変数には「Name」と「Label」をつけることができます．Name はいわゆる変
数名で，大文字小文字のアルファベット，数字，下線で８文字以内でつけます7 ．変数名は
この条件内でなんでもいいのですが，以下の名前は先約があるので使いません．
_all _b byte _coef _cons double float if in int
long _n _N _pi _pred _rc _se _skip using with
6 希望する範囲だけ印刷する方法がよく分からないのですが…．
7 Stata
7 以降は８文字制限がなくなったようです．
5
一方，Label のほうにはそのような制限はありません．スペースが入ってもかまいませ
ん．左下の Variables ウィンドウには，左側に Name，右側に Label が表示されます．ある
いは，出力結果によっては，Label が表示されることがあります．
変数はここでは小数点表示の実数ですが，欠損値は「.」で表現されます．
ではデータを入力してみましょう．
5.2
手入力するっ
メニューアイコンの右から４個目，虫眼鏡がついていないほうのアイコンをクリックす
ると Data editor が開きます．ここに直接入力するのも，もちろん，データ入力方法のひ
とつです．でも，ぼくは使ったことがありません．
5.3
Excel などから貼り付ける
よく用いられる方法のひとつは，ある程度 Excel 等で加工したデータを貼り付けるとい
う方法です．このとき，一番上の行に変数名を書いておくと，自動的にこれが Name にな
ります．
あらかじめ Excel 等で範囲を指定してコピーしておいて，Data editor を開きます．一番
左上をクリックし，メニューバーの「Edit」→「Paste」とすれば，これでデータの入力
は完了です．データを追加するときにも，データの一番左上になる部分をクリックして同
様の操作をすることで，データを入力することができます．
5.4
タブ，コンマ区切りのテキストデータを読み込む
insheet コマンドを使うそうです．
5.5
固定長のテキストデータを読み込む
この場合，ディクショナリファイルを作って読み込むという作業が必要になります．ディ
クショナリファイルは拡張子を.dct とします．これを次のコマンドで呼び出します．
infix using "C:\My Documents\sick.dct"
一方，ディクショナリファイルの中身は次のようになっています．
infix dictionary using "C:\My Documents\sick.dat" {
hhold
1-5
...
}
これでまとめて，固定長データであるsick.dat の１文字目から５文字目までを変数
名hhold として読み込む，という作業が実行されることになります．
6
5.6
Label をつけよう
ふつうは，データを読み込んだ段階では Name はついていても Label はついていません．
そこで Label をつけておきましょう．作業のときになにかと便利ですが，つけなければな
らないというものではありません．変数名が分かりやすければなくてもかまいません．
一番簡単な方法は，Data Editor の変数名の部分をダブルクリックして，
「Stata Variable
Infomation」のウィンドウを開き，Label のところに入力する方法です．しかしこれだと，
Do ファイルなどに対応できません．
コマンドをつかってラベルをつけるには以下のようにします．
label var hhold "Household ID"
これで，変数名hhold である変数に，ラベルHousehold ID がつくことになります．
5.7
保存しよう・開こう
このようにして作成したデータを Stata 形式で保存しておきましょう．Stata 形式での保
存では，Name，Label なども保存されます．ただし上書きすることになるので，データに
変更を加えている場合には気をつけましょう．
保存された Stata 形式のデータは，一般には Stata を関連付けされているので，ダブル
クリックすることによって Stata で開くことができます．しかし，Stata はデフォルトで
は 3MB しかメモリを確保しないので，大きいデータセットの時には Stata だけ起動して，
ファイルが開けないことがあります．そのような場合は，set mem コマンドでメモリを確
保した上であらためて開きましょう（「File」→「Open」）．あるいは，デフォルトで確保
するメモリを拡張しておくという方法もあります8 ．
さて，いよいよデータを操作しましょう．といきたいところですが．
コマンドを使う前に
6
6.1
Help を利用しよう
Stata は Help 機能が充実しています．使いたい統計手法があれば，メニューバーの「Help」
→「Search」を用いて，コマンドの名前や使い方を調べることができます．
また，エラーメッセージが出てきたときには，その文頭にあるエラーコードを使って，エ
ラーの内容を調べることができます．たとえばエラーコードが 191 だったとすると，直後に
search 191
と打てば，エラーの原因が出力されます．
8 あるだろうな，というだけで僕はよく知りません．後述するように，Do ファイルに，メモリを確保する，ファ
イルを開くという一連の作業を書いておく方法もあります．
7
6.2
Data Browser を利用しよう
Data Editor の右隣に，虫眼鏡が描いてある似たようなアイコンがあります．これをク
リックすると，Data Browser が起動します．Data Browser は，その名の通り，データを
眺めるウィンドウです．Data Editor では，行った操作（とくにソート）等が保持されます
が，Data Browser ならそんなことはありません．つまり，間違ってデータを並び替えてし
まうとか，データを消してしまう心配はありません．ただ眺めてみたい，というばあいに
は Data Browser を使いましょう．
6.3
Do ファイルを利用しよう
これまで時々登場してきた Do ファイルというのも，慣れてくれば使うとよいでしょう．
Stata は基本的には対話型のパッケージですが，Do ファイルを使えばバッチ形式で処理を
行うことができます．最尤法のように一連のコマンドが長くなってしまう場合には，非常
に有効です．
Do ファイルは，単にコマンドを並べていったテキストファイルです．Stata には専用の
Do file editor がついていて，Do file を部分的に実行させるなどの機能を持っているので，
使ってみるとよいでしょう．ただし，editor としてはそれほど多機能であるわけではない
ので，好みによるところでしょう．
Do ファイルはコマンドを並べたものに過ぎませんが，Do ファイルから Do ファイルを
呼び出すこともできます．いくつかのコマンドをまとめて実行するという行為をまとめて
おこなうことができるようになります．データの加工・推定をべつべつの Do ファイルに
しておき，それらを他の Do ファイルから呼び出すようにすれば，ひとつひとつの Do ファ
イルが小さくなるため，編集が楽になるかもしれません．
しかし，単にコマンドを並べておくだけでは，あとから見直したときに解読できないか
もしれません．そのときには「/* */」でコメントをはさんでおきましょう．Stata は「/*
*/」で囲まれた部分は無視してくれます．これを応用すると，長いコマンドを改行しなが
ら書いていくことができます．下の例では，最後の２行はひとつのコマンドとして認識さ
れることになります．
set mem 30m
set more off
/* Converting sick report to Stata format */
quietly infix using "C:\My Documents\2001\Ii\Intage\sicka.dct"
save "C:\My Documents\2001\Ii\Intage\sicka.dta", replace
/* generating data about 1st time to go to a doctor */
generate godoc1st = 1 if cope1n1 == 1
quietly replace godoc1st = d1d + 1 /*
*/ if cope1n1 == . & cope2n1 == 1
ではいよいよ．
8
簡単なコマンドを実行してみよう
7
どんな分析を行うときにも，必ずといっていいほど使うコマンドをいくつか紹介しましょ
う9 ．
標本統計量を求めよう
7.1
変数の一覧表にはいくつかの形式があります．たとえばdescribe コマンド，codebook
コマンド，list コマンド，inspect コマンドなどがあります．それぞれに出力形式や内容
が異なるので，必要に応じて使い分けるとよいでしょう．
標本統計量といってもいくつかあります．平均・分散・最小値・最大値を出力するには，
sum コマンドを使います．
sum hosd1-hosd30
とすれば，hosd1, hosd2, ..., hosd30 までの平均・分散・最小値・最大値を出力し
ます．中間値などを求めるときには，centile コマンドを使います．
変数が離散的（いくつかの決まった値しか取らない）ときには，table コマンドを用い
ると，各値の頻度を得ることができます．table のあとに，変数名を２つ並べるとクロス
表を書くことができます．また，相対頻度・累積頻度を求めるには，tabulate コマンドを，
相関係数はcorrelate コマンドを利用します．
サンプルを限定しよう
7.2
標本統計量を求めるコマンドに限った話ではありませんが，サンプルを限定して分析を
行いたい場合があります．分析にまったく必要ないことが確信を持って言える場合には，そ
のようなサンプルを落とす（drop コマンド）方法もありますが，そうでなければ，ある条
件を満たすサンプルに限ってさしあたりコマンドを実行するという方法があります．たと
えば，
table
duration godoc1st if
age >= 40 & duration < 11
とすれば，変数age の値が 40 以上，変数duration の値が 11 未満のサンプルに限って，
変数duration, godoc1st のクロス表を書け，ということになります．if 以後で用いる演
算子としては，==（等しい）， >=（以上）， <=（以下），>（より大きい）， <（より小さ
い）， ~=（等しくない）， &（かつ）， |（または）， +（和）， -（差）， /（商）， *（積），
^（べき乗）等があります．また，欠損値. は非常に大きな数とみなされるので，大小関係
で条件文を作るときには注意しましょう．
さて，条件を決めてコマンドを実行するのもよいのですが，これを繰り返したいとおも
うことがあります．たとえば，離散変数の値に応じて平均値を求めたいとおもうこともあ
9 ぼく自身があまり使用しないので，グラフを描く話はしません．Help かマニュアルでgraph を調べてくださ
い．
9
るでしょう．そのようなとき，if 条件文だけだと，次のように長々と書かなければなりま
せん．
sum duration if sickname == 1
sum duration if sickname == 2
...
これはできないことはありませんが，手がかかります．このように，離散変数に応じて
コマンドを実行するときには，by オプションを用います．上の例なら，
gsort sickname
by sickname : sum duration
とすれば，変数sickname の値ごとにコマンドを実行してくれます．もし，by のあとに
連続変数が来ると，その値ごとにコマンドを実行するので非常に長い出力になってしまい，
往々にして意味を持ちません．
by オプションは，データがその離散変数によってソートされていなければ実行してくれ
ません．先ほどの例の１行目は，データを並べ替えるコマンドです．
条件文を連続変数で定義して，コマンドを次々繰り返したいこともあるでしょう．たと
えば，年齢を 10 歳幅でとって，それぞれに平均値を出してみたいとします．するとコマン
ドラインは，
sum duration if age <= 9 & age >=0
sum duration if age <= 19 & age >=10
...
となります．このような繰り返しもいちいち書いていくのはしんどいだけです．このよ
うな場合には，for オプションを使います．上の例なら，
for num 0/9 : sum duration if age <= X9 & age >=X0
と書くことで同じコマンドを実行してくれます．これは，文字X のところに，そのまま，
０から９までの整数を順次代入してコマンドを実行する，というオプションです．整数の
代入先はデフォルトではX ですが，指定することもできます10 ．また，代入先が２つ以上
ある場合には，順にX, Y, Z となります．
by やfor オプションで複数のコマンドを実行したいこともあるでしょう．そのときには，
コマンドを円マークで区切っておきます．これらを応用すると次のようなコマンドを書く
こともできます．
for num 2/29 \ num 3/30 : gsort hosdcX \ by hosdcX : sum hosdY
このコマンドは以下を実行するのと同じことになります．
10 このようなケースを考えると，変数名にうかつに大文字を使うのは危険かもしれません．
10
gsort hosdc2
by hosdc2 : sum hosd3
gsort hosdc3
by hosdc3 : sum hosd4
...
gsort hosdc29
by hosdc29 : sum hosd30
by やfor オプションを使うと，出力結果が長くなってしまうので，ログを取っておくの
を忘れないようにしましょう．
変数を加工しよう
7.3
sum 等のコマンドの引数は変数名です．したがって，ある変数を対数変換したものの平
均，とか，２つの変数の比の平均，とかいったものはそのままでは求めることができませ
ん．そこで，変数を加工する必要があります．たとえば，２つの変数の和を求めたいとき
には，
generate hosdc2 = hosd1 + hosd2
とします．このとき，generate コマンドの引数に，すでにある変数名を指定するとエ
ラーがでます．また，変数を作るもとになる変数のどれかが欠損値であるサンプルや，ゼロ
で割ることになるサンプルには，欠損値を返します．自然対数はlog で求めます．generate
コマンドには使えないいくつかの関数を使うことができるのがegen コマンドです．egen
コマンドはmax, min, sum, median といった関数を含んだ表現を右辺に使うことができ
ます．
パネルデータで成長率データを作りたいといった場合には，
「前」の観測値を含んだ計算
が必要になります．この場合には，変数名に [ n − 1] といった表現を追加します．このばあ
い，順序はそのときの並び順に依存するので注意してください11 ．たとえば，パネルデー
タである主体についてはあいだが 1 回分抜けているといったことには，とくに注意を払わ
ないと自動的に対応してくれるわけではありません．具体的なコマンドラインは以下のよ
うになります．
sort id time
by id: generate dhos = hos[_n] + hos[_n-1]
変数hos の 1 階の差分系列を，id というグループごとに，dhos として作成しています．
差分の将来値がほしいときには，
sort id time
by id: generate dhos = hos[_n+1] + hos[_n]
11 ts
で始まるコマンド群は，時系列に対応しているので，最初に宣言していれば気にしないでも構わないよう
ですが，ここでは触れません．
11
となります．
変数を加工するとき，他の変数の状況によって作り方を変えたいときがあります．ある
いは，取り込んだデータでは欠損値になっているところをゼロで代理させたいなんてこと
もあるかもしれません．欠損値があると，計算から外れてしまいます．
そのようなときには，replace コマンドを使います．
generate avedoc = 0 if godoc1st ==.
for num 2/30 : quietly replace avedoc = hosdcX / duration if duration == X
quietly replace avedoc = 1 if duration == 1 & godoc1st == 1
変数avedoc により，病気が治るまでの日数のうちに病院に行った日数の比率を求めよう
としています．最初に病院に行くまでの日数godoc1st が欠損値なら，その比率はゼロです
（１行目）．治るまでの日数duration がX なら，X 日目までに通院した累積日数hosdcX を，
治るまでの日数duration で除したものがその比率になります（２行目）．また，１日目に
病院に行ってその日のうちに治ってしまえば，比率は１になります（３行目）．
既存の変数を用いて新たに変数を作る場合，しばしば登場するのがダミー変数です．条件
を満たせば１，満たさなければ０の値を取る変数で，たとえば以下のようにして作ります．
generate inc700d = 1 if income >= 700 & income ~=.
replace inc700d = 0 if income < 700
このようなダミー変数を，離散変数に対応して作りたいばあいがあります．たとえば居
住県がコード化された変数pref があるとき，その絶対値にはそれほど意味はありません．
このようなとき，ある県に住んでいるかどうかでダミー変数を県の数だけ作ります．もち
ろん，generate コマンドとreplace コマンドを組み合わせてもよいのですが，この場合
には以下のようにして作成できます．
tabulate pref, generate(pref)
こうすると，ダミー変数pref1, pref2, ... が取りうる値の数だけ発生します．このと
き，作成されたダミー変数の label は「pref == 1」等の注釈がつきます．
データセットを加工しよう
7.4
異なるデータセットにある変数を一度に推定に取り込んだりすることは，Stata ではな
かなかできません．あらかじめデータを合体させておくなどの処理が必要になります．あ
るいは，データの並べ替えなどの際にも，単なるソート以上のことをする必要があるかも
しれません．個人ベースの元データを家計ごとに作り変えるなどのケースがこれにあたり
ます．
どちらにしてもよく使うのは，変数や観測値を除去するコマンドです．
keep 変数名
drop 変数名
12
keep if 条件式
drop if 条件式
keep は変数・観測値で，条件に合うもの・並べられた変数以外の変数を全て除去します．
drop は逆に，変数・観測値で，条件に合うもの・並べられた変数を全て除去します．これ
らのコマンドで除去されたものは復元はできませんが，大量のデータを処理する場合など，
不必要な変数や観測値を残しておくと，処理速度が遅くなることがあります．
2 つの異なるデータセットを合体させるには，append コマンドかmerge コマンドを利用
します．いずれも，片方のデータセットを開いておいて，追加するデータセットをusing
オプションで指定します．追加するデータセットも Stata 形式にしておいたほうが無難で
す．append コマンドは観測値を追加しますから，同じ変数名のものがあれば，その下に
データが追加されます．どちらか片方にしかない変数については，もう片方の該当部分に
は missing variable が入ります．
append using newdata
といった形式で使います．using オプションのあとは，ダブルクォーテーション"でくくっ
てフルパスを指定することもできます．
merge コマンドは変数を追加します．開いてあるデータセットとの合体には，観測値を
マッチングさせる変数を指定する必要があります．いくつかの変数を組み合わせることで
マッチングさせることもできます．
merge recid using ds2
この場合，マッチングの鍵となる変数はrecid，追加されるデータセットはds2 となりま
す．マッチングの鍵となる変数で両方のデータセットがソートされていないと，merge コ
マンドによるマッチングが行われないので注意しましょう．merge コマンドを使うと，自
動的に merge という名前の変数が作成されます．この変数は，それぞれの観測値がどちら
のデータに属していたか，すなわち，元データのみに入っていた，新データのみに入って
いた，両方に入っていた，といった情報を示します．この変数はdrop コマンドで削除しな
い限り残っているので，データセットを次々に追加していくときには，新たな merge を定
義できなくなってエラーが返ってくることになってしまいます．
これらのコマンド群を使うと，個人ベースのデータを世帯ベースに直す，といったこと
が可能になります．たとえばいま，個人ベースの夫婦のデータセットが手に入ったとしま
しょう．すると，もとのデータから夫のデータだけ抽出（keep）→保存（save）→もとの
データを呼び出して（use）妻のデータだけ抽出（keep）→保存（save）→再び夫のデー
タを呼び出し（use）妻のデータを追加（merge）という手順を取れば世帯ベースのデータ
セットに作り直すことができます．こういった手続きを踏むときの注意点としては，merge
コマンドを実行するときに，元データと新データに同じ名前の変数名があると，新データ
のほうの変数が取り込まれないという点があります．そのため，変数名をちょっと変える
（rename var）などの手続きを途中で加えておく必要があります．もちろん，ほぼ逆の手
続きをとることで，世帯ベースを個人ベースに直す，といった手続きも行うことができる
でしょう．
13
一発コマンドで推定しよう
8
では，さらに分析らしい分析に取り掛かりましょう．とはいっても，よくある分析手法を
用いる限り，コマンドも簡単です．たとえば，平均値の差の検定はttest, by(variables)
と打てばできてしまいます．また，ほとんどの回帰分析は，
コマンド名・被説明変数・説明変数
の順に並べて書いてリターンキーを押せば結果が出力されます．それゆえ，われわれの
やるべきことは，
「Help →search」で行いたい分析手法を検索し，コマンドのフォーマッ
トを調べることだけです．分析手法によっては，オプションを指定する必要がありますが，
それも Reference Manual に丁寧に書いてあります．なので，ここからは，使いたい分析手
法に応じて Reference Manual と首っ引き，ということになります．
8.1
OLS をやってみよう
さきに述べたとおり，コマンド名・被説明変数・説明変数を順に並べれば回帰分析を行う
ことができます．最も簡単な通常の最小二乗法で，以下の式を推定することを考えましょう．
durationi = β0 + β1 gender + β2 inc700d + εi
このとき，Gauss-Markov の定理が成立して，係数 β0 , β1 , β2 の推定量が BLUE となるた
めには，誤差項が強外生性を満たしていなくてはなりません．つまり，一般化して書いて
おけば以下のとおりです．
定理 1 いま，推定式が真のモデルであり，標本（サンプル）が以下の仮定を満たすとする．
1. 推定式の線形性：推定式は次式のように線形である．
yi = β1 xi1 + β2 xi2 + · · · + βK xiK + εi
2. 誤差項の強外生性：誤差項は全ての観測値という条件付きの期待値がゼロ．
E[εi |X] = E[εi |x1 , x2 , ..., xn ] = 0
3. 多重共線性がない．すなわち，n × K のデータマトリックスのランクは確率１で K ．
4. 誤差項の分散が一定で，かつ，共分散がゼロ．
E[ε2i |x1 , x2 , ..., xn ] = σ 2 > 0
E[εi εj |x1 , x2 , ..., xn ] = 0
(i = 1, 2, ..., n)
(i = j)
このとき，OLS 推定量 b = (X X)−1 X y と真の係数ベクトル β は以下の性質を持つ．
1. 不偏性（unbiasedness）12 ：E[b|X] = β
2. 推定値の分散：V ar(b|X) = σ 2 (X X)−1
12 この性質が成立するためには誤差項の分散の仮定は必要ありません．
14
3. 最良線形不偏推定量（BLUE : Best Linear Unbiased Estimator）：b は線形不偏推定
量のなかで最も有効（efficient），すなわち最も小さい分散を持つ．
4. 残差との共分散がゼロ：cov(b, e|X) = cov(b, (y − Xb)|X) = 0
5. 推定分散の不偏性：E[e e/(n − K)|X] = σ 2
証明の一部は次の小節にとりあえず載せておきます．さて，一般には誤差項が上の定理の
仮定を十分に満たしていることはほとんどないわけですが，それでも OLS の手法をデータ
にとりあえず当てはめてみるということは，理由はともかくとして，しばしば行われます．
Stata で OLS 推定を行う場合には，reg コマンドを使います．さきほどの式を推定する場
合のコマンドラインは，
reg duration gender inc700d
となります．ここではduration が被説明変数，それ以降が説明変数です．reg コマンド
にもいくつかのオプションがあります．説明変数には自動的に定数項が含まれてしまうの
で，それを外すときにはnocons オプションを指定し，
reg duration gender inc700d, nocons
とします．標準偏差の推定に White の頑健推定量を用いるばあいには，
reg duration gender inc700d, robust
とオプションをつければよいことになります．不均一分散の修正方法としては，Davidson
and MacKinnon (1993, pp.554-556)13にある方法も用いることができ，オプションhc2 やhc3
を指定して使うことができます．誤差項に 1 階の系列相関を仮定する，Cochrane-Orcutt
による方法を用いるときは，prais コマンドにcorc オプションをつけて推定します．係数
制約の検定や非線形最小 2 乗法については後で述べます．
8.1.1
Gauss-Markov の定理の証明
説明変数が k − 1 個，サンプルの大きさ14 が t のモデルを考えてみましょう．さしあたっ
て説明変数の選択の問題は解決されているとして，定数項を明示的に取り扱うと次の線形
連立方程式で表現されます．
Y1 = β1 + β2 X21 + β3 X31 + · · · + βk Xk1 + u1
Y2 = β2 + β2 X22 + β3 X32 + · · · + βk Xk2 + u2
..
.
Yt = β1 + β2 X2t + β3 X3t + · · · + βk Xkt + ut
(1)
13 Davidson, Russell and James G. MacKinnon. 1993. Estimation and Inference in Econometrics.
Oxford University Press.
14 それぞれのデータのことは「観測値」といい，その集合のことを「サンプル」というらしいので，
「サンプル
の数」という表現はヘンだ，という妙なこだわりを持つ人も世の中にはいるらしいです．
15
これを行列で表現したいので，式をじーっと見ると，右辺の第 k 項目までは係数が共通で
あることに気がつくでしょう．そこで，
Y1
1 X21 X31
Y2
1 X22 X32
y=
..
..
.. , X = ..
.
.
.
.
Yt
1
X2t
X3t
···
Xk1
···
..
.
Xk2
..
.
β2
,β = .
.
.
···
Xkt
βk
β1
,u =
u1
u2
..
.
(2)
ut
とおいて，
y = Xβ + u
(3)
とすれば，もとの式になります．ここで，縦ベクトル y は t × 1，X は t × k ，縦ベクトル
β は k × 1，縦ベクトル u は t × 1，ですから，行列の積である Xβ の次元15 は t × 1 となり
ます．行列のふつうの積は「ヨコ×タテ」で定義されますから，1 行目がちゃんと
Y1 = β1 + β2 X21 + β3 X31 + · · · + βk Xk1 + u1
(4)
となっていることを確認しましょう．
さてここで，
（3）についていくつか確認しておきましょう．X は確率変数ではないが，u
のそれぞれの要素は実現値を観測できない確率変数です．観測できない β は定数ですが，
u が確率変数ベクトルなので y のそれぞれの要素も確率変数となります．ここで，
（3）と
いう式の形と，u についていくつかの仮定を置けば，それっぽい β を決めることができる，
というのが推定の手続きに他なりません．最小 2 乗法とは，u のそれぞれの要素はあまり
大きくないと想定して，その 2 乗和を最小にするような β をもって推定値にする，という
方法です．さてそこで，
（3）を用いて最小にされるべき 2 乗和を計算してましょう．u を移
項すると
u = y − Xβ
(5)
となります．2 乗和は
u21 + u22 + · · · + u2t = u u =
u1
u2
···
ut
u1
u2
..
.
ut
(6)
と書けます16 から，右辺も同じように書いてみると，
u u = (y − Xβ) (y − Xβ)
(7)
15 行列の縦にいくつ，横にいくつの要素が並んでいるか，ということを次元といいます（ルパン 3 世とは関係
ない）．ランクとは違うものです．次元の違う行列は足したり引いたりできないですし，一般の積の場合は左の
行列の横の次元と，右の行列の縦の次元が同じでないと計算できません．だから，行列のふつうの掛け算では交
換則が成り立ちません．掛けられる行列のそれぞれの要素に，掛ける行列を掛けるという，積となる行列の次元
がやたらと大きくなる掛け算も考えることがあって，そういうのをクロネッカー積といいますが，さしあたって
使わないので気にしないことにしましょう．行列では割り算は定義されません．
16 u は行列のタテとヨコをひっくりかえしたもので転置行列といいます．数学の人たちは uT とも書くことが
あるようです．どっちにしても微分や累乗とまぎらわしい表現ですね．
16
和をとったものの転置行列は転置行列の和になるから17 ，
u u = (y − (Xβ) )(y − Xβ)
(8)
積の転置行列は転置したものの積になるから18 ，
u u = (y − β X )(y − Xβ)
(9)
カッコを外す分配則はふつうに成り立つから19 ，
u u = y y − y Xβ − β X y − β X Xβ
(10)
3 項目を β X と y の積とみると，積の転置行列は転置したものの積だから，
u u = y y − y Xβ − y Xβ − β X Xβ
(11)
2 項目と 3 項目をまとめると
u u = y y − 2y Xβ − β X Xβ
(12)
となります．適当に β の推定値 b を決めてやれば，それに応じた e = y − Xb を計算する
ことができるので，そこから計算された e e のことを残差平方和とよびます．つまり，残
差平方和を最小にする β を見つけるのが最小 2 乗法ということになります．残差平方和は
推定値 b の関数ですから，残差平方和を最小にする b を見つけるには，残差平方和を b で
微分してゼロと置いて解けばよいわけです20 ．
「残差平方和を b で微分して」といっても，b はベクトルですからそのまま微分を考える
ことができるわけではありません．ベクトルで微分するというと難しそうですが，ベクト
ルを構成するそれぞれの要素でそれぞれ微分して，それを重ねて書き直したもののことを，
ベクトルで微分するといっているだけのことです．つまり，関数 f (x) が (x1 , x2 , · · · , xk )
の関数であるとき21 ，ベクトル x = (x1 , x2 , · · · , xk ) で微分したものは，
∂f
∂f
=
∂x
∂x1
∂f
∂x2
..
.
(13)
∂f
∂xk
です．もし，関数 f (x) が (x1 , x2 , · · · , xk ) の線形和であれば，
x1
x2
f (x) = a x = a1 a2 · · · ak .
..
(14)
xk
17 ただの足し算で位置が違うだけだからです．
18 積を取るときの縦横をひっくりかえさないといけないからです．確認してみましょう．
19 (a
+ b)(c + d) = ac + ad + bc + bd が成り立つということ．
20 最大値を見つけることになるかもしれないので，ほんとうはもうちょっとチェックが必要．
21 関数の値がベクトルではなくてただの数になるケースのみを考えます．ベクトル値関数の場合は横にも並べ
ればよいだけですが．
17
と表現できて，その微分は
a1
a2
∂f
=
.
=a
∂x
..
(15)
ak
となるので，まとめて書けば，
∂(a x)
=a
∂x
(16)
となります．2 次関数のばあいも同じようにこちゃこちゃと計算すれば，行列 A が対称行
列22 なら，
∂(x Ax)
= 2Ax
∂x
(17)
となることが分かっています．さて，これだけの準備をして
u u = y y − 2y Xβ − β X Xβ
(18)
を β で微分することを考えましょう．最小 2 乗法で求められる推定値 b は微分をゼロにす
る23 から，
∂(u u)
= −2y X + 2X Xb = 0
∂b
(19)
が成り立ちます．転置を書き直して移項して 2 で割れば，正規方程式を得ます．
X Xb = X y
(20)
さて，行列では割り算が定義されていませんから，正規方程式から b を求めるには，
X Xb = X y
(21)
の両辺に左から X X の逆行列を掛ければよいことになります．ところが，逆行列という
のはどんな行列に対しても定義されているわけではないことに注意しましょう．行列を構
成するベクトルが 1 次独立でなければ逆行列は存在せず，その場合にその行列はランク落
ちしている…のですが，さしあたっては，データを入れた行列 X が多重共線性をもってい
れば逆行列が定義できず，あるいは多重共線に近い状況であれば逆行列の値が大きくぶれ
る，ということだけをおさえておきましょう．つまり，完全な多重共線が発生していれば，
OLS 推定量は定義できないことになります．
さて，多重共線が発生していなければ，両辺に左から X X の逆行列を掛けて，OLS 推
定量を計算することができます．すなわち，
b = (X X)−1 X y
(22)
22 タテとヨコの次元が同じ行列（正方行列）のうち，転置しても同じ行列のことを対称行列といいます．分散
共分散行列は対称行列ですし，任意の行列 X に対して，X X は対称行列になります．
23 各要素がすべてゼロになるということ．
18
です．
このように求められた OLS 推定量と真の値 β との差を考えてみましょう．まず定義に
より，
b − β = (X X)−1 X y − β
(23)
y は定義により Xβ + u に等しいから，
b − β = (X X)−1 X (Xβ + u) − β
(24)
分配則をつかってばらばらにすると，
b − β = (X X)−1 (X X)β + (X X)−1 X u − β
(25)
逆行列の定義 X −1 X = I より，
b − β = β + (X X)−1 X u − β
(26)
最初の項と最後の項はおなじだから，
b − β = (X X)−1 X u
(27)
となります．これを sampling error と呼びます．
ここで，両辺の期待値をとってみましょ．
E[b − β] = E[(X X)−1 X u]
(28)
期待値には線形性があるから，確率変数ではない係数は期待値の外に出て24 ，
E[b] − β = (X X)−1 X E[u]
(29)
もし，真の誤差項 u の全ての要素の期待値がゼロならば，右辺はゼロになるから，移項して，
E[b] = β
(30)
が成り立つことになります．つまり，OLS 推定量は誤差項の期待値がゼロでありさえすれ
ば不偏性を持つことがわかります．
次に，係数推定値の検定のために，OLS 推定量 b の分散共分散行列25 を求めてみましょ
う．まず，確率変数の分散は，定数を足しても引いても変化しないから，
var(b) = var(b − β)
(31)
これは sampling error に他ならないから，
var(b) = var((X X)−1 X u)
(32)
定数行列 A と確率変数ベクトル x が与えられたときに，共分散行列について
var(Ax) = Avar(x)A
24 OLS
(33)
推定量 b は確率変数ですが，真の値 β は確率変数ではないことに注意しよう．
25 確率変数のベクトルが与えられたとき，それぞれの要素の確率変数の間の共分散を要素とする行列のことを
共分散行列といいます．k 個の確率変数があれば，共分散行列は k × k であり，対称行列となります．
19
が成り立つことが知られている26 から，
var(b) = (X X)−1 X var(u)((X X)−1 X )
(34)
もし，真の誤差項 u の全ての要素の期待値がゼロならば，共分散行列は
var(u) = E(uu )
(35)
と表現できるから，
var(b) = (X X)−1 X E(uu )((X X)−1 X )
(36)
共分散行列が単位行列の σ 2 倍で表現できる，つまり，分散が均一で相関がないとすれば，
var(b) = (X X)−1 X (σ2 In )((X X)−1 X )
(37)
となります．σ 2 はただの実数だから，
var(b) = σ 2 (X X)−1 X ((X X)−1 X )
(38)
転置行列の性質より27 ，
var(b) = σ 2 (X X)−1 X X(X X)−1 = σ 2 (X X)−1
(39)
それゆえ，ここまでをまとめてみると，
定理 2 推定式が線形で，誤差項の期待値がゼロで，データが多重共線を起こしていないと
き，OLS 推定量は不偏性を持つ．
定理 3 推定式が線形で，誤差項の期待値がゼロで，データが多重共線を起こしておらず，
系列相関も分散不均一もないとき，OLS 推定量の共分散行列は σ 2 (X X)−1 で求められる．
2 番目の定理の仮定が成り立つとき，OLS 推定量は，データを線形に組み合わせたよう
な推定量の中でもっとも分散が小さいことがわかっています．これが Gauss-Markov の定
理です．また，OLS 残差 e と OLS 推定量 b の共分散はゼロになります．
8.2
IV を使ってみよう
誤差項が強外生性を満たさないとき，操作変数法が用いられることがあります．操作変
数法は一般化積率法（GMM: Generalized Method of Moments）の特殊形と位置づけられ
ます（Hayashi 2000, pp.226-231 ）28 が，推定の考え方は次のようなものになります．いま，
推定式を
yi = δ1 z1i + δ2 z2i + · · · + δL zLi + εi
26 定数倍された確率変数の分散は，もとの確率変数の分散に対して，定数の
2 乗倍になるから．
= B−1 A−1 が成り立ちます．ABB−1 A−1 = I だから．
28 Hayashi, Fumio. 2000. Econometrics. Princeton University Press.
27 ちなみに，積の逆行列については，(AB)−1
20
としましょう．推定したい係数は δ1 , δ2 , · · · , δL であり，全部で L 個あります．説明変数
z1i , z2i , · · · , zLi は誤差項 εi と相関を持つかも知れず，強外生性を満たさないとしましょう．
いま，誤差項とは相関を持たず，しかし説明変数と相関を持つような変数が K 個あると
し，それらを x1i , x2i , · · · , xKi とします．このような変数を操作変数（IV: Instrumental
Variables）と呼びます．操作変数は誤差項と相関を持たないので，次の式が成り立ちます．
0
x1i
x2i 0
. εi = .
E
. .
. .
xKi
0
誤差項 εi は εi = yi − δ1 z1i − δ2 z2i − · · · − δL zLi と書くことができるので，代入すると，
x1i
0
x2i
0
. [yi − δ1 z1i − δ2 z2i − · · · − δL zLi ] = .
E
.
.
.
.
xKi
0
ちょっと変形して，
x1i (yi − δ1 z1i − δ2 z2i − · · · − δL zLi )
x2i (yi − δ1 z1i − δ2 z2i − · · · − δL zLi )
E
..
.
xKi (yi − δ1 z1i − δ2 z2i − · · · − δL zLi )
=
0
0
..
.
0
このような条件のことを直交条件（orthogonality condition）と呼びます．
この式を標本平均で書き直したものを考えると，xki , zli , yi はデータとして手許にある
わけですから，未知数 (δ1 , δ2 , · · · , δL ) が L 個，式が K 個あることになります．したがっ
て，代数の基本定理からわかるように，式の数（=操作変数の数）が未知数の数（=説明変
数の数）より多いとき，すなわち K < L のときには，未知数を求めることができません
（Not identified）．操作変数の数と説明変数の数が同じとき（Just identified）のときには未
知数を 1 つだけ求めることができます．操作変数の数が説明変数の数より多いとき（Over
identified）のときには，式の数のほうが多いので，うまくやらないと求めることができま
せん．Over identified のときに，直交条件が成り立っているかどうかを検定する統計量を，
GMM ではHansen の J 統計量（Hansen’s J-statistic）（Hayashi 2000, pp.217-218 ）
といい，特定化検定（specification test）の 1 種となります．
したがって，GMM 推定を行うためには，すくなくとも説明変数の数よりも多くの操作
変数を用意することが必要になります（order condition）．もちろん，説明変数のなかに
直交条件を満たすものがあれば，それを操作変数として用いることができます．説明変数
がすべて直交条件を満たすならば，わざわざ IV を用いなくても，OLS で推定できること
になります．すなわち，説明変数のみを操作変数とし，誤差項が均一分散の仮定を満たす
とき，GMM 推定量と OLS 推定量は同値です．逆に，説明変数のうちに直交条件を満たさ
21
ないものがあるときには，その変数（しばしば内生変数とよばれます29 ）の数と同じかそ
れ以上の，説明変数に含まれないような操作変数を用意する必要があります．このような
操作変数を除外変数（excluded variables）と呼んだり，除外変数の直交条件のことを除外
制約（excluded restrictions）と呼ぶこともあります．
また，OLS 推定の場合と同様，操作変数は多重共線を起こしてはいけません（rank con-
dition）．
誤差項が操作変数の条件付きで分散一定で系列相関を持たない均一分散性を満たすとき
（conditional homoskedasticity），efficient GMM 推定量は2 段階最小二乗（2SLS）推
定量と同値になることが知られています（Hayashi 2000. pp.226-231 ）．2SLS 推定とは，内
生変数を操作変数で OLS 推定し（1 段階目），その当てはめ値 zˆli（predicted values, fitted
values）を内生変数に代入して再度推定を行う（2 段階目）方法をいいます．2 段階目の推
定で求められる残差は yi − ˆ
zi δˆ ですが，統計的検定を行うときにもちいる残差は yi − zi δˆ
でなければならない点には注意が必要です．
さて，
「操作変数法」という場合には，狭く 2SLS 推定をさすことが多いようです．ここ
で説明した GMM も操作変数を用いるので操作変数法の 1 種ではあるのですが，GMM と
いった場合には，より複雑な直交条件を利用して推定を行う方法を広く指すようです．こ
れまでに述べてきた線形の 2SLS 推定は，Stata ではivreg コマンドで行うことができま
す．コマンドラインは，
ivreg 被説明変数 外生説明変数 ( 内生変数 = 除外された操作変数)
求められた残差に White の修正を加える場合には，OLS 推定のときと同様，robust オ
プションを使いますし，定数項を外すnocons オプションも有効です．また，1 段階目の推
定結果を表示させるにはfirst オプションを指定します．
2SLS 推定量は一致性（consistency）を持つのですが，有限標本（finite-sample）では必
ずしも望ましい性質を持つとは限りませんし，有限標本では不偏推定量ではありません．ま
た，操作変数が内生変数と弱い相関しか持たないとき，有限標本での近似がよくないこと
が知られています（weak instruments の問題; Hayashi 2000, p.229, Staiger and Stock
1997）30 ．そのため，1 段階目の F 統計量や，R2 を報告する傾向があるようです．
ユーザ作成の ado ファイルを使うと，線形 GMM は以下のようにして推定できます．
ivgmm0 被説明変数 外生説明変数 ( 内生変数 = 除外された操作変数)
8.3
Probit 分析をしてみよう
Stata は個票の処理に強みを持っているのですが，被説明変数が連続でないケースを検
討することもあるでしょう．被説明変数が幾つかの限られた値を取り，またその値自身に
29 計量経済学的な意味での内生性とは，単に，誤差項と説明変数が相関を持つことをさします（Wooldridge
2001, pp.49-51.）．もちろんそのなかには経済学的な意味での同時性（simultaneity）を原因とするものもあり
ますが，適切な説明変数の省略（omitted variables）や，測定誤差（measurement error）に起因するも
のもあり，経済学的な意味で同時性がないからといって，計量経済学的な意味での内生性がないとは限りません．
30 Staiger, D. and J. Stock. 1997. Instrumental variables regression with weak instruments. Econometrica 65, 557-586.
22
はとくに意味のないような回帰モデルのことを質的変量モデルといいます．たとえば被説
明変数が「Yes」
「No」の２種類の値を取るようなケースがこれに当たり，とりうる値が２
個の被説明変数をとくに２値変数とも呼びます．
質的変量モデルの推定にあたっては，連続変数としての潜在変数（latent variable）を
想定します．すなわち，潜在変数がある範囲の値を取れば質的変数がある値を取ると考え
ます．潜在変数と説明変数の関係には，最も簡単なケースでは線形の関係を仮定します．
いま，説明変数が定数とひとつの変数 xi であるケースを考えてみましょう．観測される
質的変数を yi ，潜在変数を yi∗ として，潜在変数と説明変数のあいだに線形の関係を仮定
すれば，
yi∗ = β0 + β1 xi + εi
と表現できます．ただし，β0 , β1 は推定されるべき係数であり，プロビットモデルでは，εi
は互いに独立で同一の正規分布 N (0, σ 2 ) に従う誤差項であると仮定します．潜在変数があ
る範囲の値を取れば質的変数がある値を取ると考えますが，ここで質的変数が２値変数で
あるとすれば，潜在変数が正の時には質的変数が１，負の時には０の値を取ると仮定して
一般性を失いません．すなわち，
yi = 1
y ∗ > 0 のとき
yi = 0
y ∗ ≤ 0 のとき
となります．潜在変数は可能性として任意の実数値を取りえますが，観測される値は離散
的（ここでは２つ）となります．
さて，このような設定のもとで，観測される変数が１あるいは０である確率を求めてみ
ましょう．それぞれの観測値について質的変数が観測される確率（＝尤度）を求めること
ができれば，最尤法によって係数の推定を行うことができるからです．まず，観測される
変数が１である確率は，
P (yi = 1)
= P r(y ∗ > 0) = P r(β0 + β1 xi + εi > 0)
= P r(εi > −(β0 + β1 xi ))
(40)
ここで，εi は正規分布 N (0, σ 2 ) に従うと仮定されているので，さらに変形することができ
ます．標準正規分布の累積分布関数を Φ(z) とおけば，
β0 + β1 xi
σ
P (yi = 1) = P r(εi < β0 + β1 xi ) = Φ
(41)
同様に，観測される変数が０である確率は，
P (yi = 0)
= P r(y ∗ ≤ 0) = P r(β0 + β1 xi + εi ≤ 0)
= P r(εi ≤ −(β0 + β1 xi ))
(42)
ここで，εi は正規分布 N (0, σ 2 ) に従うと仮定されているから，同様に，
P (yi = 0)
= P r(εi ≤ −(β0 + β1 xi ))
β0 + β1 xi
= Φ −
=1−Φ
σ
23
β0 + β1 xi
σ
(43)
観測値間の独立性を仮定すれば，全体の尤度は各観測値の尤度の積で表現することができ
るので，尤度 L は，
n
L =
Φ
i=1
n
ln L =
i=1
β0 + β1 xi
σ
yi∗ ln Φ
yi∗
β0 + β1 xi
σ
1−Φ
β0 + β1 xi
σ
1−yi∗
+ (1 − yi∗ ) ln 1 − Φ
β0 + β1 xi
σ
となります．ここで，関数形から分かるとおり，対数尤度を最大化する係数推定値 b0 , b1 の
値は誤差項の分散 σ に依存しないことになります．そこで，σ = 1 と置いて最大化を行い
ます．
さて，推定された係数 b0 , b1 の値はどのような意味を持っているのでしょうか．プロビッ
トモデルでは，
P (yi = 1) = Φ (b0 + b1 xi )
であるから，b0 , b1 の値が示す意味は分かりやすいとはいえないでしょう．むしろ，説明変
数が 1 単位大きくなったときの確率 P (yi = 1) の変化分を示したほうが分かりやすいとき
もあります．そこで，データの平均値 x
¯ でのこの確率の変分の大きさ，すなわち，
∂
P (yi = 1) = φ (b0 + b1 x¯) b1
∂x
¯
を推定結果として報告することもあります．この値のことをマージナル効果といい，マージ
ナル効果は確率の変化分であるので，その大きさは％ポイントで表します．たとえば，説明
変数に平均値を代入したときの確率 P (yi = 1) が 70 ％であるときにある説明変数のマージ
ナル効果が 5 ％ポイントである，とは，説明変数の値が 1 増加したときに，確率 P (yi = 1)
が 70 ％から 75 ％に変化する，ということをいいます．
さて，ここでは 2 値選択モデルのうち，プロビットモデルのみを取り上げてきました．プ
ロビットモデルでは確率分布として正規分布を用いてきましたが，ロジスティック分布を
用いる推定法もしばしば用いられ，ロジットモデルと呼ばれます．また，被説明変数が取
りうる値の数が 2 よりも大きいようなモデルを考えることもできます．たとえば，被説明
変数が「病院へ通う・市販薬を飲む・何もしない」という 3 つの選択枝を表しているケー
スがこれに当たります．3 つ以上の値を取りうる場合，それぞれの値に順序をつけること
が可能である場合と，そうでない場合があるでしょう．先程の「病院へ通う・市販薬を飲
む・何もしない」という選択の場合には，明確な順序を定めることは難しいと思われます．
このようなモデルに確率分布として正規分布を用いたものを多項プロビットモデルといい
ます．また，
「年間所得が 500 万円以下なら 0，500∼1000 万円なら 1，1000 万円以上なら
2」といった場合は，明確な順序が存在する質的変数モデルといえます．このようなモデル
で正規分布を利用するものは順序プロビットモデルと呼ばれます．これらのモデルについ
ても潜在変数を想定し，最尤法で推定を行います．
容易に想像がつくとはおもいますが，プロビットモデル・ロジットモデルは以下のよう
なコマンドラインで推定することができます．
24
probit 被説明変数 説明変数
logit 被説明変数 説明変数
マージナル効果を求める場合，プロビットモデルではprobit のところをdprobit とすれ
ば求めることができます．また，順序プロビットのばあいは．
oprobit depvar variables
となります．
8.4
Heckman の 2 段階推定をしてみよう
個票の場合，サンプル設計のためにいくつかの特定の主体の被説明変数を観測できない
ことがあります．古典的な例では，女性の賃金です．働いている女性の賃金は観測するこ
とができますが，働いていない女性の賃金は観測することができません．そのような場合
に，Heckman の 2 段階推定という手法を用いることがあります．いま被説明変数 yi は外
生性を満たす説明変数ベクトル xi と誤差項 ui の線形関数として表現できると仮定します．
すなわち，
yi = β xi + ui
と書けるとします．ただし，被説明変数 yi は条件
π zi + εi ≥ 0
が満たされなければ観測できないとしましょう．ここで，説明変数ベクトル zi は外生で，
その係数ベクトルを π ，誤差項を εi とおきます．また，誤差項 εi の分散を σε2 ，もとの推
定式の誤差項 ui との共分散を σuε とし，εi と ui が 2 変量正規分布に従っていると仮定し
ましょう．すなわち，
ui
εi
0
=N
0
,
σu2
σuε
σuε
.
σε2
さて，被説明変数 yi が入手できた観測値についてだけ最小２乗推定を行うと，Gauss-Markov
の定理が成立せず，推定量が一致性を持たない可能性があります．そこでまず，観測でき
たという条件付きの期待値を取ってみましょう．
E[yi |xi , zi , π zi + εi ≥ 0] = β xi + E[ui |εi ≥ −π zi ]
この条件付き期待値の期待誤差を ηi とおくと期待値をはずすことができて，π zi + εi ≥ 0
なる観測されたサンプルについて
yi = β xi + E[ui |εi ≥ −π zi ] + ηi .
ここで，E[ui |εi ≥ −π zi ] を求めるために，εi と ui の分散共分散行列をコレスキ分解
（Choleski decomposition）しましょう（Greene 1997, p.44, p.179）．
σε2
σuε
σuε
σu2
=
σε
σuε /σε
0
σu
σε
0
25
σuε /σε
σu
このコレスキ分解を用いると，互いに独立に標準正規分布 N (0, 1) に従う確率変数 x1 , x2
を使って εi と ui を表現することができ，
εi
ui
=
0
0
+
σε
0
x1
σuε /σε
σu
x2
=
σε x1
(σuε /σε )x1 + σu x2
これを代入すれば，期待値の線形性より，
σuε
x1 + σu x2 σε x1 ≥ −π zi
σε
E[ui |εi ≥ −π zi ] = E
=
σuε
−π zi
E x1 x1 ≥
σε
σε
標準正規分布の密度関数を φ(·)，分布関数を Φ(·) とおけば，密度関数は φ (m) = −mφ(m),
φ(m) = φ(−m) を満たしますから，
E[ui |εi < −π zi ] =
σuε
σε
∞
−π zi
σε
x1 φ(x1 )
dx1
1 − Φ(−π zi /σε )
= −
1
σuε
σε Φ(π zi /σε )
= −
1
σuε
σε Φ(π zi /σε )
1
σuε
φ
σε Φ(π zi /σε )
φ (π zi /σε )
= λ
Φ(π zi /σε )
= −
∞
−π zi
σε
−x1 φ(x1 ) dx1
∞
−π zi
σε
φ (x1 ) dx1
π zi
σε
ここで，λ = −σuε /σε とおいています．もとの式に代入すると，観測されたサンプルにつ
いて
yi = β xi + λ
φ (π zi /σε )
+ ηi
Φ(π zi /σε )
となります．期待誤差 ηi が well-behaved だと仮定すれば，２段階推定を用いることができ
ます．Heckman の２段階推定法（Heckit）とは，第１段階として Probit モデルを用いて
π を推定し，第２段階として π の推定量から φ(π zi /σε )/Φ(π zi /σε )（逆ミルズ比（Inverse
Mill’s Ratio））の一致推定量を求めて OLS 推定を行う方法です．この推定量は一致性
を持ちますが効率性を持たないので，最近では全体を最尤推定することのほうが多いよう
です．
ここまで扱ってきたモデルは，Type 2 Tobit モデルと呼ばれますが，これはType 5
Tobit モデルと呼ばれるより一般的なモデルの特殊形とみなすことができます（Amemiya
1995）．
yi
=
β xi + ui
β = β1
if
Pi = π zi + εi < 0
β = β2
if
Pi = π zi + εi ≥ 0
26
このモデルでは，Pi の値によって β がスイッチするので，switching regression とも呼
ばれます．あるいは，Pi によってレジームが選択されるとも解釈できるのでself-selection
model と呼ばれることもあります．さらに，Pi が観察できないときには未知レジームのス
イッチング回帰と呼ばれるようです．
さて，最初のほうに述べた Heckman の 2 段階推定法は，以下のようなコマンドライン
で推定することができます．
heckman 被説明変数 説明変数, select( 選択変数 )
上の表現を使えば，説明変数は xi ，選択変数は zi のことです．heckman コマンドは，デ
フォルトでは最尤法を用いた推定をするので，あえて Heckman の 2 段階推定を行う場合
には，twostep オプションを，select( ) のあとに加える必要があります．
8.5
パネルデータを使ってみよう
同一主体を時系列的に追ったデータセットをパネルデータといいます31 ．そのようなデー
タセットを用いるのも簡単です．パネルデータを扱う場合には，主体・時点ごとにひとつ
の観測値としてデータセットを構築しておくのが普通ですが，random effect model の複雑
なものなどを扱う場合には，主体ごとに１つの観測値として「横に長い」データセットに
しておくほうが扱いやすい場合もあるようです．ここでは，主体時点ごとにひとつの観測
値としてデータセットを構築しておくとしましょう．その場合，
xtdes, i(主体変数) t(時間変数)
xtsum, i(主体変数) t(時間変数)
とすることで，バランスパネルになっているかどうか，時点別平均値など，データのお
おまかな分布をみることができます．説明変数が全て強外生性を満たすときの固定効果モ
デル・変量効果モデルの推定は以下のようになります．
xtreg depvar variables, fe i(主体変数)
xtreg depvar variables, re i(主体変数)
また，変量効果モデルの推定の直後にxthaus と入力すると，固定効果モデルと変量効果
モデルについての Hausman 検定を行うことができます．
パネルデータを扱っているときに，説明変数に内生性の疑いがあり（説明変数と誤差項の
あいだに直交条件が成り立たない），操作変数法を併用したい場合には，xtivreg というコ
マンドを使います32 ．xtivreg は連立方程式アプローチを採っているわけではありません．
説明変数に強外生性（strict exogeneity）が仮定されないときの推定方法についてはさまざ
まな方法があるようですが33 ，xtivreg では，GLS random-effect model，Between-effects
31 Stata 8 の Reference には「Cross-Sectional Time-Series」と名づけられた分冊があり，ここでパネルデー
タを対象とするコマンド群の解説が行われています．
32 Stata 8 から導入されたものかもしれません．
33 パネルデータの扱い方については，Baltagi, Badi H. 2001. Econometric Analysis of Panel Data 2nd
edition. John Wiley and Sons. などを参照のこと．また，より実践的なものとして，応用計量経済学の観点か
らの良書，Wooldridge, Jeffrey M. 2001. Econometric Analysis of Cross Section and Panel Data. MIT
Press. の第 10 章，11 章を参照のこと．
27
model，Fixed-effects model，First-differenced estimator の 4 種類の推定を行うことがで
きます．コマンドラインはほぼ共通で，
xtivreg 被説明変数
強外生変数 ( 内生変数 = 操作変数 ), 推定方法オプション
となります．推定方法オプションの部分は，先に述べた順に，re be fe fd と書き，First-
differenced estimator 以外の推定を行うときには，xtreg コマンドと同じく，i(主体変数)
による変数指定が必要です．First-differenced estimator を用いるときは，tsset コマンド
による時間変数・主体変数指定を事前に行う必要があります．また，いずれの推定方法を
用いるにしても，first オプションを指定することで，1 段階目の推定結果を表示させるこ
とができます．さらに，GLS random-effect model を用いる場合，Baltagi による EC2SLS
をオプション指定することもできます．
説明変数に被説明変数の過去値が入っているモデルをダイナミックパネルと呼ぶことが
ありますが，Arellano and Bond (1991, RES) によるダイナミックパネルの推定はxtabond
コマンドで行うことができますし，説明変数と個別効果が相関を持つ場合の Amemiya and
McCurdy (1986, Econometrica) や Hausman and Taylor (1981, Econometrica) による推
定はxthtaylor コマンドに組み込まれています．さらに，パネルデータによる質的変数モ
デル（xtlogit ，xtprobit ）といったコマンド群も用意されています．Stata 8 以降は，い
くつかのパネルデータコマンドは，tsset コマンドによる事前の時間変数・主体変数指定
を必要としています．いずれにせよ，これらの詳細については，マニュアルを参照してく
ださい．
8.6
サバイバル分析をやってみよう
もっとも簡単なかたちでのサバイバル分析を行うためには，期間を表す変数と，打ち切
られたかどうかのダミー変数，時間に依存しないかたちの説明変数を用意します．まず，期
間を表す変数（ここではdurac），打ち切りダミーを表す変数（ここではfaildoc ）を指定
し，サバイバル分析を行うのだ，ということを Stata に認識させます．
stset durac, failure( faildoc)
あとはさきほどとほとんど同じです．ただしここでは，被説明変数を認識させているの
でこれを指定する必要はありません．たとえば Cox 比例ハザードモデルを推定するばあい
には
stcox gender hcon2-hcon6 inc700d copay2-copay7
ワイブル分布を仮定したパラメトリックな比例ハザードモデルでは，
streg gender hcon2-hcon6 inc700d copay2-copay7, d(weibull)
ノンパラメトリックに Kaplan-Meier 分析のグラフを描くためには，
sts graph
sts graph , na
２行目のコマンドでは累積ハザード関数を書くことができます．
28
8.7
ノンパラメトリック推定に挑戦しよう
パラメトリック推定は関数形に強い制約を置くので，情報処理技術の進展にともない，ノ
ンパラメトリック推定もしばしば利用されるようになってきています．とくに多くのデータ
が手許にあるときには，ノンパラメトリック推定は視覚的にも興味深い結果をもたらすこ
とがあります．ノンパラメトリック推定は，(1) データセットの特徴を示し，パラメトリッ
ク推定への手がかりを与える，(2) 推定されたパラメトリックモデルの診断を行う，(3) 弱
い制約のもとでのなんらかの統計的推測を行う，(4) セミパラメトリックモデルの推定を行
う，等の目的で使用されます．一般に，大量のマイクロデータを必要とします．
Density Estimation
ここでは，サンプルが独立に同一の未知の確率分布から発生している（i.i.d., identically
and independently distributed）と仮定しましょう．まず，一変数 X が密度関数 f (X) から
発生しているとします．f (x) が存在してなめらか（smooth）であるとすれば，中間値の定
理より，十分小さな実数 h > 0 に対して以下の近似が成立するとしてよいことになります．
2h · f (x) =
x+h
x−h
f (u) du = P r(X ∈ [x − h, x + h])
これを十分多くの x について行うことにより，密度関数 f (X) を推定することができます34 ．
つまり，密度関数の推定量 fˆh (X) は，
1
P r(X ∈ [x − h, x + h])
fˆh (x) =
2h
で与えられます．ところで，P r(X ∈ [x− h, x+ h]) すなわち確率変数 X が区間 [x− h, x+ h]
の値を取る確率はどのように推定すればよいのでしょうか？最も簡単な方法としては，幅
2h のこの区間 [x − h, x + h] に入っているサンプルの数を数え，全サンプル数 N で除し
た相対頻度を用いることが考えられます．定義域に入る条件が真なら１，偽なら０を取る
Indicator function I を用いれば以下のように表現できることになります．
fˆh (x) =
=
1 1
2h N
1
N
N
i=1
N
I(|Xi − x| ≤ h)
i=1
11
1
I(|Xi − x|/h ≤ 1) =
h2
N
この推定量のことを，kernel を K(u) =
1
2 I(|u|
N
i=1
1
1
K((x − Xi )/h) =
h
N
N
Kh (x − Xi )
i=1
≤ 1)， bindwidth を h とする f (x)
のkernel density estimator と呼びます．
このような kernel をstep function kernel と呼びますが，ここでは区間 [x − h, x + h]
に入っている観測値には全て同じウェイトが与えられているし，そのために推定量が不連続
になりやすい傾向があります．これらの問題を解決するためには，点 x より離れた観測値
に与える重みを連続的に小さくした kernel を用いる方法があります．このようなwindow
34 大げさに言えばこういうことですが，棒グラフの高級なもの，といったものとも言えましょうか．
29
function K を kernel とする方法として，２次・３次関数を用いる方法のほか，Gaussian
kernel と呼ばれるものがあります．具体的には以下の関数を用いることになります．
K(u) = (2π)−1/2 exp(−u2 /2)
一般に，kernel は連続（piecewise continuous）で，有界，０について対称，積分値が１に
なる関数が用いられ，kernel を表す関数 K(u) は以下を満たします．
K(u) = K(−u),
K(u) du = 1
さて，大量のデータがあるときには，ある変数の密度関数をノンパラメトリックに推定
することがあるかもしれません．その場合，
kdensity 変数名, kernel オプション
とします．Kernel としては，Gaussian kernel や triangle kernel などを用いることがで
きます．Kernel の種類によっては，さらにオプション指定をする必要があることがありま
す．デフォルトでは，コマンドを実行するとグラフが描かれますが，それぞれの点と推定
された密度の数値情報のデータgenerate( , ) オプションで作成することもできます．
Regression Estimation
ここまで述べてきたのは 1 変数の密度関数をノンパラメトリックに推定する方法でした
が，2 変数の関係についてノンパラメトリックに検討する方法について考えてみましょう．
つまり，条件付き期待値の関数 m(x) = E(Y |X = x) を推定することを考えてみましょう．
先の例と同様，２変数の確率変数の観測値ベクトル (X, Y ) は独立に同一の未知の確率分布
f (x, y) から発生していると仮定します．条件付き期待値の定義より，
yf (x, y) dy
f (x, y) dy
m(x) = E(Y |X = x) =
２変数の確率分布についての kernel density estimator も１変数と同様に考えることがで
き，次のように与えられます．
1
fˆh (x, y) =
N
N
Kh (x − Xi )Kh (y − Yi )
i=1
両辺を y について積分して，kernel の性質を考慮すると，
fˆh (x, y) dy
=
y fˆh (x, y) dy
1
N
=
1
N
=
1
N
Kh (x − Xi )Kh (y − Yi ) dy
i=1
N
Kh (x − Xi )
Kh (y − Yi ) dy =
i=1
1
N
=
N
1
N
N
Kh (x − Xi )
i=1
N
Yi Kh (x − Xi )Kh (y − Yi ) dy
i=1
N
Kh (x − Xi )Yi
i=1
30
Kh (y − Yi ) dy =
1
N
N
Kh (x − Xi )Yi
i=1
代入すれば，この推定量から条件付き期待値の推定量を得ることができます．
m
ˆ h (x) =
y fˆh (x, y) dy
=
fˆh (x, y) dy
1
N
1
N
N
i=1 Kh (x − Xi )Yi
N
i=1 Kh (x − Xi )
この推定量はNadaraya-Watson kernel estimate と呼ばれ，一致性を持ち，漸近的に
正規分布に従うことが知られています．Nadaraya-Watson kernel estimate は，以下のコマ
ンドで推定することができます．
kernreg Y 軸変数 X 軸変数, bwidth( 数値) kercode( 数値) npoint( 数値)
bwidth オプションは smoothing parameter と呼ばれるもので，kernel のパラメタを決
定します，kercode は kernel の種類を指定し，1∼7 の値をとります．npoint は，推定を
行うポイントの数を指定します．
検定してみよう
9
計量経済学で重要なのは推定ですが，しばしば統計的有意性を問う仮説検定を行う必要
があることがあります．ある変数の係数についての単純な検定，たとえば係数がゼロと等
しいかどうかといった検定は，出力結果にある係数推定値と標準偏差を用いれば t 検定す
ることができます．しかし，その他の形を採る場合には係数推定値の分散共分散行列が必
要となるためもう１ステップ必要になります．
9.1
推定のあとに
推定を行って，ちょっとややこしい仮説検定を行ったり，あるいは，係数推定値を使って
新たな変数を定義したりしたいこともあるでしょう．Stata は推定を行うと，いくつかの関
連する推定値などを一時変数として保存しています．たとえば，係数推定値は b[変数名] で
すし，とくに定数項の推定値は b[ cons] といったぐあいです．これらの表現とgenerate,
replace コマンドを併用すれば，新たな変数を定義することもできます．とくに当てはめ
値（predicted value, fitted value）の系列を作る場合には，
predict 新変数名, xb
とすれば，直前の推定から得られる当てはめ値を計算できます．predict はほとんどの
推定コマンドに対応する post estimation command ですが，プロビットやロジットのよう
な非線形な推定を行う場合には，何が生成されるのか注意してください．一般には，説明
変数の線形結合としての潜在変数あるいは被説明変数の当てはめ値が計算されます．また，
オプションstdp をxb のかわりに指定すると，標準偏差系列を得ることができます．当て
はめ値の計算は，すべての観測値について行われますが，推定に用いられた観測値につい
てのみ当てはめ値を計算したいときには，直前の推定に用いられたことを示すe(sample)
変数を用いて，
31
predict 新変数名 if e(sample), xb
とすればよいことになります．また，新変数名が既存の変数名と同じになると，エラー
が返ってきますので注意しましょう．
線形制約のばあい
9.2
線形制約の仮説検定では F 検定か χ2 分布を使った検定を行います．Stata では，係数
は b[変数名] という一時変数に保存されていますから，これを使って線形制約を表現して
やれば仮説検定を行うことができます．たとえば，次のような推定式と線形の帰無仮説を
考えてみましょう．
hoursi
β1
= β0 + β1 lnwi + β2 asseti + εi
= β2
Stata は線形の仮説検定のためにtest コマンドを用意しています．test のあとに線形制約
式を書けば検定をしてくれます．OLS を使った場合には F 検定をしますが，そのコマンド
ラインは以下のようになります．
reg hours lnw asset
test _b[lnw] = _b[asset]
係数の和が１に等しい，といった規模の経済性の検定などで使われそうな線形仮説の場
合は，
test _b[lnw] + _b[asset] = 1
となります．この場合，推定された係数の和自体は出力されないことに注意しましょう．
出力させたいときにはdisplay コマンドを使えばよいので，まとめると次のようなコマン
ドラインになるはずです．
reg hours lnw asset
display _b[lnw] + _b[asset]
test _b[lnw] + _b[asset] = 1
仮説が幾つかの式から成る場合にはaccumulate オプションを使うことになります．直
前の例で，さらに β1 = 0.5 という仮説を組み合わせる場合には以下のようになります．
test _b[lnw] + _b[asset] = 1
test _b[lnw] = 0.5, accumulate
もちろん，2 行目でaccumulate オプションを使わなければ，別々の検定をしてくれるこ
とになります．
32
9.3
非線形制約のばあい
test コマンドは便利ですが，非線形制約の検定はしてくれません．Non-Linear な仮説
検定をしたい場合には，testnl コマンドを使います．使い方はtest コマンドとほとんど
同じです（test コマンドではいくつかの表現を省略できますが，testnl コマンドではで
きません，など）．たとえば，先ほどの推定式で以下のような非線形制約を検定したいと
します．
1
= 0.3
1 + exp(β0 /β1 )
左辺の値をいったん出力させて検定するには次のようにすればいいことになります．
display 1/[1+exp(_b[_cons]/_b[lnw])]
testnl 1/[1+exp(_b[_cons]/_b[lnw])] =0.3
この場合，Wald 検定が行われます．2 つ以上の仮説をまとめて検定するときには，式を(
) でくくって並べればいいだけです．
10
ちょっと進んで
いつもいつも単純な推定式になるとは限りません．制約つきの式を推定することもある
でしょうし，非線形になっているかもしれません．最小 2 乗法だけを使うとも限りません．
ここではそのような場面に対応する最も単純なケースを扱います．
非線形最小 2 乗推定をしよう
10.1
まずは非線形最小 2 乗推定です．非線形最小 2 乗法を使う場合には，まず式を定義して
から，被説明変数を指定して推定するという順序を踏みます．推定式の右辺に当たる部分
を定義するにはprogram define ではじまり，end で終わる形式を使います．一般的には
program define nl 式名
version 6.0
if "‘1’" == "?"
{
global S_1 " パラメタ名 "
パラメタの初期化
exit
}
replace ‘1’ = 数式
end
で定義します．式の名前は必ずnl で始めます．推定は被説明変数を指定するだけです．
nl
式名
被説明変数名
33
たとえば，推定式が
hoursi = − ln
1−α
wi
α
+ ui
であるとしましょう．このとき推定したいパラメタは と α で，説明変数は定数項と wi ，
被説明変数は hoursi です．さきほどの形式にしたがって式を定義して推定するには以下の
ようにします．ここでは式の名前をnlsbst としましょう．
program define nlsbst
version 6.0
if "‘1’" == "?" {
global S_1 "alpha epsilon"
global alpha = 0.001
global epsilon = 0.4
exit
}
replace ‘1’ = -$epsilon * ln((1-$alpha)*w/$alpha)
end
nl sbst hours
非線形最小 2 乗法は収束計算をともなうので，初期値の設定が必要です．また，replace
のあとの表現では，パラメタの前には$をつける必要があり，この時点で説明変数を指定し
ておく必要があります．なので，いったん式を定義したあとで説明変数を変形したりする
ときには注意が必要です．非線形最小 2 乗法の場合も b[w] などの表現を使って検定を行
うことができます．
10.2
最尤法に挑戦してみよう
Stata はマイクロデータの処理に優れているのですが，マイクロデータを使うとさまざ
まな推定を行うことが可能になります．サンプルが多い場合には自分で尤度関数を書いて
最尤法を使って推定することもあるかもしれません．もちろん，Stata でも自分で最尤法を
使って推定することができます（収束計算の方法がすくないとの指摘もありますが…）．こ
こでは，全体の対数尤度がそれぞれの観測値の対数尤度の和で表現できるような，もっと
も簡単な場合の最尤法の，もっとも初歩的な使い方を紹介します．それゆえ，パネルデー
タの random effect のような，ちょっとややこしいケースは除きます．
Stata の最尤法プログラムとしては，lf, d0, d1, d2 の 4 通りが用意されていますが，
使いやすさと計算の正確さ・速さから，lf の使用が推奨されています．もっとも，lf は全
体の対数尤度が観測値の対数尤度の和で表現できる場合にしか対応していません．いくつ
かの観測値間に相関があるようなケースでは，それらをまとめてひとつの観測値としてし
まい（つまり横に並べてしまい），対数尤度の定義の仕方を工夫したほうが楽かもしれま
せん．
34
最尤法の使い方はさきほどの非線形最小 2 乗法とよく似ています．ただし，ここでは説
明変数の指定はあとで，関数の指定・被説明変数と説明変数の指定・推定，という順序を
採ります．ここで扱うような，観測値間の独立性が仮定されているケースではlf オプショ
ンを使います．まず一般形を書いておくと次のようになります．
program define 式名
args lnf theta1 theta2 ...
tempvar 一時変数名
対数尤度関数の定義
quietly replace ‘lnf’ = ...
end
ml model lf 式名 (被説明変数 = 説明変数名)...
ml check
ml maximize
これらの途中で初期値の設定などのオプション指定もできますがここでは省略します．
最尤法では式の名前に特に制限はありませんが，コマンド名になっているものは避けるべ
きです．lnf は対数尤度です．program の定義文の最後でquietly replace ‘lnf’ = の
形で締めくくることになります．theta1, theta2,... は説明変数（と定数項）の線形和を
表します．推定されるのは，theta たちの係数になります．標準偏差など，定数を推定す
るばあいにもこれを使います．推定は収束計算をともなうので，推定される値の域が実数
全体になるように工夫します．被説明変数とは，ここでは係数が付かない変数を指します．
対数尤度関数を定義するときには$ML y1, $ML y2,... といった形で表現されます．一時変
数とは，対数尤度関数を定義するためだけに一時的に定義する変数です．対数尤度関数が
複雑な形をとるとき，一時変数を使って見通しをよくすることができます．対数尤度関数
を定義したあと，被説明変数と説明変数を指定します． = でつながれていますが，これに
特に意味はない（最も簡単なケースでは意味がありそうに見える）だけなので気にしない
でください．
最尤法の使い方はやってみないとわかりにくいところがあるのですが，ここでひとつ例
をだします（Stata での最尤法の入門書，Gould and Sribney 1999, 27 ページ35 ）．Weibull
分布を仮定したサバイバル分析を行うケースです．この場合，対数尤度関数は以下のよう
にかけます．
ln fi = (ti exp(xi β))exp(s) + di (s − xi β + (exp(s) − 1)(ln ti − xi β))
ここで推定したいのは，説明変数の線形結合の係数 β ，Weibull 分布のパラメタ s です．ま
た，ti , di は係数が付かないのでここでは被説明変数と呼べます．いま，説明変数の線形結
合 xi β をtheta1 ，パラメタ s をtheta2 とし，ti を１つめの被説明変数，di を２つめの被
説明変数と考えてやれば，Stata で対数尤度関数を定義することができます．推定される係
数はいずれも（s も），すべての実数値を取りうることに注意してください．exp(s) に一時
35 Gould, William and William Sribney. 1999. Maximum likelihood estimation with STATA. Stata
Press, Texas.
35
変数p，右辺第 1 項に一時変数M，ln ti − xi β に一時変数R を割り当てて，式を読みやすく
してやると，Stata による尤度関数の定義はつぎのようになります．
program define myweib
version 6.0
args lnf theta1 theta2
tempvar p M R
quietly generate double ‘p’ = exp(‘theta2’)
quietly generate double ‘M’ = ($ML_y1*exp(-‘theta1’))^‘p’
quietly generate double ‘R’ = ln($ML_y1)-‘theta1’
quietly replace ‘lnf’ = -‘M’+$ML_y2*(‘theta2’-‘theta1’+(‘p’-1)*‘R’)
end
対数尤度関数を定義するときにはgenerate, replace はよく使いますが，いずれもオプ
ションquietly とdouble を使っていることに注意してください．これらがないと出力がひ
どくなったり，収束計算の精度が非常に落ちたりします．また，パラメタや一時変数はシ
ングルクォーテーションでくくっていることにも注意しましょう．ここではまだ説明変数
はなんでもいいことになっています．
次に説明変数と被説明変数を指定します．変数名は ti がtime，di がdied，xi として定
数項とage であるとしましょう．すると，
ml model lf myweib (time died = age) ()
となります．2 個目のカッコはtheta2 に対応していますので省略はできません．コロン
を使ってそれぞれに式の名前をつけることもでき，その場合には，
ml model lf myweib (timeeq: time died = age) (sigma: )
などとなります．被説明変数・説明変数の順序は，対数尤度関数の定義のときの$ML y1,
$ML y2,..., theta1, theta2,..., の順序に従います．被説明変数は = の左側，説明変数は
右側に来ますが，説明変数だけをカッコのなかにいれることもできます．被説明変数だけ
をカッコのなかにいれることはできません．被説明変数をいれるカッコは順序をまちがわ
なければどこでもいいことにはなっていますが，よくわからなければ，すべて最初のカッ
コの等号の左側に書いておけばよいでしょう．もし，説明変数の線形結合が定数項を含ま
ない場合には次のようになります．
ml model lf myweib (timeeq: time died = age, nocons) (sigma: )
ここまで書けば対数尤度の定義は終わりです．あとは最大化する命令を書けばいいこと
になります．
ml check
ml maximize
36
ここは決め打ちなので，深く考えないほうがいいです．最尤法では，最大化の収束計算
を始める前に，対数尤度がきちんと定義されているかどうかのチェックをしてくれます．そ
のチェックを通れば，とりあえずプログラムとしては間違いがないということになるよう
です．
11
おわりに
これまで，いくつかの分析手法についてのもっともベーシックな形でのコマンドライン
を紹介してきました．しかしもちろん，いろいろなオプションを指定できますし，これ以
上のさまざまな分析手法を Stata は用意しています．あとは自分で Reference Manual と
（計量経済学や統計学の教科書と）格闘するしかないのかもしれません．
37
© Copyright 2024 ExpyDoc
