●【課題 4】の解説 1.理 論(地域システム論 A 講義資料 No.3 より) ハフモデルの式に多重積形の誤差項をくわえ、その対数をとって、線形回帰モデル に適用する。 OBS ハフの原モデルを取り上げよう。行先選択確率の観測値を Pij PRD る予測値を Pij 、ハフモデルによ ε 、誤差項を e ij とおく。m 居住地、n 商業地の場合を取り上げよう。 Sj PijOBS = PijPRD e ε ij = Tijλ εij e n S j (1.1) ∑Tλ j =1 ij 式(1.1)の両辺の対数をとると、 Sj log PijOBS Tijλ εij e = log n Sj ∑Tλ j =1 ij log PijOBS log n = − log(∑ Sj λ j =1 Tij PijOBS ) + log S j − λ log Tij + ε ij = −α i − λ log Tij + ε ij , i = 1,..., m, j = 1,..., n Sj n ここで、 i のみに依存する定数 log( (1.2) (*) S ∑ T λj ) を、 α i とおいている。 j =1 ij 式(1.2)の最後の式(*)は、以下のようにして(単)線形回帰モデルに帰着できる。 (a) 自由度の考慮 OBS 確率 Pij は、 j に関する総和が 1 であるので、自由に動ける数(自由度)は n − 1 である。そこで、最後のデータを一つ取り除く。 (b) 居住地 i を固定し、居住地ごとの回帰を実行する方法 居住地を i = 1 に固定する。変数を次のように対応させる。 Y j = log P1OBS j Sj , X j = log T1 j パラメータ a1 , b1 を次のように対応させる。 a1 = −α1 , b1 = −λ 以上のもとで、居住地 1 に対する式(1.2)は次の単回帰モデルに帰着する。 Y j = a1 + b1 X j + ε1 j , j = 1,..., n − 1 同様に、他の居住地 i に対しても、個別の単回帰モデルを以下のように適用で きる。 Y j = ai + bi X j + ε ij , j = 1,..., n − 1, for i (1.3) 居住地別の単回帰モデル適用の問題点 • 距離抵抗の係数 λ は、理論上はすべての居住地で同一であるのに、居住 地ごとに異なってしまう。 次のモデルはこの問題点を解決したものである。 (c) ハフモデルの推定 (居住地を一括した単回帰モデルの適用) 居住地別の単回帰モデルの式(1.3)において、m コの居住地のサンプルをプール することを考える。これは次のように、サンプルの番号付けを変えればよい。 (i, j ) → k = (n − 1)(i − 1) + j , i = 1,..., m, j = 1,..., n − 1 (1.4) すなわち、変数の対応を次のようにつける。 Y( n −1)(i −1)+ j = log PijOBS Sj , X ( n −1)(i −1)+ j = log Tij , i = 1,..., m, j = 1,..., n − 1 (1.5) また、パラメータ bi は、居住地によらず一定で、 bi = b = −λ とおく。パラメー タ ai は、居住地ごとの定数項に対する係数であったから、1,0 の値をとるダミ i ー変数と呼ばれる、次の変数 E , i = 1,..., m を導入する。 1 if [(k − 1) (n − 1)] + 1 = i Eki = (1.6) otherwise 0 以上のもとで、式(1.2)の最後の式である次式(1.7)は log PijOBS Sj = −α i − λ log Tij + ε ij , i = 1,..., m, j = 1,..., n − 1 (1.7) 式(1.4)、(1.5)、(1.6)を用いて、次の線形重回帰モデル(1.8)に帰着する。 Yk = a1Ek1 + " + am Ekm + bX k + ε k , k = 1,..., m(n − 1) (1.8) 2.プログラム 課題 4 では、(1.8)式から最小 2 乗法でパラメータを推定するプログラムを示している。 PROC IML; RESET NOPRINT; RESET LOG; /** INPUT DATA **/ S={5000 15000 20000}; T={10 20 15,5 5 5,30 10 15}; POBS={0.3 0.22 0.48,0.13 0.36 0.51,0.03 0.61 0.36}; *仮設例の数値を入力している。① *地区 1 ダミー; E1={1 1 1, 0 0 0, 0 0 0}; *地区 2 ダミー; E2={0 0 0, 1 1 1, 0 0 0}; *地区 3 ダミー; E3={0 0 0, 0 0 0, 1 1 1}; *E1~E3 までは地区ダミー変数を定義している。② YY=POBS[,1:2]; SS=REPEAT(S,3,1); SS=SS[,1:2]; TT=T; TT=TT[,1:2]; E1=E1[,1:2]; E2=E2[,1:2]; E3=E3[,1:2]; *ここまでは、選択確率の和が 1 であることを考慮して、最後のデータを 1 つ取り除いてい る。③ YY=SHAPE(YY,1); YY=YY`; SS=SHAPE(SS,1); SS=SS`; TT=SHAPE(TT,1); TT=TT`; *ここでは、データを縦に並べている。④ Y1=LOG(YY/SS); *上では、 PijOBS Sj の log をとっている。⑤ X1=LOG(TT); *上では、 T の log をとっている。⑥ E1=SHAPE(E1,1); E1=E1`; E2=SHAPE(E2,1); E2=E2`; E3=SHAPE(E3,1); E3=E3`; *ダミー変数を縦に並べている。⑦ XX=E1||E2||E3||X1; PijOBS Sj *最小二乗法の公式に合うようにダミー変数と log を結合して、X を作っている。⑧ Y=Y1; X=XX; BETA =INV(X`*X)*X`*Y; −1 *最小 2 乗法の公式 βˆ = ( X ′X ) X ′Y より、βを推定する。⑨ PRINT BETA; 3.プログラムの流れ プログラムの流れを行列を使ってみていこう。 下の番号は、上のプログラムの解説の部分にある番号と対応している。 ① 仮設例の数値を入力している。 A B 1 0.3 0.22 P OBS = 2 0.13 0.36 3 0.03 0.61 A B S = [5000 15000 A B C 0.48 0.51 0.36 C 20000] C 1 10 20 15 T = 2 5 5 5 3 30 10 15 ② ダミー変数を定義している。 A B C 1 1 1 1 E1 = 2 0 0 0 3 0 0 0 ③ A B C 1 0 0 0 E2 = 2 1 1 1 3 0 0 0 A B C 1 0 0 0 E3 = 2 0 0 0 3 1 1 1 選択確率の和が 1 であることを考慮して、最後のデータを 1 つ取り除く。 A B A B A B P OBS 1 0.3 0.22 = 2 0.13 0.36 3 0.03 0.61 5000 15000 S = 5000 15000 5000 15000 1 10 20 T = 2 5 5 3 30 10 A B 1 1 1 E1 = 2 0 0 3 0 0 ④ A B 1 0 0 E2 = 2 1 1 3 0 0 データを縦に並べておく。 0.3 0.22 A B 1 0.3 0.22 0.13 P OBS = 2 0.13 0.36 → 0.36 3 0.03 0.61 0.03 0.61 5000 15000 A B 5000 15000 5000 S = 5000 15000 → 15000 5000 15000 5000 15000 10 20 A B 1 10 20 5 T = 2 5 5 → 5 3 30 10 30 10 ⑤ PijOBS Sj の log をとって、Y と定義する。 Y = log P OBS Sj A B 1 0 0 E3 = 2 0 0 3 1 1 −9.7212 −11.1299 −10.5574 = −10.6375 −12.0238 −10.1101 ⑥ T の log をとる 2.3026 2.9957 1.6094 log T = 1.6094 3.4012 2.3026 ⑦ ダミー変数を縦に並べている。 1 1 1 1 1 0 E1 = 2 0 0 → 0 3 0 0 0 0 A B ⑧ 0 0 1 0 0 1 E2 = 2 1 1 → 1 3 0 0 0 0 A B 0 0 1 0 0 0 E3 = 2 0 0 → 0 3 1 1 1 1 A B 重回帰モデルの公式にあうように X と Y を作る。 1 0 0 2.3026 −9.7212 1 0 0 2.9957 −11.1299 OBS 10.5574 0 1 0 1.6094 − P Y = log = X = E1 E2 E3 log T = Sj −10.6375 0 1 0 1.6094 −12.0238 0 0 1 3.4012 −10.1101 0 0 1 2.3026 −1 ⑨ 最小 2 乗法の公式 βˆ = ( X ′X ) X ′Y より、βを推定する。出力される結果は、X の並べ て順に出てくるので、βは最後の-1.82 がβの値である。 −5.59 −7.66 = ( X `X ) −1 X `Y −5.86 −1.82
© Copyright 2025 ExpyDoc