Statalist


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: st: data points within radius?


From   Austin Nichols <[email protected]>
To   [email protected]
Subject   Re: st: data points within radius?
Date   Tue, 7 Jul 2009 09:24:45 -0400

Here is another example along the lines of
http://www.stata.com/statalist/archive/2007-01/msg00098.html for
calculating properties of each observation depending on its distance
from all points on another dataset:

clear all
ssc inst vincenty, replace
tempfile p h
range pid 1 10 10
set seed 123
g p_lat=32+uniform()*10
g p_lon=-120+uniform()*40
save `p'
clear
input id cntrl VA h_lat h_lon
0785 45 1 42.103 -80.065
0787 45 1 42.103 -80.065001
0150 45 1 39.739 -75.605
0960 45 1 37.499 -77.470
0510 45 1 39.464 -77.962
0030 45 1 27.844 -82.796
0180 45 1 46.612 -112.048
0060 45 1 43.618 -116.194
0085 33 0 41.242 -110.991
0540 45 1 34.568 -112.456
0075 33 0 40.698 -111.990
0190 33 0 40.044 -111.716
1180 45 1 46.053 -118.356
0009 33 0 33.860 -118.149
1538 33 0 34.093 -118.344
end
isid id
save `h'
desc using `h', sh
local nh=r(N)
use `p'
local np=_N
g long n1k=.
g double mindist=.
g ids=""
merge using `h'
forv i=1/`np' {
qui vincenty `=p_lat[`i']' `=p_lon[`i']' h_lat h_lon, hav(d)
g within1000=d<1000
su within1000, meanonly
qui replace n1k=r(sum) in `i'
su d, meanonly
qui replace mindist=r(min) in `i'
qui levelsof id if d==mindist[`i'], loc(vs)
foreach v of loc vs {
 qui replace ids="`v' "+ids in `i'
 }
drop d within1000
}
drop if _m==2
drop id cntrl VA h_lat h_lon _m
compress
la var n1k "Number of VA hosp w/in 100mi"
la var mindist "Distance to nearest VA hosp"
la var ids "ID number(s) of nearest hospital(s)"
describe
list
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/



© Copyright 1996–2024 StataCorp LLC   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index