GEE代码学习 day9

7.2 超参数优化

评估随机森林分类器中的树数量如何影响分类准确性

// Hyperparameter tuning. 
var numTrees = ee.List.sequence(5, 100, 5);  
var accuracies = numTrees.map(function(t) { 
                   var classifier = ee.Classifier.smileRandomForest(t) 
                   .train({ features: trainingSet, 
                    classProperty: 'class', 
                    inputProperties: predictionBands }); 
return testingSet .classify(classifier) .errorMatrix('class', 'classification') .accuracy(); });  
print(ui.Chart.array.values({ array: ee.Array(accuracies), 
                              axis: 0, 
                            xLabels: numTrees }) 
                             .setOptions({ hAxis: { title: 'Number of trees' }, 
                                           vAxis: { title: 'Accuracy' }, 
                                           title: 'Accuracy per number of trees' }));

8 解释图像：回归

回归是可用于从查看影像转变为分析影像的基本工具之一。

确定因变量与一个或多个自变量之间关系的强度和特征，以更好地理解或预测因变量。在遥感环境中，这些因变量可能是自然灾害或物种发生的可能性，而自变量可能是波段、指数或其他栅格数据集，如树高或土壤类型。我们对数值因变量使用回归。分类用于分类因变量，例如土地覆被或土地利用。

过程

1. 收集有关因变量的已知实例的数据。这可能是已知的滑坡区域、已知的物种出现区域等。

2. 定义并选择自变量。

3. 运行回归以生成描述因变量和自变量之间关系的方程。

4. 可选：使用此方程和自变量图层，您可以创建整个区域的因变量地图。

使用遥感数据的主要目标是汇总空间、时间和波段信息以查找关系。Earth Engine 为遥感数据提供了大量汇总技术，这些技术属于统称为 reducers。

一组数字的缩减器的一个示例是 median，它在集合中找到中间数字，从而将大型集合的复杂性降低为具有代表性的估计值。该减少操作在 Earth Engine 中使用 ee 实现。Reducer.median. 的 Reducer.median. 中，或者对于这个特别常见的操作，使用简写操作 median。

Reducer 可以逐个像素运行，也可以在感知像素的周围环境的情况下运行。空间对象的 reducer 的一个示例是计算邻域的中位数高程，这可以通过嵌入对 ee 的调用在 Earth Engine 中实现。Reducer.median 函数中。reduceRegion 调用指示 Earth Engine 分析跨区域分组的像素，从而允许将城市街区或公园边界等区域视为统一的兴趣单元。Reducer 可用于汇总 Earth Engine 中的多种形式的数据，包括图像集合、图像、列表和特征集合。

Reducer 可以相对简单，如中位数，也可以更复杂。想象一个 1000 人的身高和体重样本：线性回归会将 2000 个组装值减少到两个：斜率和截距。当然，在给定的集合中，这种关系可能是弱的或强的;更大的想法是，许多简化和总结大量数据的操作都可以使用 reducer 的这个想法和术语来概念化。

线性拟合

Earth Engine 中可用的最简单的回归是通过 reducer linearFit 实现的。此函数是具有一个自变量和一个因变量的线性函数的最小二乘估计值。这个回归方程写成 Y = α + β X + ，其中 α 是线的截距，β 是斜率，Y 是因变量，X 是自变量，是误差。假设目标是根据有关意大利都灵的已知信息估计每个 Landsat 像素中的树木覆盖百分比。下面，我们将定义一个围绕都灵的矩形几何，并将其命名为都灵。

// Define a Turin polygon. 
var Turin = ee.Geometry.Polygon(  [ 
                                    [  
                                       [7.455553918110218, 45.258245019259036], 
                                       [7.455553918110218, 44.71237367431335], 
                                       [8.573412804828967, 44.71237367431335], 
                                        [8.573412804828967, 45.258245019259036] 
                                                                                  ] 
                                                                                ], 
                                     null, 
                                    false);  
// Center on Turin 
Map.centerObject(Turin, 9);

访问使用“MOD44B.006 Terra Vegetation Continuous Fields Yearly Global 250 m”数据集以用作我们的因变量。

var mod44b = ee.ImageCollection('MODIS/006/MOD44B');

访问 ImageCollection 的那部分，该部分表示 2020 年分辨率为 250 m 的全球树木覆盖率。我们将使用 filterDate 命令访问该时间段的图像，使用first命令将该单个图像转换为 Image 类型，将其剪辑到 Turin，然后选择适当的波段。然后，我们将有关新映像 percentTree2020 的信息打印到控制台。

///// 
// Start Linear Fit /////  
// Put together the dependent variable by filtering the 
// ImageCollection to just the 2020 image near Turin and 
// selecting the percent tree cover band. 
var percentTree2020 = mod44b  .filterDate('2020-01-01', '2021-01-01') .first() .clip(Turin) .select('Percent_Tree_Cover');  
// You can print information to the console for inspection.  
print('2020 Image', percentTree2020);  
Map.addLayer(percentTree2020, { max: 100 }, 'Percent Tree Cover');

使用“USGS Landsat 8 Collection 2 Tier 1 and Real-Time Data Raw Scenes”用作我们的自变量。

var landsat8_raw =  
ee.ImageCollection('LANDSAT/LC08/C02/T1_RT');

按日期和位置筛选此集合。我们会将集合筛选到 2020 年的明确（无云）日期，然后按位置进行筛选。

// Put together the independent variable. 
var landsat8filtered = landsat8_raw  
                      .filterBounds(Turin.centroid({ 'maxError': 1 })) 
                      .filterDate('2020-04-01', '2020-4-30') .first();  
print('Landsat8 filtered', landsat8filtered);  
// Display the L8 image. 
var visParams = {  bands: ['B4', 'B3', 'B2'], max: 16000 }; 
Map.addLayer(landsat8filtered, visParams, 'Landsat 8 Image');

使用 centroid 函数将选择与 Turin 几何中心相交的图像，而不是几何中的任何位置。（注意，我们在本练习中使用单个图像，但在实践中，您几乎总是需要将图像集合过滤到您感兴趣区域的边界，然后创建合成。在这种情况下，您将使用合成 Earth Engine 算法来获得 Landsat 影像的无云合成。）

使用 Landsat 影像的波段计算 NDVI，我们将将其用作自变量。

// Calculate NDVI which will be the independent variable. 
var ndvi = landsat8filtered.normalizedDifference(['B5',  'B4']);

现在，我们开始将 Earth Engine 中的数据组合成正确的格式。首先，使用 addBands 函数创建具有两个波段的图像：首先，使用NDVI，它将充当自变量;第二个是之前创建的映象 percentTree2020。

// Create the training image. 
var trainingImage = ndvi.addBands(percentTree2020);  
print('training image for linear fit', trainingImage);

现在我们可以设置和运行回归，在我们的几何上使用线性拟合 reducer 并打印结果。由于构建回归模型需要从图像周围组装点以供考虑，因此它是使用 reduceRegion 而不是 reduce 实现的。我们需要包括比例变量（此处为 30 m，这是 Landsat 的分辨率）。

// Independent variable first, dependent variable second. 
// You need to include the scale variable. 
var linearFit = trainingImage.reduceRegion({ 
                                      reducer: ee.Reducer.linearFit(),
                                      geometry: Turin, 
                                      scale: 30, 
                                      bestEffort: true  });  
// Inspect the results. 
print('OLS estimates:', linearFit); 
print('y-intercept:', linearFit.get('offset')); 
print('Slope:', linearFit.get('scale'));

如果您尝试在没有 bestEffort： true 参数的情况下运行此操作，则很可能会收到错误：“Image.reduceRegion： Too many pixels in the region.”这意味着所涉及的像素数已超过 Earth Engine 的默认 maxPixels 限制 1000 万。原因是我们请求的回归的复杂性。如果您曾经使用 100 个点进行过回归，您可能已经看到或制作了一个包含 100 个点的散点图;现在，想象一个包含超过 1000 万个点的散点图，以设想所请求内容的规模。在这里，限制不是 Earth Engine 的计算能力，而更像是一个通知，表明代码正在调用特别大的计算，并且新手用户可能并不完全打算执行这种复杂的操作。错误文本指出了解决此问题的几种方法。首先，我们可以增加 maxPixels。其次，我们可以以较低的分辨率进行聚合（例如，将比例尺从 30 m 增加到 50 m）。第三，我们可以像这里一样将 bestEffort 参数设置为 true，这会指示 reducer 使用最高分辨率扩展而不超过maxPixels。第四，我们可以减少该区域的面积（即，缩小都灵周围的研究区域）。最后，我们可以对图像堆栈进行随机采样，并将该样本用于回归。

最后，让我们将回归应用于整个 NDVI 区域，它比都灵边界多边形大。

// Create a prediction based on the linearFit model. 
var predictedTree = ndvi.expression(  'intercept + slope * ndvi', { 'ndvi': ndvi.select('nd'), 'intercept': ee.Number(linearFit.get('offset')), 'slope': ee.Number(linearFit.get('scale')) });  
print('predictedTree', predictedTree);  
// Display the results. 
Map.addLayer(predictedTree, { max: 100 }, 'Predicted Percent Tree Cover');