实用小程序记录1-计算地球两点之间的距离

 

修改：20191212。

关于geodist 命令，文中的代码写法有问题。重新修改了一个可用的例子代码如下

 

 

*firstly,This process calls for a use written command named "geodist"
*type in codes here below and install this command
*ssc install geodist


*Simulating source data files
*variable meaning
*x-latitude 
*y-longitude

clear
set more off

//simulating latitude and longitude information.
set obs 10
gen x = uniform()*90
gen y = uniform()*180


***********************
* display simulated address
list 

// matrix used to save data 
matrix Mat_dist = I(10)
//compute distance between every two points, store the distance in macros
forval i = 1/10{
	forval j = 1/10{
		preserve
		gen Ax = x[`i']
		gen Ay = y[`i']
		gen Bx = x[`j']
		gen By = y[`j']
		keep in 1
		geodist Ax Ay Bx By, gen(dist)
		matrix Mat_dist[`i',`j'] = dist[1]
		restore 
	}
}
matrix list Mat_dist

 

 

 

 

 

分割线：

--------------------------------------------------------------------

 

 

 

 

问题描述：

对一定区域范围内医院面临的竞争程度进行度量。

 

方案：

1.以每家医院周边一定范围内的医院数量作为竞争程度的衡量指标。

2.实现思路：

 

• 得到每家医院的坐标位置（经度、维度）
• 计算每两家医院的距离
• 判断距离是否小于临界值
• 计数，符合判定条件的医院总数。

 

 

 

 

3. Stata 实现代码

*firstly,This process calls for a use written command named "geodist"
*type in codes here below and install this command
*ssc install geodist


*Simulating source data files
*variable meaning
*x-latitude 
*y-longitude

clear
set more off

//simulating 50 points with latitude and longitude information.
set obs 10
gen x = uniform()*100
gen y = uniform()*100

//compute distance between every two points, store the distance in macros
forval i = 1/10{
	forval j = 1/10{
		local Ax = x[`i']
		local Ay = y[`i']
		local Bx = x[`j']
		local By = y[`j']
		local dist_`i'_`j' = geodist(Ax,Ay,Bx,By)
	}
}

//generate a data file containing distance information
forval i = 1/10{
	gen v`i' =.
	forval j = 1/10{
		replace v`i' = `dist_`i'_`j'' if _n == `j'
	}
}

//judge if distance between two points is less than 30.
forval i = 1/10{
	replace v`i' = . if v`i'>=30|v`i'==0  
}
//cout numbers 
forval i = 1/10{
	egen num`i' = count(v`i')
}
keep num*
keep if _n==1

 

4.Stata代码的总结

首先通过模拟的方法生成了10个医院的经度、维度数据，分别存放在两个变量中。然后利用双重循环计算出了每家医院与其他医院的距离，并存储在local中，然后将这些距离的数据转存为一个矩阵形状的数据中。效果如下：

随后需要做如下处理：对于每家医院，距离小于临界值（这里选择30）的医院有多少个？对应到上述数据，在每一列中，变量取值小于30且不等于0的又几个？
这里用到了Stata的egen命令下的count函数:

 

        count(exp)                                                                             (allows by varlist:)  
            creates a constant (within varlist) containing the number of nonmissing observations of exp.  Also see  rownonmiss() and rowmiss().

可以利用count统计变量中非missing得record的数目，因此，可先将v*的取值大于30活等于0的替换为missing，然后计算non-missing的数量。

 

missing 替换的结果：

计算出的复合条件的number结果