【发布时间】：2019-12-09 01:32:08
【问题描述】：
我需要对存储为 xarray.Dataset 并包含 nan 值（出于技术原因给定像素的每个波段都是 nan）的一些 Landsat（卫星图像）场景应用 PCA 转换。
这是创建示例数据集的代码：
import numpy as np
import xarray as xr
# Create a demo xarray.Dataset
ncols = 25
nrows = 50
lon = [50 + x * 0.2 for x in range(nrows)]
lat = [30 + x * 0.2 for x in range(ncols)]
red = np.random.rand(nrows, ncols) * 10000
green = np.random.rand(nrows, ncols) * 10000
blue = np.random.rand(nrows, ncols) * 10000
nir = np.random.rand(nrows, ncols) * 10000
swir1 = np.random.rand(nrows, ncols) * 10000
swir2 = np.random.rand(nrows, ncols) * 10000
ds = xr.Dataset({'red': (['longitude', 'latitude'], red),
'green': (['longitude', 'latitude'], green),
'blue': (['longitude', 'latitude'], blue),
'nir': (['longitude', 'latitude'], nir),
'swir1': (['longitude', 'latitude'], swir1),
'swir2': (['longitude', 'latitude'], swir2)},
coords = {'longitude': (['longitude'], lon),
'latitude': (['latitude'], lat)})
# To keep example realistic let's add some nodata
ds = ds.where(ds.latitude + ds.longitude < 90)
print(ds)
<xarray.Dataset> Dimensions: (latitude: 25, longitude: 50) Coordinates: * longitude (longitude) float64 50.0 50.2 50.4 50.6
50.8 51.0 51.2 51.4 ... * latitude (latitude) float64 30.0 30.2 30.4 30.6 30.8 31.0 31.2 31.4 ... Data variables:
red (longitude, latitude) float64 6.07e+03 13.8 9.682e+03 ...
green (longitude, latitude) float64 5.476e+03 350.4 7.556e+03 ...
blue (longitude, latitude) float64 4.306e+03 2.104e+03 9.267e+03 ...
nir (longitude, latitude) float64 1.445e+03 8.633e+03 6.388e+03 ...
swir1 (longitude, latitude) float64 6.005e+03 7.692e+03 4.004e+03 ...
swir2 (longitude, latitude) float64 8.235e+03 3.127e+03 674.6 ...
在网上搜索了一下，尝试实现sklearn.decomposition PCA函数，没有成功。
我首先将每个二维带转换为一个维度：
# flatten dataset
tmp_list = []
for b in ['red', 'green', 'blue','nir','swir1','swir2']:
tmp_list.append(ds[b].values.flatten().astype('float64'))
flat_ds = np.array(tmp_list)
然后我尝试计算 PCA 并在没有 nan 的位置转换原始数据。我成功生成了一些输出，但与使用 ArcGIS 或 Grass 生成的输出完全不同。
当我更改我的位置时，似乎 sklearn 函数无法处理包含 nan 的数据。所以我从扁平数据集中删除了 nan 值，当我对扁平 PCA 结果进行缩减时，这是有问题的，因为它不包含原始数据集维度的倍数。
# deflate PCAs
dims = ds.dims['longitude'], ds.dims['latitude']
pcas = xr.Dataset()
for i in range(flat_pcas.shape[0]):
pcas['PCA_%i' % (i + 1)] = xr.DataArray(np.reshape(flat_pcas[i], dims),
coords=[ds.longitude.values, ds.latitude.values],
dims=['longitude','latitude'])
恢复情况：
是否存在另一种更简单的方法来在 xarray.Dataset 上实现 PCA 转换？
nan怎么处理？
【问题讨论】：
标签： python scikit-learn pca python-xarray