Google Earth Engine（GEE）——影像数组的转换与运算.matrixSolve

Earth Engine 支持转置、逆和伪逆等数组变换。例如，考虑一个时间序列图像的普通最小二乘 (OLS) 回归。在以下示例中，具有预测变量和响应的带的图像被转换为​​数组图像，然后“求解”以获得最小二乘系数估计三种方式。首先，组装图像数据并转换为数组：

这里对于计算系数有三种方式：下面代码中有所展示，具体就是更好的方法是使用该pseudoInverse()方法（matrixPseudoInverse()对于数组图像）：

matrixPseudoInverse()

计算矩阵的 Moore-Penrose 伪逆。
 

Computes the Moore-Penrose pseudoinverse of the matrix.

Arguments:

this:value (Image):

应用操作的图像。
 

The image to which the operation is applied.

Returns: Image

matrixSolve(image2)

求解矩阵方程 A*x=B 中的 x，如果对于 image1 和 image2 中每个匹配的波段对 A 是超过估计的，则找到最小二乘解。如果 image1 或 image2 只有 1 个波段，则将其用于另一个图像中的所有波段。如果图像具有相同数量的波段，但名称不同，则它们按自然顺序成对使用。输出波段以两个输入中较长的命名，或者如果它们的长度相等，则按 image1 的顺序命名。输出像素的类型是输入类型的并集。
 

Solves for x in the matrix equation A*x=B, finding a least-squares solution if A is overdetermined for each matched pair of bands in image1 and image2. If either image1 or image2 has only 1 band, then it is used against all the bands in the other image. If the images have the same number of bands, but not the same names, they're used pairwise in the natural order. The output bands are named for the longer of the two inputs, or if they're equal in length, in image1's order. The type of the output pixels is the union of the input types.

Arguments:

  

   

    

     

    
    

     

      图像1（图像）：
     
    
   

    

     

    
    

     

      从中获取左操作数带的图像。
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      图像2（图像）：
     
    
   

    

     

    
    

     

      从中获取正确操作数带的图像。
     
    
  
 

this:image1 (Image):

The image from which the left operand bands are taken.

image2 (Image):

The image from which the right operand bands are taken.

Returns: Image

或者

从可读性和计算效率的角度来看，获得 OLS 系数的最佳方法是solve()（matrixSolve()对于阵列图像）。该 solve()函数确定如何从输入的特征中最好地求解系统，使用超定系统的伪逆、方阵的逆和近似奇异矩阵的特殊技术：

代码：

  

   

    

     

    
    

     

      // 此函数使用简单云分数的阈值屏蔽输入。
     
    
   

    

     

    
    

     

      var cloudMask = function(img) {
     
    
   

    

     

    
    

     
  var cloudscore = ee.Algorithms.Landsat.simpleCloudScore(img).select('cloud');
     
    
   

    

     

    
    

     
  return img.updateMask(cloudscore.lt(50));
     
    
   

    

     

    
    

     

      };
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 加载Landsat5影像，对于以后影像操作基本上下面的4-5步逗得用
     
    
   

    

     

    
    

     

      var collection = ee.ImageCollection('LANDSAT/LT05/C01/T1_TOA')
     
    
   

    

     

    
    

     
  // 筛选08-10年的影像
     
    
   

    

     

    
    

     

        .filterDate('2008-04-01', '2010-04-01')
     
    
   

    

     

    
    

     
  // 过滤以仅获取感兴趣点的图像。
     
    
   

    

     

    
    

     

        .filterBounds(ee.Geometry.Point(-122.2627, 37.8735))
     
    
   

    

     

    
    

     
  // 通过在集合上遍历 cloudMask 函数来屏蔽云。
     
    
   

    

     

    
    

     

        .map(cloudMask)
     
    
   

    

     

    
    

     
  // 仅选择 NIR 和红色波段。
     
    
   

    

     

    
    

     

        .select(['B4', 'B3'])
     
    
   

    

     

    
    

     
  // 按时间顺序对集合进行排序。
     
    
   

    

     

    
    

     

        .sort('system:time_start', true);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 此函数根据输入计算预测变量和响应。
     
    
   

    

     

    
    

     

      var makeVariables = function(image) {
     
    
   

    

     

    
    

     
  // 计算图像相对于 Epoch 的小数年时间。
     
    
   

    

     

    
    

     
  var year = ee.Image(image.date().difference(ee.Date('1970-01-01'), 'year'));
     
    
   

    

     

    
    

     
  // 以弧度计算季节，每年一个周期。
     
    
   

    

     

    
    

     
  var season = year.multiply(2 * Math.PI);
     
    
   

    

     

    
    

     
  // 返回预测变量的图像，然后是响应。
     
    
   

    

     

    
    

     
  return image.select()
     
    
   

    

     

    
    

     

          .addBands(ee.Image(1))                                  // 0. constant
     
    
   

    

     

    
    

     

          .addBands(year.rename('t'))                             // 1. linear trend
     
    
   

    

     

    
    

     

          .addBands(season.sin().rename('sin'))                   // 2. seasonal
     
    
   

    

     

    
    

     

          .addBands(season.cos().rename('cos'))                   // 3. seasonal
     
    
   

    

     

    
    

     

          .addBands(image.normalizedDifference().rename('NDVI'))  // 4. response
     
    
   

    

     

    
    

     

          .toFloat();
     
    
   

    

     

    
    

     

      };
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 定义集合数组中的变化轴。
     
    
   

    

     

    
    

     

      var imageAxis = 0;
     
    
   

    

     

    
    

     

      var bandAxis = 1;
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 将集合转换为数组。
     
    
   

    

     

    
    

     

      var array = collection.map(makeVariables).toArray();
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 检查图像轴的长度（图像数）。
     
    
   

    

     

    
    

     

      var arrayLength = array.arrayLength(imageAxis);
     
    
   

    

     

    
    

     

      // 更新msak以确保图像数量大于或
     
    
   

    

     

    
    

     

      // 等于预测变量的数量（线性模型是可解的）
     
    
   

    

     

    
    

     

      array = array.updateMask(arrayLength.gt(4));
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 根据沿带轴的位置获取数组的切片。
     
    
   

    

     

    
    

     

      var predictors = array.arraySlice(bandAxis, 0, 4);
     
    
   

    

     

    
    

     

      var response = array.arraySlice(bandAxis, 4);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 以困难的方式计算系数。不推荐
     
    
   

    

     

    
    

     

      var coefficients1 = predictors.arrayTranspose().matrixMultiply(predictors)
     
    
   

    

     

    
    

     

        .matrixInverse().matrixMultiply(predictors.arrayTranspose())
     
    
   

    

     

    
    

     

          .matrixMultiply(response);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      //以简单的方式计算
     
    
   

    

     

    
    

     

      // 推荐使用这种犯法
     
    
   

    

     

    
    

     

      var coefficients2 = predictors.matrixPseudoInverse()
     
    
   

    

     

    
    

     

        .matrixMultiply(response);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 计算系数最简单的方法。强烈推荐
     
    
   

    

     

    
    

     

      var coefficients3 = predictors.matrixSolve(response);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      // 将结果转换为多波段图像。
     
    
   

    

     

    
    

     

      var coefficientsImage = coefficients3
     
    
   

    

     

    
    

     
  // Get rid of the extra dimensions.
     
    
   

    

     

    
    

     

        .arrayProject([0])
     
    
   

    

     

    
    

     

        .arrayFlatten([
     
    
   

    

     

    
    

     

          ['constant', 'trend', 'sin', 'cos']
     
    
   

    

     

    
    

     

      ]);
     
    
   

    

     

    
    

     
 
     
    
   

    

     

    
    

     

      //自己可以答应输出结果看看，发现结果一致，所以我就不再展示图片了
     
    
  
 

文章来源: blog.csdn.net，作者：此星光明2021年博客之星云计算Top3，版权归原作者所有，如需转载，请联系作者。

原文链接：blog.csdn.net/qq_31988139/article/details/119764306