matlab 误差椭圆,求3倍标准差误差椭圆分析的程序
根據《白話空間統計之九:方向分布(標準差橢圓)修正版》(有些地方沒有理解清楚),寫了下面的程序。但是好像結果不對
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個標準差的結果,生成的橢圓如上,具體數據如下:
總結
以上是生活随笔為你收集整理的matlab 误差椭圆,求3倍标准差误差椭圆分析的程序的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: AI公司奥创光年Mogic AI获千万美
- 下一篇: 用matlab做纹理合成,关于图像纹理合