根据《白话空间统计之九：方向分布(标准差椭圆)修正版》(有些地方没有理解清楚)，写了下面的程序。但是好像结果不对

Z=mvnrnd([0.5 1.5], [0.025 0.03 ; 0.03 0.16], 50);

X=Z(:,1);  Y=Z(:,2);

mean_X=nanmean(X);  mean_Y=nanmean(Y);   %椭圆圆心

%确定长短半轴

SDE_X=nanstd(X,1);

SDE_Y=nanstd(Y,1);

%确定方向

N=size(X,1)*size(X,2)-size(find(isnan(X)),1);  %非NaN元素个数

X2=nanvar(X(:),1)*N;

Y2=nanvar(Y(:),1)*N;

A1=nancov(X,Y,1);  % 是2*2矩阵，含4个元素，元素1是S(X(:))，及X的方差；元素4是Y的方差；元素2与3相等，是X与Y的协方差；

X1=A1(1)*N; % (X-mean_X)平方求和

Y1=A1(4)*N; % (Y-mean_Y) 平方求和

XY=A1(2)*N; %(X-mean_X)(Y-mean_Y) 求和

A=X1-Y1;

B=sqrt(A^2+4*XY^2);

C=2*XY;

theta=atan((A+B)/C);  % 椭圆方向，狐度，以X轴为准，正北方(12点方向)为0度，顺时针旋转    *180/pi=角度

%确定x、y轴的标准差

XStd=sqrt((cos(theta)^2*X1+sin(theta)^2*Y1-sin(2*theta)*XY)/N)*sqrt(2);  %X轴标准差，不知道是否有*sqrt(2)

YStd=sqrt((sin(theta)^2*X1+cos(theta)^2*Y1+sin(2*theta)*XY)/N)*sqrt(2);  %Y轴标准差，不知道是否有*sqrt(2)

% XStd=sqrt((cos(theta)^2*X1+sin(theta)^2*Y1-sin(2*theta)*XY)/N);  %X轴标准差，不知道是否有*sqrt(2)

% YStd=sqrt((sin(theta)^2*X1+cos(theta)^2*Y1+sin(2*theta)*XY)/N);  %Y轴标准差，不知道是否有*sqrt(2)

STD = 3;                     %# 3 standard deviations

conf = 2*normcdf(STD)-1;     %# covers around 99% of population

a=XStd*sqrt(conf);   b=YStd*sqrt(conf);

% a=0.504770;  b=13.867618;

A=a^2*sin(theta)^2+b^2*cos(theta)^2;

B=-(a^2-b^2)*sin(2*theta);

C=a^2*cos(theta)^2+b^2*sin(theta)^2;

f=-a^2*b^2;

X0=mean_X;  Y0=mean_Y;

plot(X,Y,'r.')

hold on

ezplot(subs('A*(x-X0)^2+B*(x-X0)*(y-Y0)+C*(y-Y0)^2+f=0'))

其中，我们前面说的好几种算法，如中心要素、中位数中心和平均中心，都是关于点方位的分析，那么今天我们要讲的这个算法，就是同时对点的方向和分布进行分析的一种经典算法——标准差椭圆。

这算法最早是由美国南加州大学(Universityof Southern California)社会学教授韦尔蒂.利菲弗(D. Welty Lefever)在1926年提出，所以有的书里面，也把这个算法称为Lefever's "Standard DeviationalEllipse"(利菲弗方向性分布)(又到每天的历史起源科普时间……)。

这个算法最大的特点，就如同他的名词一样，是用来度量一组数据的方向和分布的，生成的结果又正如他的别名一样，会输出一个椭圆，如下：

红色的点是伤寒发病的案例，蓝色的河流是长江太湖流域段，从计算的结果来看，发病的数据方向与长江的流向方向基本一致，而范围较大。

从上图，我们基本上就可以看出方向分布工具的主要作用了，它可以识别一组数据的方向以及分布的趋势，并且了解到这份数据是否具有一些特性，至于有哪些特性，我们后面再说。

我们先来看看这个标准差椭圆的生成算法。

其实算法很简单，要画出一个椭圆，虽然比画圆麻烦点，但是也麻烦不了多少，关键的参数如下：

1、确定圆心。

2、确定旋转角度。

3、确定XY轴的长度。

首先是确定圆心，方向分布工具的圆心，直接利用的是算数平均中心来计算椭圆的圆心(算术平均中心请查看我在2015年8月17日写的《空间统计之八：平均中心和中位数中心》一文)

然后就确定椭圆的形式了，公式如下：

其中，Xi和Yi是每个要素的空间位置坐标，X和Y是算数平均中心。

SDEx和SDEy就是计算出来的椭圆的方差，总所周知，椭圆的大小取决于方差大小，长半轴表示最大方差，短半轴表示最小方差，在空间统计上面，用X、Y的方差进行计算，得到长短半轴。

然后确定椭圆的方向，以X轴为准，正北方(12点方向)为0度，顺时针旋转，计算公式如下：

最后确定XY轴的标准差，公式如下：

标准差的作用是确定椭圆的方程，一般椭圆方程如下：

S是置信度的值，可以根据数据量来查询卡方概率表(Table:Chi-Square Probabilities)，这个大家有兴趣去百度一下就有了。

把所有的数据都带入到公式中，就很容易的把所有的参数都计算出来，接下去只需要再地图上画出结果就行。

那么这个椭圆揭示了一些什么意义呢？

使用ArcGIS的话，方向分布工具除了生成这样一个椭圆以外，还会给出如下结果：

其中，Shape_Leng和Shape_Area是生成的椭圆的周长和面积，单位与你数据的单位相同，这里我的数据是经纬度的，所以生成的结果只能作为相对参考结果。

CenterX和CenterY表示的是椭圆的中心点。

XstdDist和YStdDist表示的X轴的长度和Y轴的长度。

Rotation表示的是椭圆的方向角度。如下：

结果解读如下：

1、椭圆的长半轴表示的是数据分布的方向，短半轴表示的是数据分布的范围，长短半轴的值差距越大(扁率越大)，表示数据的方向性越明显。反之，如果长短半轴越接近，表示方向性越不明显。如果长短半轴完全相等，就等于是一个圆了，圆的话就表示没有任何的方向特征。

2、短半轴表示数据分布的范围，短半轴越短，表示数据呈现的向心力越明显；反之，短半轴越长，表示数据的离散程度越大。同样，如果短半轴与长半轴完全相等了，就表示数据没有任何的分布特征。

3、中心点表示了整个数据的中心位置，一般来说，只要数据的变异程度不是很大的话，这个中心点的位置大约与算数平均数的位置基本上是一致的，至于数据变异是什么情况，请看下面第4点。

4、有的同学会很疑惑，为什么你画的这个椭圆，还有很多的点都在外面，没有把所有的点都包含进去？那么就是就是“标准差椭圆”这个名词里面的“标准差”的含义所在了。

在ArcGIS工具里面(其他的工具也都差不多)，提供了“椭圆大小”(Ellipse_Size)这个参数，这个参数表示你生成的椭圆的级别，一共有三个，如下表：

三个级别的椭圆，分别表示了你生成的椭圆，能够包含68%，95%和99%三个级别的数据，我们通过可以指定要表示的标准差数(1、2 或 3)来决定你生成的椭圆包含的数据比例。

当要素具有空间正态分布时(即这些要素在中心处最为密集，而在接近外围时会逐渐变得稀疏)，第一级标准差(默认值)范围可将约占总数 68%的输入要素的质心包含在内。第二级标准差范围会将约占总数 95%的要素包含在内，而第三级标准差范围则会覆盖约占总数 99%的要素的质心。

所以，当你选择不同标准差等级的时候，你发现你的中心点的位置也可能不同。

当然，作为空间分析工具，方向分布一样可以进行加权计算，这个计算主要还是与中心点的位置确定以及椭圆标准差等级生成的椭圆大小有关系。

下面我们来通过一个实例来了解方向分布工具的应用：

一共有两年的伤寒病数据，如下，红色的是2000年的，蓝色是2001年的：

使用1个标准差的结果，生成的椭圆如上，具体数据如下：