 |
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: st: spatial analysis - DISTANCE
Austin,
This is a great start. In all of my tinkering with the code
yesterday, I don't know why I didn't think of that. Thanks, as always.
--------------------------------------------
David Diego Torres
Quoting Austin Nichols <[email protected]>:
David Torres<[email protected]> :
The key question to answer yourself before you adopt any of the many
strategies elucidated on the list is what you are going to do with a
separate list of identifiers for each observation on your original
dataset. Having this answer in hand, think about how the data should
be organized to do that efficiently. Then get the lists together.
Here is one way of doing the last in your original example (note I
changed only one thing, an equality test to an inequality), but I
cannot see the utility of this exercise, so it is unlikely to be what
you actually want:
g long rad10=.
g double dist1=.
g id1=""
merge using univpoints00_0
forv i=1/`np' {
qui vincenty `=latitude[`i']' `=longitude[`i']' sch_lat sch_lon, hav(dis)
g within10=dis<10
su within10, meanonly
qui replace rad10=r(sum) in `i'
su dis, meanonly
qui replace dist1=r(min) in `i'
qui levelsof unitid if dis<10, loc(vs)
foreach v of loc vs {
qui replace id1="`v' "+id1 in `i'
}
drop dis within10
}
On Sun, Aug 30, 2009 at 10:15 AM, David Torres<[email protected]> wrote:
Austin, I've looked over the code you offered and it seems that what I'm
getting is just the xy coordinates rather than what I need. There's got to
be a way, somehow, to have stata return a list of all the unitids that are
within the specified radius. When I run the code I have, all I get is a
return of the unitid on the nearest school to my 65000 tracts. However,
within the rad10 variable, several tracts have multiple schools within my
radius of 10 miles. I need to get all of these unitids regardless of
whether I end up with wide or long data; I can always reshape later.
Is my question clearer?
Well, any help you or anyone else can give is greatly appreciated. In the
meantime, I'll keep working at it.
Thanks
David D. Torres
Quoting Austin Nichols <[email protected]>:
David Torres<[email protected]> :
This is in reference to e.g.
http://stata.com/statalist/archive/2009-07/msg00083.html
http://stata.com/statalist/archive/2007-01/msg00098.html
I believe.
Note that standard Stata commands like -summarize- and -replace-
inside the loop do all the work of saving characteristics of the
nearest, or k nearest, observations (and whatever else, such as mean
of Z over those with distance less than 10 units). No program like
-vincenty- (on SSC) does that; a -replace- command does it.
-vincenty- does not return "the min, max, mean, or whatever" and it
appears that -globdist- (which I found with -findit-) does only part
of what -vincenty- does (and the -globdist- help file contains the
caveat that "Testing of this program has been limited"). In other
words, add a -replace- to put the result you want in a variable.
Maybe you want something like:
g long rad10=.
g double dist1=.
g id1=""
g xy1=""
merge using univpoints00_0
forv i=1/`np' {
qui vincenty `=latitude[`i']' `=longitude[`i']' sch_lat sch_lon, hav(dis)
g within10=dis<10
su within10, meanonly
qui replace rad10=r(sum) in `i'
su dis, meanonly
qui replace dist1=r(min) in `i'
qui levelsof unitid if dis==dist1[`i'], loc(vs)
foreach v of loc vs {
qui replace id1="`v' "+id1 in `i'
su sch_lat if unitid==`v', meanonly
loc y1=r(mean)
su sch_lon if unitid==`v', meanonly
loc x1=r(mean)
qui replace xy1="(`x1',`y1') "+id1 in `i'
}
drop dis within10
}
On Sat, Aug 29, 2009 at 4:26 PM, David Torres<[email protected]> wrote:
I am trying to calculate multiple distances from one set of xy
coordinates
to another
using vincenty. In addition to getting those distances, however, I need
to
somehow list
the ids for the xy coordinates used to calculate the distance.
Here is the piece of code I'm working with:
code excerpt:
g long rad10=.
g double dist1=.
g id1=""
merge using univpoints00_0
forv i=1/`np' {
qui vincenty `=latitude[`i']' `=longitude[`i']' sch_lat sch_lon, hav(dis)
g within10=dis<10
su within10, meanonly
qui replace rad10=r(sum) in `i'
su dis, meanonly
qui replace dist1=r(min) in `i'
qui levelsof unitid if dis==dist1[`i'], loc(vs)
foreach v of loc vs {
qui replace id1="`v' "+id1 in `i'
}
drop dis within10
}
end excerpt
As you can see, I'm getting a count of all sch_lat and sch_lon
coordinates
within a 10
mile radius of the other set of latitude and longitude points. Then I
return only the
minimum distance. Will vincenty allow for returning more than just the
min,
max, mean,
or whatever? Or is globdist a better choice for returning multiple
distances and the ids
attached to the point the distance to which is being calculated?
Thanks
David Torres
UMich
*
* 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/
*
* 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/