地學(xué)計算第四章課件

2022/10/111第四章變異函數(shù)和擬合模型資源與環(huán)境學(xué)院楊勇2022/10/112華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院目錄半方差函數(shù)影響半方差函數(shù)的因素變異函數(shù)的理論擬合模型變異函數(shù)（半方差函數(shù)）r(h)是地統(tǒng)計學(xué)的基石2022/10/113華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.1半方差函數(shù)—半方差云圖表示兩個實測值之間差的一半，稱為半方差值，一般意義上，半方差值會隨著兩點(diǎn)間距離的加大而加大，這是因為距離相近的樣品點(diǎn)的性質(zhì)較為相似。根據(jù)點(diǎn)對之間的空間滯后距離h將所有半方差值繪制成的散點(diǎn)圖稱為半方差云圖（如下圖，紅色圈內(nèi)為特異值）2022/10/114華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.1半方差函數(shù)—經(jīng)驗半方差函數(shù)在實踐中，一般是將滯后距（h）劃分為不同的級別，計算每個滯后級別上的平均半方差值，進(jìn)而得到經(jīng)驗半方差值和經(jīng)驗半方差圖來分析數(shù)據(jù)的空間關(guān)系。表示滯后級別，更為常見的表示方法為：2022/10/115華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院一般地，在實際計算時，假設(shè)步長為lag，當(dāng)前滯后級別為n(n為正整數(shù))，則h=n*lag，應(yīng)該這樣處理：（1）研究區(qū)所有點(diǎn)，找到點(diǎn)對（Pi,Pj），其符合條件：(n-1)*lag<dis(Pi,Pj)<=n*lag，它們之間的距離記為DISi（2）計算[z(pi)-z(pj)]2,記為Si.（3）設(shè)找到N(h)個這樣的點(diǎn)對，計算平均距離（4）計算，為n滯后級別上的經(jīng)驗半方差值。（5）將各個級別的（havg，r*(havg)），繪制在圖上，形成經(jīng)驗半方差圖2022/10/117華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院課堂練習(xí)假設(shè)給出某地土壤元素含量采樣數(shù)據(jù)，其格式為（x,y,k），共N個采樣點(diǎn)，其中，x和y分別表示采樣點(diǎn)的橫縱坐標(biāo)值，k為實驗室測得的元素含量值，步長為h，試給出相應(yīng)的程序步驟，計算該區(qū)域內(nèi)土壤元素的經(jīng)驗半方差圖。（需計算的最大滯后據(jù)為樣點(diǎn)間最大距離的一半）2022/10/118華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院2022/10/1110華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院變異函數(shù)的一些現(xiàn)象2、塊金效應(yīng)當(dāng)h=0時，變異函數(shù)r(h)不等于0，而是等于1個常數(shù)C0，這種現(xiàn)象稱為“塊金效應(yīng)”2022/10/1111華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院變異函數(shù)的一些現(xiàn)象3、各向異性2022/10/1112華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.2影響半方差函數(shù)的因素一、統(tǒng)計方面1、統(tǒng)計分布：特異值的影響，極大或極小值越靠近區(qū)域中間，影響越大。采樣數(shù)據(jù)的不均勻分布（h-散點(diǎn)圖）。2、屏蔽：特異值可能是錯誤的值，也可能是真實的值，可區(qū)分對待。如果檢查發(fā)現(xiàn)時錯誤。的值，必須去除。但如果是真實的值，可以去除，也可以保留（如污染區(qū)的污染物含量值）3、偏斜：數(shù)據(jù)的偏態(tài)分布，也會造成方差的不穩(wěn)定

異常值對變異函數(shù)的影響2022/10/1114華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.2影響半方差函數(shù)的因素三、采樣間距：應(yīng)大致地了解變異的空間尺度以選擇一個合適的采樣密度。（1）對于完全未知的領(lǐng)域，不知道變異的尺度或模式，應(yīng)預(yù)采樣，找到合適的采樣間距（2）若先前的采樣間距造成經(jīng)驗半方差函數(shù)呈純塊金效應(yīng)，沒有明顯的空間自相關(guān)性，則先前的采樣間距過大，應(yīng)縮小。（3）若依先前的采樣間距得到的半方差函數(shù)有明顯的結(jié)構(gòu)，但研究區(qū)內(nèi)樣點(diǎn)分布不均勻，則可進(jìn)一步優(yōu)化其分布。對于第（1）（2）情況，可以采取套合采樣解決。2022/10/1115華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.3變異函數(shù)的理論擬合模型理論變異函數(shù)用來擬合一些列經(jīng)驗變異函數(shù)值，供后續(xù)進(jìn)行插值估計時使用。選用理論變異函數(shù)模型是，要根據(jù)經(jīng)驗半方差圖的形狀來選取合適的模型2022/10/1117華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院4.3變異函數(shù)的理論擬合模型變異函數(shù)的理論模型：有基臺值模型無基臺值模型2022/10/1118華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院有基臺值模型—球狀模型C0：塊金常數(shù)C0+C：基臺值C：拱高a：變程應(yīng)用最廣的模型2022/10/1119華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院有基臺值模型—指數(shù)模型C0：塊金常數(shù)C0+C：基臺值C：拱高3a：變程當(dāng)C0=0,C=1時，稱為標(biāo)準(zhǔn)指數(shù)函數(shù)模型2022/10/1120華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院有基臺值模型—高斯模型C0：塊金常數(shù)C0+C：基臺值C：拱高：變程當(dāng)C0=0,C=1時，稱為標(biāo)準(zhǔn)高斯函數(shù)模型2022/10/1121華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院三種常用模型比較0.952022/10/1122華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院有基臺值模型—線性有基臺值模型C0：塊金常數(shù)C0+C：基臺值C：拱高A：常數(shù)，表示直線斜率當(dāng)C0=0,C=1時，稱為標(biāo)準(zhǔn)指數(shù)函數(shù)模型2022/10/1124華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院無基臺值模型——線性無基臺值模型2022/10/1125華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院無基臺值模型——冪函數(shù)值模型2022/10/1127華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型在實際中，有時區(qū)域化隨機(jī)變量Z(x)的變化相當(dāng)復(fù)雜，往往包含各種尺度及各種層次的變化，反映在變異函數(shù)r(h)上，就是單一的模型結(jié)構(gòu)不能將其合理表達(dá)，而是多層次的結(jié)構(gòu)相互疊加在一起，地統(tǒng)計學(xué)上稱為套合。所謂套合結(jié)構(gòu)，就是把分別出現(xiàn)在不同距離h上或不同方向上同時起作用的變異性組合起來，對全部有效的結(jié)構(gòu)信息，作定量化的概括，以表示區(qū)域化變量的主要特征。2022/10/1128華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型土壤是一個不均與、具有高度空間異質(zhì)性的復(fù)合體，它與土壤母質(zhì)、氣候、水文、地形和生物等因素有關(guān)，分析土壤空間變異的因素，可將其變異分為系統(tǒng)變異（土壤形成因素相互作用造成）和隨機(jī)變異（可以觀測到的，但與土壤形成印務(wù)無關(guān)且不能直接分析的）兩大類。如由h分開的兩個點(diǎn)x和x+h的土壤某一性質(zhì)Z(x)和Z(x+h)。當(dāng)h趨近于0時，可以認(rèn)為兩點(diǎn)間的差異完全是由取樣和測定誤差造成，當(dāng)h逐步增大，如h<1m，差異可能還要加上諸如水分等因素，當(dāng)h<100m時，在新的變異要考慮地形的作用。2022/10/1129華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型當(dāng)h一定時，變異函數(shù)r(h)應(yīng)包含小于h的所有影響因素，因此，絕大多數(shù)變異函數(shù)都由下面兩個變異函數(shù)組成：

r(h)=r0(h)+r1(h)，即一個代表純塊金方差，一個代表空間相關(guān)的方差。一般情況下，套合模型可以用放映各種不同尺度變化的多個變異函數(shù)之和表示，即：ri(h)可以是相同的或不同的理論模型2022/10/1130華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型如，區(qū)域化變量Z(x)的變異性由r0(h),r1(h)和r2(h)組成，其中2022/10/1131華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型三者組成的套合模型為：2022/10/1132華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型2022/10/1133華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型實例某區(qū)域鎳蘊(yùn)藏量的變異函數(shù)計算值2022/10/1134華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型實例鎳含量的變異函數(shù)曲線圖可以看出基臺值在2.55左右對頭兩點(diǎn)作切線與縱坐標(biāo)相交的交點(diǎn)就是塊金值，大約在0.4左右0.4/2.55=0.156則大約16%的誤差是隨機(jī)現(xiàn)象造成的，而空間自相關(guān)現(xiàn)象造成的誤差為84%2022/10/1135華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型實例C0=0.4C1=1.55a1=14C2=0.6a2=50從圖中可看出，理論值與實際值差異較大，尤其是在15到40m之間，因此，需進(jìn)行反復(fù)修改2022/10/1136華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院套合模型實例C0=0.4C1=1.15a1=12C2=1a2=60從圖中可看出，理論值與實際值差異擬合較好2022/10/1137華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院變異函數(shù)理論模型的最優(yōu)擬合要使變異函數(shù)的理論模型真實地描述變量的變化規(guī)律，在建立理論模型過程中，要對模型進(jìn)行最優(yōu)擬合，即確定模型中各個參數(shù)的值。地統(tǒng)計學(xué)中，擬合過程主要包括三個步驟：確定曲線類型、參數(shù)最優(yōu)估計、最優(yōu)曲線的確定。2022/10/1138華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院最優(yōu)擬合—確定曲線類型一般來講，通過(h,r(h))散點(diǎn)圖確定曲線的大致類型或初步類型，對這個大致的初步類型進(jìn)行參數(shù)最優(yōu)估計，確定是否為最優(yōu)曲線。2022/10/1139華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院最優(yōu)擬合—參數(shù)最優(yōu)估計變異函數(shù)的理論模型主要是曲線模型，將曲線模型經(jīng)過適當(dāng)?shù)淖儞Q，化為線性模型，然后用最小二乘法原理求未知參數(shù)的估計。2022/10/1140華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院最優(yōu)擬合—參數(shù)最優(yōu)估計上表式地統(tǒng)計學(xué)中常用的變異函數(shù)理論模型經(jīng)過適當(dāng)變換后，化為的線性模型，對于球狀模型、指數(shù)模型和高斯模型，只討論0<h<=a時的擬合問題，對變換后的變異函數(shù)理論模型，除球狀模型為二元線性回歸模型外，其余均為一元線性回歸模型，根據(jù)最小二乘法原理對這兩類線性回歸模型進(jìn)行參數(shù)估計。2022/10/1141華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院一元線性回歸模型的參數(shù)估計指數(shù)模型，高斯模型，冪函數(shù)模型和對數(shù)模型，均可以變換成一元線性回歸模型根據(jù)最小二乘法：通過得到的b0和b1，可以反算理論模型中的C0，C和a等參數(shù)值2022/10/1142華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院一元線性回歸模型的參數(shù)估計但上述一般方法沒有注意到變異函數(shù)曲線上的頭幾個點(diǎn)的重要性遠(yuǎn)大于其他點(diǎn)的重要性（因為后續(xù)預(yù)測算法中，與預(yù)測點(diǎn)距離較近的已知點(diǎn)對預(yù)測的結(jié)果影響較大），因此，不應(yīng)把它們與其他實際變異函數(shù)曲線上的點(diǎn)平均對待，為了克服這個問題，采用加權(quán)回歸的方法比較合適，擬合度較高，權(quán)重系數(shù)主要是采用每一距離上的樣本對數(shù)N(hi)2022/10/1143華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院加權(quán)回歸法其中：2022/10/1144華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院二元線性回歸模型的參數(shù)估計在地統(tǒng)計學(xué)的理論模型中，只有球狀模型線性化后成為二元線性回歸模型，共有三個參數(shù)：，三個參數(shù)加權(quán)最小二乘法的參數(shù)最優(yōu)估計為：2022/10/1145華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院二元線性回歸模型的參數(shù)估計2022/10/1146華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院計算出的b0,b1,b2后，還要分三種情況討論（1）b0>0,b1>0,b2<0,此時球狀模型中的三個參數(shù)C0,C,a為2022/10/1147華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院（2）如b0<0,b1>0,b2<0，此時b0<0，即C0<0，顯然不符合球狀模型的要求，可設(shè)b0=0，這是線性公式為y=b1x1+b2x2，重新根據(jù)最小二乘法求出參數(shù)b1和b2，在b0=0的條件下，仍可求出C0,C和a2022/10/1148華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院（3）如b0>0,b1>0,b2>=0,此時應(yīng)分兩種情況，一是b2=0，則線性模型變?yōu)閥=b0+b1x1，為一元線性模型，不是球狀模型，可按一元線性回歸模型參數(shù)估計方法求解參數(shù)，另一是b2>0，這時對原始數(shù)據(jù)進(jìn)行調(diào)整，增加或刪除一些不重要的實際變異函數(shù)點(diǎn)的數(shù)據(jù)，反復(fù)多次調(diào)整，直到b2<0為止，然后代入下述公式中進(jìn)行計算。2022/10/1149華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院回歸模型的檢驗請回顧第二章的相關(guān)內(nèi)容2022/10/1150華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院實例100個樣方土壤表層有機(jī)質(zhì)的含量

2022/10/1151華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院四種變異函數(shù)理論模型擬合參數(shù)線性有基臺值模型的殘差最小，決定系數(shù)最大，其次是球狀模型，基臺值基本相同，但線性有基臺值模型的變程最小，其次是高斯模型和球狀模型。塊金值是指數(shù)模型的最小，其次是球狀模型和線性有基臺值模型。對這5個參數(shù)，顯然最重要的是考慮決定系數(shù)R2的大小，其次是考慮殘差RSS的大小，然后再考慮變程和塊金值的大小，根據(jù)這個原則，選擇球狀模型作為本實例的變異函數(shù)理論模型是比較合適的，這個理論模型除了具有較高的擬合精度外，對變程內(nèi)的模擬可以得到滿意的結(jié)果。2022/10/1152華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院2022/10/1153華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于優(yōu)化搜索算法的參數(shù)擬合對于結(jié)構(gòu)復(fù)雜的變異函數(shù)理論模型，特別是套合結(jié)構(gòu)模型，參數(shù)復(fù)雜，難以用一般的通用方法求解出模型中的參數(shù)。但一些智能優(yōu)化算法，如遺傳算法、模擬退火算法、蟻群算法能夠使用統(tǒng)一的流程求解出接近最優(yōu)的參數(shù)。我們以遺傳算法為例，介紹該算法在求解套合結(jié)構(gòu)模型時的流程。2022/10/1154華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于遺傳算法的變異函數(shù)理論模型參數(shù)估計1、多尺度套合模型的規(guī)范表達(dá)2022/10/1155華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于遺傳算法的變異函數(shù)理論模型參數(shù)估計從上式可以看出，需求解的參數(shù)為2n+1個（因為第一個模型總是純塊金模型）。而在實際計算時，可以令，這樣方便從經(jīng)驗半方差圖中識別ci取值區(qū)間。并有以下約束：2022/10/1156華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于遺傳算法的變異函數(shù)理論模型參數(shù)估計編碼策略及初始群體產(chǎn)生假設(shè)需要顧及m(m<=2n+1)個參數(shù)，每個參數(shù)的取值范圍和估值精度分別是Umin,Umax和Qi，則將m個參數(shù)分別以L1,L2,……,Lm為長度進(jìn)行二進(jìn)制編碼，其中則每條染色體長度為，染色體中每個參數(shù)編碼對應(yīng)的解碼公式為：

以這種編碼方式隨機(jī)產(chǎn)生T組染色體2022/10/1157華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于遺傳算法的變異函數(shù)理論模型參數(shù)估計確定個體適應(yīng)度評價函數(shù)2022/10/1158華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院基于遺傳算法的變異函數(shù)理論模型參數(shù)估計遺傳操作遺傳算法主要包括3個基本算子，即選擇、交叉和變異，為此，需確定交叉概率Pc和變異概率Pm，3個過程執(zhí)行以后，將產(chǎn)生新一代種群，并記錄適應(yīng)度最高的染色體選擇—輪盤賭選擇法如下表表示11個個體的適應(yīng)度函數(shù)值、選擇概率和累計概率。每一輪產(chǎn)生1個[0,1]隨機(jī)數(shù)，將該隨機(jī)數(shù)作為選擇指針來確定被選個體。如隨機(jī)數(shù)為0.81，則第6個個體被選中。2022/10/1159華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院交叉算法單點(diǎn)交叉：考慮如下兩個11位變量的父個體父個體1：01110011010父個體2：10101100101交叉點(diǎn)位置為5，則交叉后兩個子個體：子個體1：01110100101子個體2：101010110102022/10/1160華中農(nóng)業(yè)大學(xué)資源與環(huán)境學(xué)院多點(diǎn)交叉M個交叉位置，可無重復(fù)隨機(jī)地選擇，交叉點(diǎn)之間的變量間續(xù)地相互交換，產(chǎn)生新個體，但在第一位置與第一個個交叉點(diǎn)之間的一段不做交換：父個體1：0111