深度学习已经占据了解决复杂问题的大多数领域，地理空间领域也不例外。文章的标题让您感兴趣，因此希望熟悉卫星数据集 ; 目前，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文件。
代码语言：javascript复制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数组读取。
打印一下读过的数据的大小。
代码语言：javascript复制print("Bangalore multispectral image shape: ", featuresBangalore.shape)
print("Bangalore binary built-up image shape: ", labelBangalore.shape)
print("Hyderabad multispectral image shape: ", featuresHyderabad.shape)输出：
代码语言：javascript复制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包的转换模块将做到这一点。
重组数据的架构
代码语言：javascript复制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)输出：
代码语言：javascript复制Bangalore multispectral image shape: 4198376, 6
Bangalore binary built-up image shape: 4198376
Hyderabad multispectral image shape: 1391808, 6在上面代码片段的第七行中，提取值为1的所有像素。这是一种故障安全措施，可以避免由于NoData像素导致的问题，这些像素通常具有极高和极低的值。
现在，将分割数据以进行训练和验证。这样做是为了确保模型没有看到测试数据，并且它对新数据的表现同样出色。否则，模型将过度拟合并且仅在训练数据上表现良好。
代码语言：javascript复制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)输出：
代码语言：javascript复制(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)输出：
代码语言：javascript复制(2519025, 1, 6) (1679351, 1, 6) (1391808, 1, 6)现在一切都已到位，使用keras构建模型。首先，将使用顺序模型，一个接一个地添加图层。有一个输入层，节点数等于nBands。使用具有14个节点和“ relu ”作为激活功能的一个隐藏层。最后一层包含两个节点，用于二进制构建类，具有' softmax '激活功能，适用于分类输出。
代码语言：javascript复制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上运行模型。根据数据大小和计算能力，安装模型需要一些时间。模型编译后可以看到以下内容：
预测单独保存的测试数据的值，并执行各种精度检查。
代码语言：javascript复制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。具有微小调整的类似模型可以应用于类似的应用。
代码语言：javascript复制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
