【GEE笔记】主成分分析（PCA）算法的实现和应用

news2025/1/9 16:53:47

前言

主成分分析（PCA）是一种常用的降维方法，它可以将多个相关的变量转换为少数几个不相关的变量，称为主成分（PC）。这些主成分可以反映原始变量的大部分信息，同时减少数据的复杂度和冗余性。在遥感领域，PCA可以用来提取影像的特征，消除噪声，增强对比度，或者进行分类和变化检测等。

本文介绍如何使用Google Earth Engine（GEE）平台实现PCA算法，并且展示一个应用案例，即利用PCA对哨兵二号（Sentinel-2）影像进行降维。

PCA算法原理

PCA算法的基本思想是通过正交变换，将原始数据投影到一个新的坐标系中，使得新坐标系的各个轴（即主成分）之间相互正交，且按照方差大小递减排序。这样，第一个主成分就是原始数据在最大方差方向上的投影，它包含了最多的信息；第二个主成分就是原始数据在次大方差方向上的投影，它包含了次多的信息；以此类推，直到所有的主成分都被提取出来。通常情况下，我们只需要保留前几个主成分，就可以达到较好的降维效果。

PCA算法的具体步骤如下：

  1. 对原始数据进行中心化处理，即每个变量减去其均值，使得新的数据集的均值为零。
  2. 计算原始数据的协方差矩阵，它反映了各个变量之间的线性相关程度。
  3. 对协方差矩阵进行特征值分解，得到特征值和特征向量。特征值表示了各个主成分的方差大小，特征向量表示了各个主成分的方向。
  4. 按照特征值的大小降序排列，选择前k个最大的特征值对应的特征向量作为新坐标系的基。
  5. 将原始数据乘以特征向量矩阵，得到新坐标系下的数据，即主成分。

GEE平台上的PCA实现

数据准备

定义了一个矩形区域作为研究区域

var geometry = 
    ee.Geometry.Polygon(
        [[[119.34885102549033, 36.519170262470254],
          [119.34885102549033, 36.32081057978827],
          [119.56583100595908, 36.32081057978827],
          [119.56583100595908, 36.519170262470254]]], null, false);

获取哨兵2影像，并对其进行一些预处理。选择2020年3月至5月期间，在哨兵2表面反射率（SR）影像集合中根据自定义的几何区域geometry进行裁剪。选择10个波段（B2, B3, B4, B5, B6, B7, B8, B8A, B11, B12），并对每个波段应用一个云掩膜函数，以去除云影响。然后均值合成为一幅图像。最后进行图像重采样，将空间分辨率转为统一的10米。以下是相关的代码：

var year = 2020
var bandlist = ['B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'B8', 'B8A', 'B11', 'B12']
var start = ee.Date(year + '-3-01');
var finish = ee.Date(year + '-5-1');
var dataset = ee.ImageCollection('COPERNICUS/S2_SR')
    .filterDate(start, finish)
    .filterBounds(geometry)
    .map(maskS2clouds)
var dataset10 = dataset.select(bandlist);
var bandsimage = dataset10.mean().clip(geometry)
print("输入影像", bandsimage)

bandsimage = bandsimage.reproject('EPSG:4326', null, 10);
var scaledataset10 = bandsimage.projection().nominalScale();
var rgbVis = {
    min: 0.0,
    max: 3000,
    bands: ['B4', 'B3', 'B2'],
};
Map.centerObject(geometry, 12)
Map.addLayer(bandsimage, rgbVis, 'bandsimage');

// 定义一个云掩膜函数
function maskS2clouds(image) {
    var qa = image.select('QA60');
    // Bits 10 and 11 are clouds and cirrus, respectively.
    var cloudBitMask = 1 << 10;
    var cirrusBitMask = 1 << 11;
    var mask = qa.bitwiseAnd(cloudBitMask).eq(0)
        .and(qa.bitwiseAnd(cirrusBitMask).eq(0));
    return image.updateMask(mask).clip(geometry);
}

运行上述代码后，我们可以在地图上看到输入影像的真彩色显示，如下图所示：
在这里插入图片描述

PCA变换

接下来，实现PCA算法并将其封装为一个函数 doPCAandCR 。该函数的输入参数是一个多波段影像，输出参数是一个包含两个属性的对象: pcImage 和 cntributionRate 。 pcImage 是一个多波段影像，每个波段对应一个主成分分量， contributionRate 是一个数组，每个元素对应一个主成分的贡献率

PCA算法的主要步骤如下:

  • 计算输入影像的每个波段的均值，并将其作为一个常数影像
  • 将输入影像减去均值影像，得到中心化后的影像
  • 将中心化后的影像转换为数组格式
  • 计算数组的协方差知阵
  • 对协方差知阵进行特征值分解，得到特征值和特征向量
  • 将特征向量知阵与数组相乘，得到主成分数组
  • 将主成分数组转换为影像格式，并按照波段顺序命名
  • 将每个波段的特征值除以所有特征值之和，得到每个波段的贡献率
var result = doPCAandCR(bandsimage);
print("PCA", result.pcImage);
Map.addLayer(result.pcImage.select(0), {}, 'resultPCA1');
print("每个波段的贡献率", result.contributionRate);

// 定义一个PCA和贡献率计算函数
function doPCAandCR(image) {
    // Get the band names of the image
    var bandNames = image.bandNames();
    var region = image.geometry();
    var meanDict = image.reduceRegion({
        reducer: ee.Reducer.mean(),
        geometry: geometry,
        scale: scaledataset10,
        tileScale: 16,
        maxPixels: 1e9
    });
    var means = ee.Image.constant(meanDict.values(bandNames));
    var centered = image.subtract(means);
    var arrays = centered.toArray();
    var covar = arrays.reduceRegion({
        reducer: ee.Reducer.centeredCovariance(),
        geometry: geometry,
        scale: scaledataset10,
        tileScale: 16,
        maxPixels: 1e9
    });
    var covarArray = ee.Array(covar.get('array'));
    var eigens = covarArray.eigen();
    var eigenValues = eigens.slice(1, 0, 1);
    var eigenVectors = eigens.slice(1, 1);
    var arrayImage = arrays.toArray(1);
    var principalComponents = ee.Image(eigenVectors).matrixMultiply(arrayImage);
    var sdImage = ee.Image(eigenValues.sqrt())
        .arrayProject([0])
        .arrayFlatten([getNewBandNames('sd', bandNames)]);
    var pcImage = principalComponents
        .arrayProject([0])
        .arrayFlatten([getNewBandNames('pc', bandNames)])
        .divide(sdImage);
    // 通过将每个波段的特征值除以所有特征值之和来计算每个波段的贡献率
    var sumEigenValues = eigenValues.reduce(ee.Reducer.sum(), [0]).get([0, 0]);
    var contributionRate = eigenValues.divide(sumEigenValues).project([0]);
    return { pcImage: pcImage, contributionRate: contributionRate };
}

function getNewBandNames(prefix, bandNames) {
    var seq = ee.List.sequence(1, bandNames.length());
    return seq.map(function (b) {
        return ee.String(prefix).cat(ee.Number(b).int());
    });
}

将准备的数据进行主成分变换，得到各主成分及其贡献率，并将第一主成分可视化
在这里插入图片描述

在这里插入图片描述

计算累积贡献率

// 对贡献率进行累加操作
var cumsum = ee.List(result.contributionRate.toList().iterate(accumulate, ee.List([])));
print("累积贡献率", cumsum);

// 定义一个累加函数
function accumulate(value, list) {
    // 将list转换为ee.List对象
    list = ee.List(list);
    // 获取上一次的累加结果，如果没有则用0代替
    var previous = ee.Algorithms.If(list.size(), ee.Number(list.get(-1)), 0)//ee.Number(list.get(-1)).or(0);
    // 将当前值与上一次的结果相加
    var added = ee.Number(value).add(previous);
    // 返回更新后的list
    return list.add(added);
};

在这里插入图片描述

本文来自互联网用户投稿，该文观点仅代表作者本人，不代表本站立场。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如若转载，请注明出处：http://www.coloradmin.cn/o/787534.html

如若内容造成侵权/违法违规/事实不符，请联系多彩编程网进行投诉反馈，一经查实，立即删除！

相关文章

结构型设计模式之组合模式【设计模式系列】

结构型设计模式之组合模式【设计模式系列】

系列文章目录 C技能系列 Linux通信架构系列 C高性能优化编程系列 深入理解软件架构设计系列 高级C并发线程编程 设计模式系列 期待你的关注哦&#xff01;&#xff01;&#xff01; 现在的一切都是为将来的梦想编织翅膀&#xff0c;让梦想在现实中展翅高飞。 Now everythi…
阅读更多...
【nginx】nginx中root与alias的区别：

【nginx】nginx中root与alias的区别：

文章目录 root与alias主要区别在于nginx如何解释location后面的uri&#xff0c;这会使两者分别以不同的方式将请求映射到服务器文件上。 root的处理结果是&#xff1a;root路径&#xff0b;location路径 alias的处理结果是&#xff1a;使用alias路径替换location路径 alias是一…
阅读更多...
前端配置nginx反向代理

前端配置nginx反向代理

1.下载nginx 2.进入nginx.conf 3.配置 4.将nginx放入前端项目根目录 5.启动后端项目和nginx服务 启动命令 start nginx ; 重新启动命令 nginx -s reload 6.访问localhost就能看到项目了&#xff08;默认80端口&#xff09; 注意&#xff1a;如果nginx配置了域名&#xff0…
阅读更多...
Android ConstraintLayout使用攻略

Android ConstraintLayout使用攻略

原文链接 Android ConstraintLayout使用攻略 ConstraintLayout是新一代的布局&#xff0c;它汲取了众家之长&#xff0c;把布局的概念进行了大统一&#xff0c;灵活且强大&#xff0c;基本上可以干掉以前所有的常用布局&#xff08;LinearLayout, RelativeLayout和FrameLayout…
阅读更多...
94.qt qml-分页Table表格组件

94.qt qml-分页Table表格组件

在我们之前学习了87.qt qml-分页组件控件(支持设置任意折叠页数等)_qt分页控件_诺谦的博客-CSDN博客 然后我们又学习了Table实现，所以本章实现一个分页Table表格组件，配合分页控件, 模拟请求服务器数据来实现数据分解效果，因为一般使用分页的时候，一般都是分页请求，避免数…
阅读更多...
Flutter的开发环境搭建-图解

Flutter的开发环境搭建-图解

前言&#xff1a;Flutter作为一个移动应用开发框架&#xff0c;具有许多优点和一些局限性。最大的优点就是-跨平台开发&#xff1a;Flutter可以在iOS和Android等多个平台上进行跨平台开发&#xff0c;使用一套代码编写应用程序&#xff0c;节省开发时间和成本。 Flutter可以编…
阅读更多...
python将大文件拆分为多个小文件

python将大文件拆分为多个小文件

如上图&#xff0c;目前采用单行不停写入的方式&#xff0c;这里是读了两次文件&#xff0c;第一次读取文件是为了获取总行数&#xff0c;第二次读取是取数据内容。 如果只读取一次文件&#xff0c;则会对内存有一定的要求&#xff0c;会需要在第一次读取数据的时候就将文件内…
阅读更多...
ONNX Runtime 加速深度学习(C++ 、python)详细介绍

ONNX Runtime 加速深度学习(C++ 、python)详细介绍

ONNX Runtime 加速深度学习(C 、python)详细介绍 本文在 https://blog.csdn.net/u013250861/article/details/127829944 基础上进行了更改&#xff0c;感谢原作&#xff01; ONNXRuntime(Open Neural Network Exchange)是微软推出的一款针对ONNX模型格式的推理框架&#xff0c…
阅读更多...
Redis常用数据类型和使用场景

Redis常用数据类型和使用场景

Redis目前支持5种数据类型&#xff0c;分别是&#xff1a; String&#xff08;字符串&#xff09; List&#xff08;列表&#xff09; Hash&#xff08;字典&#xff09; Set&#xff08;集合&#xff09; Sorted Set&#xff08;有序集合&#xff09; 下面就分别介绍这五…
阅读更多...
Qt Core学习日记——第六天QMetaMethod

Qt Core学习日记——第六天QMetaMethod

Qt子类会将每一个函数封装成QMetaMethod存储在对应的QMetaObject中&#xff0c;包括信号、槽函数、普通函数、构造函数、析构函数 函数解析 QMetaMethod::methodSignature 获取方法的签名 比如函数slot2&#xff0c;对应签名是“slot2(int*)” QMetaMethod::name 方法名称。…
阅读更多...
13.2.3 【Linux】新增与移除群组

13.2.3 【Linux】新增与移除群组

基本上&#xff0c;群组的内容都与这两个文件有关&#xff1a;/etc/group, /etc/gshadow。 群组的内容其实很简单&#xff0c;都是上面两个文件的新增、修改与移除而已。 groupadd 为了让使用者的 UID/GID 成对&#xff0c;建议新建的与使用者私有群组无关的其他群组时&#x…
阅读更多...
RabbitMQ入门，springboot整合RabbitMQ

RabbitMQ入门，springboot整合RabbitMQ

周末的两天没有写文章&#xff0c;因为项目分离出来了一个权限管理平台&#xff0c;花了一点时间整理项目&#xff0c;同时完成了一些功能的开发。 今天这篇文章介绍一下RabbitMQ这个消息中间件&#xff0c;以及通过springboot整合RabbitMQ。 目录 一、初步了解RabbitMQ 二、…
阅读更多...
学Java有哪些就业方向？

学Java有哪些就业方向？

俗话说&#xff1a;男怕入错行&#xff0c;女怕嫁错郎。众所周知&#xff0c;选工作就是选行业&#xff0c;行业和方向选对了&#xff0c;个人的发展就会随着行业风向青云直上&#xff0c;比同龄人更快的积累到财富。那究竟未来什么会是热门行业呢?这个真的很难预测&#xff0…
阅读更多...
【1++的C++初阶】之模板（二）

【1++的C++初阶】之模板（二）

&#x1f44d;作者主页&#xff1a;进击的1 &#x1f929; 专栏链接&#xff1a;【1的C初阶】 文章目录 一&#xff0c;非类型模板参数二&#xff0c;模板特化三&#xff0c;模板分离编译 一&#xff0c;非类型模板参数 模板参数分为类类型模板参数与非类型模板参数。 类类型形…
阅读更多...
【雕爷学编程】Arduino动手做（167）---MG996R金属齿轮舵机2

【雕爷学编程】Arduino动手做（167）---MG996R金属齿轮舵机2

37款传感器与模块的提法&#xff0c;在网络上广泛流传&#xff0c;其实Arduino能够兼容的传感器模块肯定是不止37种的。鉴于本人手头积累了一些传感器和执行器模块&#xff0c;依照实践出真知&#xff08;一定要动手做&#xff09;的理念&#xff0c;以学习和交流为目的&#x…
阅读更多...
苹果“空间音频导航”专利曝光，提供导航指引，跟声音走就对啦?

苹果“空间音频导航”专利曝光，提供导航指引，跟声音走就对啦?

近日&#xff0c;苹果公司成功申请一项专利&#xff0c;该专利名为“空间音频导航”&#xff0c;该专利详细说明了如何利用双耳音频设备&#xff08;AirPods或Apple Vision Pro&#xff09;为用户提供导航指引。 “空间音频导航”是一种模拟声音来源方向和距离的技术&#xff0…
阅读更多...
STM32MP157驱动开发——按键驱动（POLL 机制）

STM32MP157驱动开发——按键驱动（POLL 机制）

文章目录 “POLL ”机制&#xff1a;APP执行过程驱动使用的函数应用使用的函数pollfd结构体poll函数事件类型实现原理 poll方式的按键驱动程序(stm32mp157)gpio_key_drv.cbutton_test.cMakefile修改设备树文件编译测试 “POLL ”机制&#xff1a; 使用休眠-唤醒的方式等待某个…
阅读更多...
c# Outlook检索设定问题

c# Outlook检索设定问题

基于c# 设定outlook约会予定&#xff0c;时间格式是YYYY-MM-DD HH:mm 的情报。 问题发生&#xff1a; 根据开始时间&#xff08;2023/01/01 7:00&#xff09;条件查询该时间是否存在outlook信息时&#xff0c;明明存在一条数据&#xff0c;就是查询不出来数据 c#代码 Strin…
阅读更多...
单源最短路的扩展应用

单源最短路的扩展应用

AcWing 1137. 选择最佳线路 有一天&#xff0c;琪琪想乘坐公交车去拜访她的一位朋友。 由于琪琪非常容易晕车&#xff0c;所以她想尽快到达朋友家。 现在给定你一张城市交通路线图&#xff0c;上面包含城市的公交站台以及公交线路的具体分布。 已知城市中共包含 n 个车站…
阅读更多...
解决 Visual Studio Code 编译器代码自动格式化

解决 Visual Studio Code 编译器代码自动格式化

首先找到.vscode下的settings.json配置文件 将vue3snippets.enable-compile-vue-file-on-did-save-code更改为false
阅读更多...
推荐文章
最新文章

Copyright @ 2022~2023