基於Python的Tensorflow衛星數據分類神經網絡

  • 2019 年 10 月 4 日
  • 筆記

深度學習已經佔據了解決複雜問題的大多數領域,地理空間領域也不例外。文章的標題讓您感興趣,因此希望熟悉衛星數據集 ; 目前,Landsat 5 TM。機器學習(ML)算法如何工作的知識很少,將幫助快速掌握這本動手教程。對於那些不熟悉ML概念的人,簡而言之,它是建立一個實體的一些特徵(特徵或X)與其他屬性(值或標籤或Y)之間的關係 – 提供了大量的例子(標記數據) )到模型,以便從中學習,然後預測新數據(未標記數據)的值/標籤。這對於機器學習來說已經足夠理論了!

衛星數據的一般問題:

衛星數據中的兩個或更多要素類(例如,建造/貧瘠/採石場)可具有相似的光譜值,這使得該分類在過去的幾十年中成為具有挑戰性的任務。由於上述問題,傳統的監督和無監督方法不能成為完美的分類器,儘管它們可以穩健地執行分類。但是總會有相關的問題。通過以下示例來理解這一點:

在上圖中,如果使用垂直線作為分類器並僅沿着x軸移動它,使其將所有圖像分類為右側作為房屋,則答案可能不是直截了當的。這是因為數據的分佈是這樣的,即不可能只用一條垂直線將它們分開。但是,這並不意味着房屋根本無法歸類!

假設使用紅線(如上圖所示)來分隔這兩個要素。在這種情況下,大多數房屋都是由分類器確定的,但房子仍被遺漏,一棵樹被誤分類為房屋。為了確保不會留下任何一個房子,可以使用藍線。在這種情況下,分類器將覆蓋所有房屋; 這被稱為高召回率。然而,並非所有的分類圖像都是真正的房屋,這被稱為低精度。同樣,如果使用綠線,所有分類為房屋的圖像都是房屋; 因此,分類器具有高精度。在這種情況下召回的次數會減少,因為還有三所房子被遺漏了。在大多數情況下,這種權衡 在精確度和召回之間保持。

上面展示的房屋和樹木問題類似於建築物,採石場和貧瘠土地的情況。衛星數據的分類優先級可能因目的而異。例如,如果想確保所有的組合單元被歸類為組合,沒有留下任何東西,並且你更少關心具有類似簽名的其他類的像素被歸類為組合,那麼一個模型與需要高召回率。相反,如果優先級是僅對純組合像素進行分類而不包括任何其他類像素,並且可以放棄混合的組合像素,則需要高精度分類器。通用模型將使用房屋和樹木的紅線來保持精確度和召回之間的平衡。

當前範圍中使用的數據

在這裡,將把Landsat 5 TM的六個波段(波段2 – 波段7)視為特徵,並嘗試預測二進制構建類。2011年為班加羅爾及其相應的二元建築層獲得的多光譜Landsat 5數據將用於訓練和測試。最後,2005年為海德拉巴收購的另一個多光譜Landsat 5數據將用於新的預測。

由於使用標記數據來訓練模型,因此這是一種受監督的ML方法。

多光譜訓練數據及其相應的二進制構建層

將在Python中使用Google的Tensorflow庫來構建神經網絡(NN)。需要以下其他庫,請確保提前安裝它們:

  • pyrsgis – 讀寫 GeoTIFF
  • scikit-learn – 用於數據預處理和準確性檢查
  • numpy – 用於基本數組操作

沒有進一步的延遲,開始編碼。

將所有三個文件放在一個目錄中 – 在腳本中分配路徑和輸入文件名,並讀取GeoTIFF文件。

import os  from pyrsgis import raster    os.chdir("E:\yourDirectoryName")  mxBangalore = 'l5_Bangalore2011_raw.tif'  builtupBangalore = 'l5_Bangalore2011_builtup.tif'  mxHyderabad = 'l5_Hyderabad2011_raw.tif'    # Read the rasters as array  ds1, featuresBangalore = raster.read(mxBangalore, bands='all')  ds2, labelBangalore = raster.read(builtupBangalore, bands=1)  ds3, featuresHyderabad = raster.read(mxHyderabad, bands='all')

pyrsgis包的柵格模塊分別將GeoTIFF的地理位置信息和數字編號(DN)值作為NumPy數組讀取。

打印一下讀過的數據的大小。

print("Bangalore multispectral image shape: ", featuresBangalore.shape)  print("Bangalore binary built-up image shape: ", labelBangalore.shape)  print("Hyderabad multispectral image shape: ", featuresHyderabad.shape)

輸出:

Bangalore multispectral image shape: 6, 2054, 2044  Bangalore binary built-up image shape: 2054, 2044  Hyderabad multispectral image shape: 6, 1318, 1056

從輸出中可以看出,班加羅爾圖像中的行數和列數是相同的,並且多光譜圖像中的層數是相同的。該模型將基於所有頻帶上的相應DN值來學習確定像素是否構建,因此,多光譜圖像應具有以相同順序堆疊的相同數量的特徵(頻帶)。

現在將數組的形狀更改為二維數組,這是大多數ML算法所期望的,其中每行代表一個像素。pyrsgis包的轉換模塊將做到這一點。

重組數據的架構

from pyrsgis.convert import changeDimension    featuresBangalore = changeDimension(featuresBangalore)  labelBangalore = changeDimension (labelBangalore)  featuresHyderabad = changeDimension(featuresHyderabad)  nBands = featuresBangalore.shape[1]  labelBangalore = (labelBangalore == 1).astype(int)    print("Bangalore multispectral image shape: ", featuresBangalore.shape)  print("Bangalore binary built-up image shape: ", labelBangalore.shape)  print("Hyderabad multispectral image shape: ", featuresHyderabad.shape)

輸出:

Bangalore multispectral image shape: 4198376, 6  Bangalore binary built-up image shape: 4198376  Hyderabad multispectral image shape: 1391808, 6

在上面代碼片段的第七行中,提取值為1的所有像素。這是一種故障安全措施,可以避免由於NoData像素導致的問題,這些像素通常具有極高和極低的值。

現在,將分割數據以進行訓練和驗證。這樣做是為了確保模型沒有看到測試數據,並且它對新數據的表現同樣出色。否則,模型將過度擬合併且僅在訓練數據上表現良好。

from sklearn.model_selection import train_test_split    xTrain, xTest, yTrain, yTest = train_test_split(featuresBangalore, labelBangalore, test_size=0.4, random_state=42)    print(xTrain.shape)  print(yTrain.shape)    print(xTest.shape)  print(yTest.shape)

輸出:

(2519025, 6)  (2519025,)  (1679351, 6)  (1679351,)

上面代碼片段中的test_size(0.4)表示訓練測試比例為60/40。

包括NN在內的許多ML算法都期望歸一化數據。這意味着直方圖在一定範圍(此處為0到1)之間被拉伸和縮放。將規範化功能以滿足此要求。可以通過減去最小值併除以範圍來實現歸一化。由於Landsat數據是8位數據,因此最小值和最大值分別為0和255(2⁸= 256個值)。

請注意,從標準化數據計算最小值和最大值始終是一個好習慣。為避免複雜性,將在此處堅持使用8位數據的默認範圍。

另一個額外的預處理步驟是將特徵從二維重塑為三維,使得每行代表單個像素。

# Normalise the data  xTrain = xTrain / 255.0  xTest = xTest / 255.0  featuresHyderabad = featuresHyderabad / 255.0    # Reshape the data  xTrain = xTrain.reshape((xTrain.shape[0], 1, xTrain.shape[1]))  xTest = xTest.reshape((xTest.shape[0], 1, xTest.shape[1]))  featuresHyderabad = featuresHyderabad.reshape((featuresHyderabad.shape[0], 1, featuresHyderabad.shape[1]))    # Print the shape of reshaped data  print(xTrain.shape, xTest.shape, featuresHyderabad.shape)

輸出:

(2519025, 1, 6) (1679351, 1, 6) (1391808, 1, 6)

現在一切都已到位,使用keras構建模型。首先,將使用順序模型,一個接一個地添加圖層。有一個輸入層,節點數等於nBands。使用具有14個節點和「 relu 」作為激活功能的一個隱藏層。最後一層包含兩個節點,用於二進制構建類,具有' softmax '激活功能,適用於分類輸出。

from tensorflow import keras    # Define the parameters of the model  model = keras.Sequential([      keras.layers.Flatten(input_shape=(1, nBands)),      keras.layers.Dense(14, activation='relu'),      keras.layers.Dense(2, activation='softmax')])    # Define the accuracy metrics and parameters  model.compile(optimizer="adam", loss="sparse_categorical_crossentropy", metrics=["accuracy"])    # Run the model  model.fit(xTrain, yTrain, epochs=2)

神經網絡架構

如第10行所述,使用' adam '優化器編譯模型。(還有其他幾個你可以檢查。)現在將使用的損失類型是分類 – 稀疏 – 交叉熵。模型性能評估的度量標準是「 準確性 」。

最後,使用兩個時期(或迭代)在xTrain和yTrain上運行模型。根據數據大小和計算能力,安裝模型需要一些時間。模型編譯後可以看到以下內容:

預測單獨保存的測試數據的值,並執行各種精度檢查。

from sklearn.metrics import confusion_matrix, precision_score, recall_score    # Predict for test data  yTestPredicted = model.predict(xTest)  yTestPredicted = yTestPredicted[:,1]    # Calculate and display the error metrics  yTestPredicted = (yTestPredicted>0.5).astype(int)  cMatrix = confusion_matrix(yTest, yTestPredicted)  pScore = precision_score(yTest, yTestPredicted)  rScore = recall_score(yTest, yTestPredicted)    print("Confusion matrix: for 14 nodesn", cMatrix)  print("nP-Score: %.3f, R-Score: %.3f" % (pScore, rScore))

所述SOFTMAX函數產生用於每個類類型單獨的列的概率值。只提取第一類(構建),如上面代碼片段中的第六行所述。用於地理空間相關分析的模型變得難以評估,因為與其他一般ML問題不同,依賴於廣義總結誤差是不公平的; 空間位置是獲勝模型的關鍵。因此,混淆矩陣,精度和召回可以更清晰地反映模型的表現。

終端中顯示的混淆矩陣,精度和召回

如上面的混淆矩陣所示,有數千個組合像素被分類為非組合,反之亦然,但與總數據大小的比例較小。在測試數據上獲得的精度和召回率大於0.8。

總是可以花一些時間並執行一些迭代來找到隱藏層的最佳數量,每個隱藏層中的節點數以及獲得準確性的時期數。一些常用的遙感指數,如NDBI或NDWI,也可以在需要時用作特徵。達到所需精度後,使用模型預測新數據並導出GeoTIFF。具有微小調整的類似模型可以應用於類似的應用。

predicted = model.predict(feature2005)  predicted = predicted[:,1]    #Export raster  prediction = np.reshape(predicted, (ds.RasterYSize, ds.RasterXSize))  outFile = 'Hyderabad_2011_BuiltupNN_predicted.tif'  raster.export(prediction, ds3, filename=outFile, dtype='float')

請注意,使用預測的概率值導出GeoTIFF,而不是其閾值化的二進制版本。總是可以在以後的GIS環境中對浮點類型圖層進行閾值處理,如下圖所示。

Hyderabad構建層由模型使用多光譜數據預測

已經精確評估了模型的準確性並進行了調用 – 還可以對新預測的柵格進行傳統檢查(例如kappa係數)。除了上述衛星數據分類的挑戰之外,其他直觀的限制包括由於光譜特徵的變化,模型無法預測在不同季節和不同區域獲得的數據。

在本文中使用的模型是NN的一個非常基本的架構,包括卷積神經網絡(CNN)在內的一些複雜模型已經被研究人員證明可以產生更好的結果。這種分類的主要優點是一旦模型被訓練就具有可擴展性。

請在此處查找使用的數據和完整腳本。

https://github.com/PratyushTripathy/Landsat-Classification-Using-Neural-Network