当前位置：首页 > news >正文

通过卫星图像追踪城市扩张

news 2025/7/13 7:00:55

到本世纪中叶，非洲城市将成为城市化的中心。最佳估计显示，到2050年，非洲人口将翻一番，其中三分之二的增长将集中在城市。因此，各国和各城市亟需以促进经济发展的方式进行城市规划。然而，即使是像尼日利亚这样幅员辽阔的中低收入国家，也严重缺乏可靠且频繁的人口普查或人口普查，对城市发展的直接管理也十分有限。

卫星图像可以在这方面帮助我们，揭示发展随时间的变化，并识别出未来可能发展的地区，而这些地区比传统数据源覆盖的时间更长。卫星图像已证明其对发展经济学的实用性，因为它可以帮助我们研究关键的经济问题，即使在数据不易获取的环境中，也能依靠从太空获取的信息。

在这篇文章中，我将展示如何结合使用遥感和监督分类来观察 2000 年至 2020 年期间尼日利亚城市的发展。随机森林分类使我能够从基本卫星图像中识别土地覆盖类型，更具体地说是建筑环境，并使用乌贾瓦尔·甘地（Ujaval Gandhi）的空间思想中的工具，我们可以使用此分类来计算在研究的时间段内，哪些区域和土地利用因城市化进程而丧失。

寻找合适的卫星图像

市面上有大量卫星图像可供使用——访问和了解这些图像的一种简单方法是通过 GEE 代码编辑器的数据目录。对于每个栅格，GEE 数据目录都会识别数据中包含的特定波段，并指出所选图像中的潜在问题和需要的校正。

这里我使用了 LandSat 7，它可以追溯到 2000 年。为了获得特定年份数据的最佳图像，我会校正上述所有问题（例如校正表面反射率问题或去除云层），然后通过取这些图像的中值来聚合特定年份的图像。这样我就能获得尽可能清晰的图像。

运行监督分类

本次练习的目标是将上面的卫星图像转化为一张能够识别图像中每个 3000 万像素具体土地用途的图像。为此，我使用了监督分类，这意味着我根据现有的分类方法，告诉机器学习算法我希望在新图像中识别的不同类别。由于我这里只关注建筑环境，所以我只有 4 个类别：裸地（未开发的土地）、植被、城市和水域。

现有的土地覆盖数据通常有更具体的类别，例如美国国家土地覆盖数据库 (NBLDB) ，该数据库包含 16 个类别，其中 4 个为城市土地；或者，数据以更精细的空间和时间分辨率收集，例如每日 10 米分辨率的 Dynamic World 数据库。遗憾的是，像Dynamic World这样的数据来源仅从最近几年（2015 年以后）开始收集数据，而内容更丰富的土地覆盖地图更有可能覆盖高收入国家。

创建类别：
我没有依赖现有的训练数据，而是通过对 2020 年 Landsat 7 原始图像中每个感兴趣的类别（裸露、城市、水和植被）进行采样，手动创建自己的分类。为此，我在 GEE 地图窗口中使用图钉和/或小多边形；对于每个类别，我创建一个新的几何体，并为其分配一个特定的值，但为这 4 个类别的所有几何体赋予相似的类名。为了确保采样的准确性，我还比较了同一对象（此处为拉各斯市）在不同季节的变化。

在构建了裸露、城市、水和植被的 4 个几何图形后，它们可以合并为一个变量，供分类器使用。

图像改进：
该分类器是一个随机森林分类器，包含 50 棵树，使用合并后的几何体作为训练样本。在这里，我没有使用上面构建的 LandSat 7 中值图像以及我为每个类别收集的观测数据，而是在图像中添加了更多信息，以改进分类效果。

在本例中，由于我关注的是区分水、植被和建筑环境，因此我添加了由有助于突出这些不同特征的指数组成的附加波段，例如 NDVI、NDBI、mNDWI 和 BSI。这些指数是图像中包含的不同波段值之间的归一化差异：例如，NDVI（归一化差异植被指数）是 NIR（近红外）波段和红光波段值之间的归一化差异。它在遥感中最常用作绿度指标。这依赖于健康的植被会吸收红光波段并反射近红外波段的事实，因此指数值较高（接近 1）表示植被。其他归一化差异使用相同原理来识别水（改进的 NDWI）和建筑环境（NDBI 和建筑表面指数，后者也是归一化差异，但结合了多个波段：一边是近红外和蓝光，另一边是红光和短波红外 1）。

此外，还可以将追踪其他参数（例如坡度）的波段添加到待分类的图像中。所有这些都添加到中值图像的原始波段上，生成的图像在下方表示为“合成”

var addIndices = function(image){var ndvi = imagenormalizedDifference(["SR_B4","SR_B3"]).rename(["ndvi"]);var ndbi = imagenormalizedDifference(["SR_B5","SR_B4"]).rename(["ndbi"]);var mndwi = imagenormalizedDifference(["SR_B2","SR_B5"]).rename(["mndwi"]);var bsi = imageexpression("(( X + Y ) - (A + B)) /(( X + Y ) + (A + B))",{'X': imageselect("SR_B5"),//swir1'Y': imageselect("SR_B3"),//red'A': imageselect("SR_B4"),// nir'B': imageselect("SR_B1"),// blue}).rename("bsi");return imageaddBands(ndvi).addBands(ndbi).addBands(mndwi).addBands(bsi);};var composite = addIndices(landsat7_2020);

运行分类器
下面的代码可以更好地解释这一点，但我们使用机器学习（下文将介绍一个随机的第一分类器）来检测土地覆盖。使用我们上面构建的训练样本，并基于添加了附加波段的起始图像，随机森林分类器将尝试在新图像中检测（可能会有一些误差）与训练样本中识别出的相同类别。对于研究人员来说，这里唯一的选择变量是树木的数量，其中每棵树都可以被视为从图像中提取的数据样本，用于应用算法。

// Create a random forest classifier with 50 trees.var classifier = eeClassifiersmileRandomForest({numberOfTrees:50}).train({// Train the classifier.// Use the sampled pixels.features: training,// This is the class that we want to predict.classProperty:"class",// The bands the classifier will use for training and classification.inputProperties: compositebandNames()});// Apply the classifier on the 2020 image.var classified20 = compositeclassify(classifier);

运行分类器后，我们可以创建一个显示分类图像的新地图。最终结果如下所示，其中红色表示建筑环境，绿色表示植被，蓝色表示来自 2020 年合成图像的水。由于我们使用 2020 年拉各斯的像素来训练分类器，因此这种分类不如观察 2010 年和 2000 年发生了什么以及建筑环境发生了哪些变化那么有趣。然后，我们可以应用相同的分类器对 2010 年和 2000 年的图像进行分类。

下图基于随机森林分类和 2020 年训练样本，展示了 2010 年（上）和 2000 年（下）的拉各斯。

为了更系统地观察变化，我们可以绘制两张地图之间的差异图。回想一下，由于分类的不同，每个像素对应一个特定的值（裸露、城市、水体、植被）。有多种方法可以做到这一点：最简单的方法是将 2010 年和 2020 年之间的差异制成表格并绘制成图。我们可以得到一张地图，只显示 2010 年至 2020 年期间分类发生变化的区域。

由于这并不能告诉我们太多关于建筑环境变化的信息，甚至不能告诉我们哪些类型的土地用途已经转变为城市用地，我们可以尝试精确地确定哪些土地随着时间的推移已经变成了城市用地（参见下面的代码）。我们可以进一步使用ee.Reducer.frequencyHistogram()函数创建一个转换矩阵，但需要注意的是，该矩阵会返回一个计数值（即从一个类别转换为另一个类别的像素数量），而这可能不是有用的信息。

## since we had previously assigned the value of 0 to one of the categories, we reassign new values herevar classified10_new = classified10remap([0,1,2,3],[1,2,3,4]);var classified20_new = classified20remap([0,1,2,3],[1,2,3,4]);## this allows us to keep track of class in 2010 and in 2020, by multiply the class value in 2010 by 100## the first digit gives us class in 2010 and the last one class in 2020var envt2020 = classified10_newmultiply(100).add(classified20_new).rename("transitions20");## we use this new variable to identify land that was not urban in 2010 (first number is not 2) and has become urban by 2020 (units digit is equal to 2)var urbandiff20_10 = envt2020eq(102).or(envt2020eq(302)).or(envt2020eq(402));## this second variable tracks the land that remained urban between 2010 and 2020var urban20102020 = envt2020eq(202);## this adds the map to visualize the changes, as seen belowMapaddLayer(urbandiff20_10clip(lagos),{min:0,max:1,palette:["white","#7a0177"]},"newly urban - in 2020");