1

この質問は、Stata からのフォローアップの質問です: replace, if, forvalues。次のデータを検討してください。

set seed 123456
set obs 5000
g firmid = "firm" + string(_n)    /* Observation (firm) id */
g nw = floor(100*runiform())      /* Number of workers in a firm */
g double lat = 39+runiform()      /* Latitude in decimal degree of a firm */
g double lon = -76+runiform()     /* Longitude in decimal degree of a firm */

最初の 10 個の観測値は次のとおりです。

     +--------------------------------------+
     | firmid   nw         lat          lon |
     |--------------------------------------|
  1. |  firm1   81   39.915526   -75.505018 |
  2. |  firm2   35   39.548523   -75.201567 |
  3. |  firm3   10   39.657866    -75.17988 |
  4. |  firm4   83   39.957938   -75.898837 |
  5. |  firm5   56   39.575881   -75.169157 |
  6. |  firm6   73   39.886184   -75.857255 |
  7. |  firm7   27    39.33288   -75.724665 |
  8. |  firm8   75   39.165549    -75.96502 |
  9. |  firm9   64   39.688819   -75.232764 |
 10. | firm10   76   39.012228   -75.166272 |
     +--------------------------------------+

会社 1 と他のすべての会社の間の距離を計算する必要があります。したがって、vincentyコマンドは次のようになります。

. scalar theLat = 39.915526
. scalar theLon = -75.505018
. vincenty lat lon theLat theLon, hav(distance_km) inkm

vincenty コマンドは、各観測値と企業 1 の間の距離を持つdistance_km変数を作成します。ここでは、39.915526 と -75.505018 の 2 つの数値を手動でコピーして貼り付けます。

質問 1 : これらの数値を抽出する構文は何ですか?

これで、distances_km <= 2 の観測を維持できます。そして、

. egen near_nw_sum = sum(nw)

会社 1 から 2 km 以内にあるワーカーの合計を作成します (または、collapseコマンドで実行できます)。

質問 2 : すべての企業に対してこれを行う必要があり、最終的なデータは次のようになります。

     +-----------------------------------------------------------------+
     | firmid   nw         lat          lon            near_nw_sum     |
     |-----------------------------------------------------------------|
  1. |  firm1   81   39.915526   -75.505018  (# workers near firm1)    |
  2. |  firm2   35   39.548523   -75.201567  (# workers near firm2)    |
  3. |  firm3   10   39.657866    -75.17988  (# workers near firm3)    |
  4. |  firm4   83   39.957938   -75.898837  (# workers near firm4)    |
  5. |  firm5   56   39.575881   -75.169157  (# workers near firm5)    |
  6. |  firm6   73   39.886184   -75.857255  (# workers near firm6)    |
  7. |  firm7   27    39.33288   -75.724665  (# workers near firm7)    |
  8. |  firm8   75   39.165549    -75.96502  (# workers near firm8)    |
  9. |  firm9   64   39.688819   -75.232764  (# workers near firm9)    |
 10. | firm10   76   39.012228   -75.166272  (# workers near firm10)   |
     +-----------------------------------------------------------------+

near_nw_sum変数を作成することが私の最終目標です。私の弱いデータ管理スキルのためにここであなたの助けが必要です.

4

1 に答える 1

2

以下は、基本的にここにあるのと同じ戦略であり、「最終目標」に基づいています。繰り返しますが、元のデータセットのサイズによっては便利です。joinby観測を作成するため、Stata の制限を超える可能性があります。しかし、私はそれがあなたが望むことをすると信じています。

clear all
set more off

set seed 123456
set obs 10
g firmid = _n   /* Observation (firm) id */
g nw = floor(100*runiform())      /* Number of workers in a firm */
g double lat = 39+runiform()      /* Latitude in decimal degree of a firm */
g double lon = -76+runiform()     /* Longitude in decimal degree of a firm */
gen dum = 1
list

* joinby procedure
tempfile main
save "`main'"

rename (firmid lat lon nw) =0
joinby dum using "`main'"
drop dum

* Pretty print
sort firmid0 firmid
order firmid0 firmid
list, sepby(firmid0)

* Uncomment if you do not want to include workers in the "base" firm.
*drop if firmid0 == firmid

* Compute distance
vincenty lat0 lon0 lat lon, hav(distance_km) inkm
keep if distance_km <= 40 // an arbitrary distance
list, sepby(firmid0)

* Compute workers of nearby-firms
collapse (sum) nw_sum=nw (mean) nw0 lat0 lon0, by(firmid0)
list

それが行うことは、企業のペアごとの組み合わせを形成して、距離を計算し、近くの企業の労働者を合計することです。ここでは、質問 1 で尋ねたようにスカラーを抽出する必要はありません。また、firmid文字列に変換する変数を複雑にする必要もありません。

以下は、観測数に関するStataの制限の問題を克服します。

clear all
set more off

* Create empty database
gen x = .
tempfile results
save "`results'", replace

* Create input for exercise
set seed 123456
set obs 500
g firmid = _n   /* Observation (firm) id */
g nw = floor(100*runiform())      /* Number of workers in a firm */
g double lat = 39+runiform()      /* Latitude in decimal degree of a firm */
g double lon = -76+runiform()     /* Longitude in decimal degree of a firm */
gen dum = 1
*list

* Save number of firms
local size = _N
display "`size'"

* joinby procedure
tempfile main
save "`main'"

timer clear 1
timer clear 2
timer clear 3
timer clear 4

quietly {
    timer on 1
    forvalues i=1/`size'{
        timer on 2
        use "`main'" in `i', clear // assumed sorted on firmid
        rename (firmid lat lon nw) =0

        joinby dum using "`main'", unmatched(using)
        drop _merge dum
        order firmid0 firmid
        timer off 2

        timer on 3
        vincenty lat0 lon0 lat lon, hav(dist) inkm
        timer off 3
        keep if dist <= 40 // an arbitrary distance

        timer on 4
        collapse (sum) nw_sum=nw (mean) nw0 lat0 lon0, by(firmid0)

        append using "`results'"
        save "`results'", replace
        timer off 4
    }
    timer off 1
}

use "`results'", clear
sort firmid0
drop x
list

timer list

しかし非効率的ですが、 を使用したいくつかのテストtimerでは、ほとんどの計算時間はvincentyエスケープできないコマンドに費やされることが示されています。以下は、Intel Core i5 プロセッサと従来のハード ドライブ (SSD ではない) を使用した場合の 10,000 回の観測時間 (秒単位) です。タイマー 1 は合計で、2、3、4 はコンポーネント (概算) です。タイマー 3 は以下に対応しvincentyます。

. timer list
   1:   1953.99 /        1 =    1953.9940
   2:    169.19 /    10000 =       0.0169
   3:   1669.95 /    10000 =       0.1670
   4:     94.47 /    10000 =       0.0094

もちろん、両方のコードで距離の重複計算が行われることに注意してください (たとえば、会社 1-会社 2 と会社 2-会社 1 の間の両方の距離が計算されます)。これはおそらく回避できます。現状では、110,000 回の観測には長い時間がかかります。良い面として、この 2 番目のセットアップでは、最初のセットアップでの同量の観測と比較して、必要な RAM が非常に少ないことに気付きました。実際、私の 4GB マシンは後者でフリーズします。

また、私はあなたと同じシードを使用していますが、異なる数の観測値 (5000 ではない) を作成しているため、データが異なり、変数作成プロセスに違いがあることに注意してください。

(ちなみに、値をスカラーとして保存したい場合は、添え字を使用できます: scalar latitude = lat[1])。

于 2013-11-05T08:02:35.937 に答える