Geocomputation （2）属性数据操作

已于 2023-12-08 16:20:11 修改
阅读量1.7k 收藏
点赞数
分类专栏： GIS空间分析 文章标签： python 开发语言
于 2023-06-21 01:10:13 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/weixin_44026026/article/details/131318689
版权

Geocomputation （2）属性数据操作

来源：https://github.com/geocompx/geocompy

1.准备

#| echo: false
import pandas as pd
import matplotlib.pyplot as plt
pd.options.display.max_rows = 6
pd.options.display.max_columns = 6
pd.options.display.max_colwidth = 35
plt.rcParams['figure.figsize'] = (5, 5)

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import geopandas as gpd
import rasterio

示例数据

#| echo: false
import os
from urllib.request import urlretrieve

data_path = 'data'
if not os.path.exists(data_path):
  os.mkdir(data_path)
  
file_path = 'data/landsat.tif'
url = 'https://github.com/geocompx/geocompy/releases/download/0.1/landsat.tif'
if not os.path.exists(file_path):
  print('Attempting to get the data')
  urlretrieve(url, file_path)

Attempting to get the data

world = gpd.read_file('data/world.gpkg')
src_elev = rasterio.open('data/elev.tif')
src_multi_rast = rasterio.open('data/landsat.tif')

2. Introduction

属性数据是与地理（几何）数据相关联的非空间信息。
公共汽车站提供了一个简单的示例：除了名称之外，它的位置通常由纬度和经度坐标（几何数据）表示。
例如，伦敦的 Elephant & Castle / New Kent Road 停靠点的坐标为经度“-0.098”和纬度“51.495”，在@中描述的简单特征表示中可以表示为“POINT (-0.098 51.495)” sec-空间类。
诸如“POINT”特征（使用简单特征术语）的名称属性等属性是本章的主题。

另一个例子是栅格数据中特定网格单元的高程值（属性）。 与矢量数据模型不同，栅格数据模型间接存储网格单元的坐标，这意味着属性和空间信息之间的区别不太明确。
为了说明这一点，请考虑光栅矩阵第 3 行和第 4 列中的像素。
它的空间位置由其在矩阵中的索引定义：从原点在 x 方向（通常在地图上向东和向右）移动四个像元，在 y 方向（通常向南和向下）移动三个像元。
光栅的分辨率定义了标题中指定的每个 x 和 y 步长的距离。
标头是栅格数据集的重要组成部分，它指定像素如何与地理坐标相关（另请参阅@spatial-operations 章）。

本章介绍如何根据属性（例如矢量数据集中的公交车站名称和栅格数据集中的像素高程）来操作地理对象。
对于矢量数据，这意味着子集化和聚合等技术（参见@sec-vector-attribute-subsetting 和@sec-vector-attribute-aggregation）。
@sec-vector-attribute-joining 和@sec-creating-attributes-and-removing-spatial-information 分别演示了如何使用共享 ID 将数据连接到简单要素对象以及如何创建新变量。
这些操作中的每一个都有一个空间等价物：例如，用于使用布尔值Series 对(Geo)DataFrame 进行子集化的[ 运算符适用于基于对象的属性和空间关系的子集化对象，这些方法使用诸如 .相交； 您还可以使用空间连接连接两个地理数据集中的属性。 这是个好消息：本章培养的技能可以交叉迁移。
@sec-spatial-operations 将此处介绍的方法扩展到空间世界。

在下一节深入探讨各种类型的矢量属性操作之后，@sec-manipulating-raster-objects 介绍了栅格属性数据操作，它演示了如何创建包含连续和分类属性的栅格图层并从其中提取像元值 或更多层（栅格子集）。
@sec-summarizing-raster-objects 提供了“全局”栅格操作的概述，可用于汇总整个栅格数据集。

3. 矢量属性操作

如@sec-vector-layers 中所述，矢量图层（GeoDataFrame，来自包geopandas）基本上是扩展表（来自包pandas 的DataFrame），不同之处在于矢量图层具有 几何列。
由于 GeoDataFrame 扩展了 DataFrame，矢量图层也支持包 pandas 中所有与表相关的普通操作，如下所示。

3.1 向量属性子集 {#sec-vector-attribute-subsetting}

pandas 支持多个子集接口，但最推荐 是：

  • .loc，它使用 pandas 索引，以及
  • .iloc，它使用（隐式）numpy 风格的数字索引。

在这两种情况下，该方法后跟方括号和两个索引，以逗号分隔。 每个索引可以包括：

  • 一个特定的值，如 1
  • 一片，如 0:3
  • 一个 list，如 [0,2,4]
  • :——表示“所有”索引

此规则的一个例外是使用列表选择列，如“df[[‘a’,‘b’]]”，而不是“df.loc[:, [‘a’,‘b’]]”， 从 df 中选择列 'a' 和 'b'。

以下是对世界各国的“GeoDataFrame”进行子集化的几个示例。

按位置子集行，例如前三行：

world.iloc[0:3, :]
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
0FJFijiOceania...69.9608222.253784MULTIPOLYGON (((-180.00000 -16....
1TZTanzaniaAfrica...64.1632402.099404MULTIPOLYGON (((33.90371 -0.950...
2EHWestern SaharaAfrica...NaNNaNMULTIPOLYGON (((-8.66559 27.656...

3 rows × 11 columns

这相当于：

world.iloc[:3]
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
0FJFijiOceania...69.9608222.253784MULTIPOLYGON (((-180.00000 -16....
1TZTanzaniaAfrica...64.1632402.099404MULTIPOLYGON (((33.90371 -0.950...
2EHWestern SaharaAfrica...NaNNaNMULTIPOLYGON (((-8.66559 27.656...

3 rows × 11 columns

world.head(3)
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
0FJFijiOceania...69.9608222.253784MULTIPOLYGON (((-180.00000 -16....
1TZTanzaniaAfrica...64.1632402.099404MULTIPOLYGON (((33.90371 -0.950...
2EHWestern SaharaAfrica...NaNNaNMULTIPOLYGON (((-8.66559 27.656...

3 rows × 11 columns

world.iloc[:, 0:3]
iso_a2name_longcontinent
0FJFijiOceania
1TZTanzaniaAfrica
2EHWestern SaharaAfrica
............
174XKKosovoEurope
175TTTrinidad and TobagoNorth America
176SSSouth SudanAfrica

177 rows × 3 columns

按位置子集行和列：

world.iloc[0:3, 0:3]
iso_a2name_longcontinent
0FJFijiOceania
1TZTanzaniaAfrica
2EHWestern SaharaAfrica

按名称子集列：

world[['name_long', 'geometry']]
name_longgeometry
0FijiMULTIPOLYGON (((-180.00000 -16....
1TanzaniaMULTIPOLYGON (((33.90371 -0.950...
2Western SaharaMULTIPOLYGON (((-8.66559 27.656...
.........
174KosovoMULTIPOLYGON (((20.59025 41.855...
175Trinidad and TobagoMULTIPOLYGON (((-61.68000 10.76...
176South SudanMULTIPOLYGON (((30.83385 3.5091...

177 rows × 2 columns

给定列之间的“切片”列：

world.loc[:, 'name_long':'pop']
name_longcontinentregion_un...typearea_km2pop
0FijiOceaniaOceania...Sovereign country19289.970733885806.0
1TanzaniaAfricaAfrica...Sovereign country932745.79235752234869.0
2Western SaharaAfricaAfrica...Indeterminate96270.601041NaN
........................
174KosovoEuropeEurope...Sovereign country11230.2616721821800.0
175Trinidad and TobagoNorth AmericaAmericas...Sovereign country7737.8098551354493.0
176South SudanAfricaAfrica...Sovereign country624909.09908611530971.0

177 rows × 7 columns

通过布尔值列表（0 和 1 或 True 和 False）进行子集化：

x = [1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0] 
world.iloc[:, x]
name_longname_longiso_a2...name_longiso_a2iso_a2
0FijiFijiFJ...FijiFJFJ
1TanzaniaTanzaniaTZ...TanzaniaTZTZ
2Western SaharaWestern SaharaEH...Western SaharaEHEH
........................
174KosovoKosovoXK...KosovoXKXK
175Trinidad and TobagoTrinidad and TobagoTT...Trinidad and TobagoTTTT
176South SudanSouth SudanSS...South SudanSSSS

177 rows × 11 columns

我们可以使用 .drop 方法通过 id 删除特定行，例如，删除第 2、3 和 5 行：

world.drop([2, 3, 5])
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
0FJFijiOceania...69.9600008222.253784MULTIPOLYGON (((-180.00000 -16....
1TZTanzaniaAfrica...64.1630002402.099404MULTIPOLYGON (((33.90371 -0.950...
4USUnited StatesNorth America...78.84146351921.984639MULTIPOLYGON (((-171.73166 63.7...
........................
174XKKosovoEurope...71.0975618698.291559MULTIPOLYGON (((20.59025 41.855...
175TTTrinidad and TobagoNorth America...70.42600031181.821196MULTIPOLYGON (((-61.68000 10.76...
176SSSouth SudanAfrica...55.8170001935.879400MULTIPOLYGON (((30.83385 3.5091...

174 rows × 11 columns

或者使用 .drop 方法和 axis=1 （即列）删除特定列：

world.drop(['name_long', 'continent'], axis=1)
iso_a2region_unsubregion...lifeExpgdpPercapgeometry
0FJOceaniaMelanesia...69.9600008222.253784MULTIPOLYGON (((-180.00000 -16....
1TZAfricaEastern Africa...64.1630002402.099404MULTIPOLYGON (((33.90371 -0.950...
2EHAfricaNorthern Africa...NaNNaNMULTIPOLYGON (((-8.66559 27.656...
........................
174XKEuropeSouthern Europe...71.0975618698.291559MULTIPOLYGON (((20.59025 41.855...
175TTAmericasCaribbean...70.42600031181.821196MULTIPOLYGON (((-61.68000 10.76...
176SSAfricaEastern Africa...55.8170001935.879400MULTIPOLYGON (((30.83385 3.5091...

177 rows × 9 columns

我们可以使用 .rename 方法重命名列：

world[['name_long', 'pop']].rename(columns={'pop': 'population'})
name_longpopulation
0Fiji885806.0
1Tanzania52234869.0
2Western SaharaNaN
.........
174Kosovo1821800.0
175Trinidad and Tobago1354493.0
176South Sudan11530971.0

177 rows × 2 columns

标准的 numpy 比较运算符可用于布尔子集化，如表 @tbl-comparison-operators 所示。

符号名字
==等于
!=不等于
>, <大于/小于
>=, <=大于/小于或等于
&, `, ~`

: 返回布尔值的比较运算符 (True/False)。 {#tbl-比较运算符}

以下示例通过创建一个名为“small_countries”的新“GeoDataFrame”对象演示了用于子集化的逻辑向量，该对象仅包含表面积小于 10,000 km2 的国家：

idx_small = world['area_km2'] < 10000  ## a logical 'Series'
small_countries = world[idx_small]
small_countries
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
45PRPuerto RicoNorth America...79.39012235066.046376MULTIPOLYGON (((-66.28243 18.51...
79PSPalestineAsia...73.1260004319.528283MULTIPOLYGON (((35.39756 31.489...
89VUVanuatuOceania...71.7090002892.341604MULTIPOLYGON (((166.79316 -15.6...
........................
160NoneNorthern CyprusAsia...NaNNaNMULTIPOLYGON (((32.73178 35.140...
161CYCyprusAsia...80.17300029786.365653MULTIPOLYGON (((32.73178 35.140...
175TTTrinidad and TobagoNorth America...70.42600031181.821196MULTIPOLYGON (((-61.68000 10.76...

7 rows × 11 columns

中介“idx_small”（代表小国家的索引的缩写）是一个布尔值“Series”，可用于按表面积划分世界上七个最小国家的子集。
省略中间对象的更简洁的命令会生成相同的结果：

small_countries = world[world['area_km2'] < 10000]
small_countries
iso_a2name_longcontinent...lifeExpgdpPercapgeometry
45PRPuerto RicoNorth America...79.39012235066.046376MULTIPOLYGON (((-66.28243 18.51...
79PSPalestineAsia...73.1260004319.528283MULTIPOLYGON (((35.39756 31.489...
89VUVanuatuOceania...71.7090002892.341604MULTIPOLYGON (((166.79316 -15.6...
........................
160NoneNorthern CyprusAsia...NaNNaNMULTIPOLYGON (((32.73178 35.140...
161CYCyprusAsia...80.17300029786.365653MULTIPOLYGON (((32.73178 35.140...
175TTTrinidad and TobagoNorth America...70.42600031181.821196MULTIPOLYGON (((-61.68000 10.76...

7 rows × 11 columns

上面显示的各种方法可以链接起来，用于与多个子集步骤的任意组合，例如：

world[world['continent'] == 'Asia']  \
    .loc[:, ['name_long', 'continent']]  \
    .iloc[0:5, :]
name_longcontinent
5KazakhstanAsia
6UzbekistanAsia
8IndonesiaAsia
24Timor-LesteAsia
76IsraelAsia

我们还可以组合索引：

idx_small = world['area_km2'] < 10000
idx_asia = world['continent'] == 'Asia'
world.loc[idx_small & idx_asia, ['name_long', 'continent', 'area_km2']]
name_longcontinentarea_km2
79PalestineAsia5037.103826
160Northern CyprusAsia3786.364506
161CyprusAsia6207.006191

3.2 矢量属性聚合

聚合涉及基于一个或多个分组变量（通常是列中的值；地理聚合将在下一章介绍）汇总数据。 这种基于属性的聚合的一个典型示例是根据国家级数据（每个国家一行）计算每个大陆的人数。
“world”数据集包含必要的成分：列“pop”和“continent”，分别是人口和分组变量。 目的是找到每个大陆的国家人口的“sum()”，从而得到更小的数据框。 （由于聚合是数据缩减的一种形式，因此在处理大型数据集时它可能是一个有用的早期步骤）。 这种聚合可以使用 .groupby 和 .sum 的组合来实现：

world_agg1 = world[['continent', 'pop']].groupby('continent').sum()
world_agg1
pop
continent
Africa1.154947e+09
Antarctica0.000000e+00
Asia4.311408e+09
......
Oceania3.775783e+07
Seven seas (open ocean)0.000000e+00
South America4.120608e+08

8 rows × 1 columns

如果您不喜欢默认使用的科学记数法来显示人口总和，您可以像这样更改浮点值的 Pandas 显示格式：

pd.set_option('display.float_format', '{:.0f}'.format)
world_agg1
pop
continent
Africa1154946633
Antarctica0
Asia4311408059
......
Oceania37757833
Seven seas (open ocean)0
South America412060811

8 rows × 1 columns

结果是一个（非空间）表，有八行，每个大陆一行，两列报告每个大陆的名称和人口。

如果我们想在聚合结果中包含几何图形，我们可以使用 .dissolve 方法。
这样，除了人口总和之外，我们还可以获得每个大陆的相关几何图形，即所有国家的联合。
请注意，我们使用“by”参数来选择用于分组的列，并使用“aggfunc”参数来选择非几何列的聚合函数：

world_agg2 = world[['continent', 'pop', 'geometry']] \
    .dissolve(by='continent', aggfunc='sum') \
    .reset_index()
world_agg2
continentgeometrypop
0AfricaMULTIPOLYGON (((32.83012 -26.74...1154946633
1AntarcticaMULTIPOLYGON (((-66.29003 -80.2...0
2AsiaMULTIPOLYGON (((120.29501 -10.2...4311408059
............
5OceaniaMULTIPOLYGON (((169.66781 -43.5...37757833
6Seven seas (open ocean)POLYGON ((68.93500 -48.62500, 6...0
7South AmericaMULTIPOLYGON (((-66.95992 -54.8...412060811

8 rows × 3 columns

@fig-spatial-aggregation 显示结果：

#| label: fig-spatial-aggregation
#| fig-cap: Continents with summed population

world_agg2.plot(column='pop', legend=True);

在这里插入图片描述

生成的“world_agg2”对象是一个“GeoDataFrame”，其中包含代表世界各大洲（和公海）的 8 个特征。

.dissolve 中 aggfunc 参数的其他选项 包括：

*'第一'
*'最后'
*'分钟'
*'最大'
*'总和'
*'平均'
*'中位数'

此外，我们可以传递自定义函数。

作为一个更复杂的例子，下面是我们如何计算每个大陆的总人口、面积和国家数量：

world_agg3 = world.dissolve(
    by='continent', aggfunc={
         "name_long": "count",
         "pop": "sum",
         'area_km2': "sum"
     }).rename(columns={'name_long': 'n'})
world_agg3
geometrynpoparea_km2
continent
AfricaMULTIPOLYGON (((32.83012 -26.74...51115494663329946198
AntarcticaMULTIPOLYGON (((-66.29003 -80.2...1012335956
AsiaMULTIPOLYGON (((120.29501 -10.2...47431140805931252459
...............
OceaniaMULTIPOLYGON (((169.66781 -43.5...7377578338504489
Seven seas (open ocean)POLYGON ((68.93500 -48.62500, 6...1011603
South AmericaMULTIPOLYGON (((-66.95992 -54.8...1341206081117762592

8 rows × 4 columns

图 @fig-spatial-aggregation-different-functions 可视化了具有三个聚合属性的大陆的结果层 (world_agg3)。

#| label: fig-spatial-aggregation-different-functions
#| fig-cap: 'Continent properties, calculated using spatial aggregation using different functions'

fig, axes = plt.subplots(2, 2, figsize=(9, 5))
world_agg3.plot(column='pop', edgecolor='black', legend=True, ax=axes[0][0])
world_agg3.plot(column='area_km2', edgecolor='black', legend=True, ax=axes[0][1])
world_agg3.plot(column='n', edgecolor='black', legend=True, ax=axes[1][0])
axes[0][0].set_title('Summed population')
axes[0][1].set_title('Summed area')
axes[1][0].set_title('Count of countries')
fig.delaxes(axes[1][1]);

在这里插入图片描述

下面继续最后的结果来演示其他表相关的操作。 鉴于“world_agg3”大陆摘要（@fig-spatial-aggregation-different-functions），我们：

  • 删除几何列，
    *计算每个大陆的人口密度，
  • 按大洲包含的国家数量排列大洲，以及
  • 只保留 3 个人口最多的大陆。
world_agg4 = world_agg3.drop(columns=['geometry'])
world_agg4['density'] = world_agg4['pop'] / world_agg4['area_km2']
world_agg4 = world_agg4.sort_values(by='n', ascending=False)
world_agg4 = world_agg4.head(3)
world_agg4
npoparea_km2density
continent
Africa5111549466332994619839
Asia47431140805931252459138
Europe396690362562306521929

3.3 矢量属性加入

合并来自不同来源的数据是数据准备中的一项常见任务。 联接通过基于共享的“键”变量组合表来实现这一点。
pandas 有一个名为 pd.merge 的函数，用于基于公共列加入 (Geo)DataFrames。
pd.merge 函数遵循数据库语言 SQL 中使用的约定（Grolemund 和 Wickham 2016）。
pd.merge 函数对 DataFrame 和 GeoDataFrame 对象的作用相同。
pd.merge 的结果可以是 DataFrame 或 GeoDataFrame 对象，具体取决于输入。

一种常见的空间数据属性连接类型是将“DataFrames”连接到“GeoDataFrames”。
为实现这一点，我们使用带有 GeoDataFrame 作为第一个参数的 pd.merge，并从指定为第二个参数的 DataFrame 向其添加列。
在以下示例中，我们将咖啡生产数据与“世界”数据集结合起来。
咖啡数据位于从主要咖啡生产国的 CSV 文件中导入的名为“coffee_data”的“DataFrame”中：

coffee_data = pd.read_csv('data/coffee_data.csv')
coffee_data
name_longcoffee_production_2016coffee_production_2017
0AngolaNaNNaN
1Bolivia34
2Brazil32772786
............
44Zambia3NaN
45Zimbabwe11
46Others2326

47 rows × 3 columns

它的三列是：

  • name_long 国家名称
  • coffee_production_2016 和 coffee_production_2017 包含以每年 60 公斤袋为单位的咖啡产量估计值。

保留第一个数据集的左连接基于常见的“name_long”列将“world”与“coffee_data”合并：

world_coffee = pd.merge(world, coffee_data, on='name_long', how='left')
world_coffee
iso_a2name_longcontinent...geometrycoffee_production_2016coffee_production_2017
0FJFijiOceania...MULTIPOLYGON (((-180.00000 -16....NaNNaN
1TZTanzaniaAfrica...MULTIPOLYGON (((33.90371 -0.950...8166
2EHWestern SaharaAfrica...MULTIPOLYGON (((-8.66559 27.656...NaNNaN
........................
174XKKosovoEurope...MULTIPOLYGON (((20.59025 41.855...NaNNaN
175TTTrinidad and TobagoNorth America...MULTIPOLYGON (((-61.68000 10.76...NaNNaN
176SSSouth SudanAfrica...MULTIPOLYGON (((30.83385 3.5091...NaNNaN

177 rows × 13 columns

结果是一个与原始 world 对象相同的 GeoDataFrame 对象，但有两个关于咖啡生产的新变量（coffee_production_2016 和 coffee_production_2017）。
这可以绘制成地图，如@fig-join-coffee-production 所示：

#| label: fig-join-coffee-production
#| fig-cap: 'World coffee production, thousand 60-kg bags by country, in 2017 (source: International Coffee Organization).'

base = world_coffee.plot(color='white', edgecolor='lightgrey')
coffee_map = world_coffee.plot(ax=base, column='coffee_production_2017')
coffee_map.set_title('Coffee production');

在这里插入图片描述

要工作，基于属性的连接需要在两个数据集中都有一个“关键变量”（“pd.merge”的“on”参数）。
在上面的示例中，world_coffee 和 world DataFrames 都包含一个名为 name_long 的列。
（默认情况下，pd.merge 使用所有具有匹配名称的列。但是，建议明确指定要用于匹配的列的名称，就像我们在上一个示例中所做的那样。）

如果列名不同，您可以使用 left_on 和 right_on 来指定各自的列。

请注意，结果“world_coffee”与原始数据集“world”具有相同的行数。
虽然 coffee_data 中只有 47 行，但 world_coffee 中所有 177 个国家/地区的记录都完好无损。 原始数据集中没有匹配项的行被分配了新咖啡生产变量的“np.nan”值。
这是左连接的一个特征（用 how='left' 指定），也是我们通常想要做的。

如果我们只想保留关键变量匹配的国家怎么办？
在这种情况下，可以使用内部连接：

pd.merge(world, coffee_data, on='name_long', how='inner')
iso_a2name_longcontinent...geometrycoffee_production_2016coffee_production_2017
0TZTanzaniaAfrica...MULTIPOLYGON (((33.90371 -0.950...8166
1PGPapua New GuineaOceania...MULTIPOLYGON (((141.00021 -2.60...11474
2IDIndonesiaAsia...MULTIPOLYGON (((104.36999 -1.08...742360
........................
42ETEthiopiaAfrica...MULTIPOLYGON (((47.78942 8.0030...215283
43UGUgandaAfrica...MULTIPOLYGON (((33.90371 -0.950...408443
44RWRwandaAfrica...MULTIPOLYGON (((30.41910 -1.134...3642

45 rows × 13 columns

连接两个 (Geo)DataFrame 的另一种方法是恰当地称为“join”函数：

world.join(coffee_data.set_index('name_long'), on='name_long', how='inner')
iso_a2name_longcontinent...geometrycoffee_production_2016coffee_production_2017
1TZTanzaniaAfrica...MULTIPOLYGON (((33.90371 -0.950...8166
7PGPapua New GuineaOceania...MULTIPOLYGON (((141.00021 -2.60...11474
8IDIndonesiaAsia...MULTIPOLYGON (((104.36999 -1.08...742360
........................
165ETEthiopiaAfrica...MULTIPOLYGON (((47.78942 8.0030...215283
168UGUgandaAfrica...MULTIPOLYGON (((33.90371 -0.950...408443
169RWRwandaAfrica...MULTIPOLYGON (((30.41910 -1.134...3642

45 rows × 13 columns

请注意，在这种情况下，我们需要将 coffee_data 的索引设置为 name_long 值以避免错误消息。

3.4 创建属性和删除空间信息

通常，我们希望在现有列的基础上创建一个新列。 例如，我们要计算每个国家的人口密度。 为此，我们需要将人口列（此处为“pop”）除以面积列（此处为“area_km2”）。 请注意，我们正在处理名为“world2”的“world”副本，因此我们不会修改原始图层：

world2 = world.copy()
world2['pop_dens'] = world2['pop'] / world2['area_km2']
world2
iso_a2name_longcontinent...gdpPercapgeometrypop_dens
0FJFijiOceania...8222MULTIPOLYGON (((-180.00000 -16....46
1TZTanzaniaAfrica...2402MULTIPOLYGON (((33.90371 -0.950...56
2EHWestern SaharaAfrica...NaNMULTIPOLYGON (((-8.66559 27.656...NaN
........................
174XKKosovoEurope...8698MULTIPOLYGON (((20.59025 41.855...162
175TTTrinidad and TobagoNorth America...31182MULTIPOLYGON (((-61.68000 10.76...175
176SSSouth SudanAfrica...1936MULTIPOLYGON (((30.83385 3.5091...18

177 rows × 12 columns

要将现有列粘贴（即连接）在一起，我们可以使用普通的 Python 字符串运算符“+”，就好像我们在处理单个字符串而不是“系列”一样。 例如，我们希望将“continent”和“region_un”列合并到名为“con_reg”的新列中，使用“:”作为分隔符。 随后，我们使用 .drop 删除原始列：

world2['con_reg'] = world['continent'] + ':' + world2['region_un']
world2 = world2.drop(['continent', 'region_un'], axis=1)
world2
iso_a2name_longsubregion...geometrypop_denscon_reg
0FJFijiMelanesia...MULTIPOLYGON (((-180.00000 -16....46Oceania:Oceania
1TZTanzaniaEastern Africa...MULTIPOLYGON (((33.90371 -0.950...56Africa:Africa
2EHWestern SaharaNorthern Africa...MULTIPOLYGON (((-8.66559 27.656...NaNAfrica:Africa
........................
174XKKosovoSouthern Europe...MULTIPOLYGON (((20.59025 41.855...162Europe:Europe
175TTTrinidad and TobagoCaribbean...MULTIPOLYGON (((-61.68000 10.76...175North America:Americas
176SSSouth SudanEastern Africa...MULTIPOLYGON (((30.83385 3.5091...18Africa:Africa

177 rows × 11 columns

生成的 sf 对象有一个名为“con_reg”的新列，代表每个国家/地区的大陆和地区，例如，“South America:Americas”代表阿根廷和其他南美国家。 相反的操作，根据分隔符字符串将一列拆分为多列，使用 .str.split 方法完成。 结果，我们回到了两个单独的“continent”和“region_un”列的先前状态（只是它们的位置现在是最后一个，因为它们是新创建的）：

world2[['continent', 'region_un']] = world2['con_reg'] \
    .str.split(':', expand=True)
world2
iso_a2name_longsubregion...con_regcontinentregion_un
0FJFijiMelanesia...Oceania:OceaniaOceaniaOceania
1TZTanzaniaEastern Africa...Africa:AfricaAfricaAfrica
2EHWestern SaharaNorthern Africa...Africa:AfricaAfricaAfrica
........................
174XKKosovoSouthern Europe...Europe:EuropeEuropeEurope
175TTTrinidad and TobagoCaribbean...North America:AmericasNorth AmericaAmericas
176SSSouth SudanEastern Africa...Africa:AfricaAfricaAfrica

177 rows × 13 columns

重命名一个或多个列可以使用 .rename 方法结合 columns 参数来完成，它应该是 old_name:new_name 形式的字典。 例如，以下命令将冗长的“name_long”列重命名为简单的“name”：

world2.rename(columns={'name_long': 'name'})
iso_a2namesubregion...con_regcontinentregion_un
0FJFijiMelanesia...Oceania:OceaniaOceaniaOceania
1TZTanzaniaEastern Africa...Africa:AfricaAfricaAfrica
2EHWestern SaharaNorthern Africa...Africa:AfricaAfricaAfrica
........................
174XKKosovoSouthern Europe...Europe:EuropeEuropeEurope
175TTTrinidad and TobagoCaribbean...North America:AmericasNorth AmericaAmericas
176SSSouth SudanEastern Africa...Africa:AfricaAfricaAfrica

177 rows × 13 columns

要一次更改所有列名，我们将“新”列名的“列表”分配给“.columns”属性。 list 的长度必须与列数相同（即 world.shape[1]）。 如下图所示，它输出相同的“world2”对象，但名称非常短：

new_names = ['a', 'b', 'c', 'd', 'e', 'f', 'g', 'h', 'geom', 'i', 'j', 'k', 'l']
world2.columns = new_names
world2
abc...jkl
0FJFijiMelanesia...Oceania:OceaniaOceaniaOceania
1TZTanzaniaEastern Africa...Africa:AfricaAfricaAfrica
2EHWestern SaharaNorthern Africa...Africa:AfricaAfricaAfrica
........................
174XKKosovoSouthern Europe...Europe:EuropeEuropeEurope
175TTTrinidad and TobagoCaribbean...North America:AmericasNorth AmericaAmericas
176SSSouth SudanEastern Africa...Africa:AfricaAfricaAfrica

177 rows × 13 columns

要对列重新排序，我们可以将修改后的列列表传递给子集运算符“[”。 例如，以下表达式按相反的字母顺序对 world2 列进行重新排序：

names = sorted(world2.columns, reverse=True)
world2 = world2[names]
world2
lkj...cba
0OceaniaOceaniaOceania:Oceania...MelanesiaFijiFJ
1AfricaAfricaAfrica:Africa...Eastern AfricaTanzaniaTZ
2AfricaAfricaAfrica:Africa...Northern AfricaWestern SaharaEH
........................
174EuropeEuropeEurope:Europe...Southern EuropeKosovoXK
175AmericasNorth AmericaNorth America:Americas...CaribbeanTrinidad and TobagoTT
176AfricaAfricaAfrica:Africa...Eastern AfricaSouth SudanSS

177 rows × 13 columns

这些属性数据操作中的每一个，即使它们在“pandas”包中定义并且适用于任何“DataFrame”，也会保留几何列和“GeoDataFrame”类。 然而，有时移除几何体是有意义的，例如为了加速聚合或仅导出属性数据以进行统计分析。 要从“GeoDataFrame”到“DataFrame”，我们需要：

1.删除几何列

  1. 从 GeoDataFrame 转换为 DataFrame

例如：

world2 = world2.drop('geom', axis=1)
world2 = pd.DataFrame(world2)
world2
lkj...cba
0OceaniaOceaniaOceania:Oceania...MelanesiaFijiFJ
1AfricaAfricaAfrica:Africa...Eastern AfricaTanzaniaTZ
2AfricaAfricaAfrica:Africa...Northern AfricaWestern SaharaEH
........................
174EuropeEuropeEurope:Europe...Southern EuropeKosovoXK
175AmericasNorth AmericaNorth America:Americas...CaribbeanTrinidad and TobagoTT
176AfricaAfricaAfrica:Africa...Eastern AfricaSouth SudanSS

177 rows × 12 columns

4. 操作栅格对象

4.1 栅格子集

使用 rasterio 时，可以通过 numpy 数组访问栅格值，该数组可以使用 .read 方法导入：

elev = src_elev.read(1)
elev

array([[ 1,  2,  3,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

然后，我们可以使用 numpy 方法访问单元格值的任何子集，例如：

elev[0, 0]  ## Value at row 1, column 1

1

可以通过结合子集操作覆盖现有值来修改单元格值，例如 将 elev 的左上角单元格设置为 0：

elev[0, 0] = 0
elev

array([[ 0,  2,  3,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

多个单元格也可以这样修改：

elev[0, 0:3] = 0
elev

array([[ 0,  0,  0,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

4.2 汇总栅格对象

可以通过在具有栅格值的数组上应用“numpy”汇总函数来计算栅格值的全局汇总，例如 np.mean：

np.mean(elev)

18.333333333333332

请注意，如果栅格包含需要忽略的“无数据”值，则应使用“无数据”安全函数（例如 np.nanmean）。 在我们证明这一点之前，我们必须将数组从“int”转换为“float”，因为“int”数组不能包含“np.nan”（由于计算机内存限制）：

elev1 = elev.copy()
elev1 = elev1.astype('float64')
elev1

array([[ 0.,  0.,  0.,  4.,  5.,  6.],
       [ 7.,  8.,  9., 10., 11., 12.],
       [13., 14., 15., 16., 17., 18.],
       [19., 20., 21., 22., 23., 24.],
       [25., 26., 27., 28., 29., 30.],
       [31., 32., 33., 34., 35., 36.]])

现在我们可以在数组中插入一个 np.nan 值。 （尝试在原始的 elev 数组中这样做会引发错误，请自己尝试看看）

elev1[0, 2] = np.nan
elev1

array([[ 0.,  0., nan,  4.,  5.,  6.],
       [ 7.,  8.,  9., 10., 11., 12.],
       [13., 14., 15., 16., 17., 18.],
       [19., 20., 21., 22., 23., 24.],
       [25., 26., 27., 28., 29., 30.],
       [31., 32., 33., 34., 35., 36.]])

使用 np.nan 值，汇总值变得未知：

np.mean(elev1)

nan

因此，我们需要忽略“无数据”值：

np.nanmean(elev1)

18.857142857142858

栅格值统计数据可以通过多种方式可视化。
一种方法是将栅格值“扁平化”为一维数组，然后使用图形函数，例如 plt.hist 或 plt.boxplot（来自 matplotlib.pyplot）。
例如：

x = elev.flatten()
plt.hist(x);

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-NUW9s6Av-1687280429016)(03-attribute-operations_files/03-attribute-operations_91_0.png)]

…


174
Europe
Europe
Europe:Europe
…
Southern Europe
Kosovo
XK


175
Americas
North America
North America:Americas
…
Caribbean
Trinidad and Tobago
TT


176
Africa
Africa
Africa:Africa
…
Eastern Africa
South Sudan
SS

177 rows × 12 columns

5.操作栅格对象 {#sec-manipulating-raster-objects}

5.1 栅格子集 {#sec-raster-subsetting}

使用 rasterio 时，可以通过 numpy 数组访问栅格值，该数组可以使用 .read 方法导入：

elev = src_elev.read(1)
elev

array([[ 1,  2,  3,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

然后，我们可以使用 numpy 方法访问单元格值的任何子集，例如：

elev[0, 0]  ## Value at row 1, column 1

1

可以通过结合子集操作覆盖现有值来修改单元格值，例如 将 elev 的左上角单元格设置为 0：

elev[0, 0] = 0
elev

array([[ 0,  2,  3,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

多个单元格也可以这样修改：

elev[0, 0:3] = 0
elev

array([[ 0,  0,  0,  4,  5,  6],
       [ 7,  8,  9, 10, 11, 12],
       [13, 14, 15, 16, 17, 18],
       [19, 20, 21, 22, 23, 24],
       [25, 26, 27, 28, 29, 30],
       [31, 32, 33, 34, 35, 36]], dtype=uint8)

5.2 汇总栅格对象

可以通过在具有栅格值的数组上应用“numpy”汇总函数来计算栅格值的全局汇总，例如 np.mean：

np.mean(elev)

18.333333333333332

请注意，如果栅格包含需要忽略的“无数据”值，则应使用“无数据”安全函数（例如 np.nanmean）。 在我们证明这一点之前，我们必须将数组从“int”转换为“float”，因为“int”数组不能包含“np.nan”（由于计算机内存限制）：

elev1 = elev.copy()
elev1 = elev1.astype('float64')
elev1

array([[ 0.,  0.,  0.,  4.,  5.,  6.],
       [ 7.,  8.,  9., 10., 11., 12.],
       [13., 14., 15., 16., 17., 18.],
       [19., 20., 21., 22., 23., 24.],
       [25., 26., 27., 28., 29., 30.],
       [31., 32., 33., 34., 35., 36.]])

现在我们可以在数组中插入一个 np.nan 值。 （尝试在原始的 elev 数组中这样做会引发错误，请自己尝试看看）

elev1[0, 2] = np.nan
elev1

array([[ 0.,  0., nan,  4.,  5.,  6.],
       [ 7.,  8.,  9., 10., 11., 12.],
       [13., 14., 15., 16., 17., 18.],
       [19., 20., 21., 22., 23., 24.],
       [25., 26., 27., 28., 29., 30.],
       [31., 32., 33., 34., 35., 36.]])

使用 np.nan 值，汇总值变得未知：

np.mean(elev1)

nan

因此，我们需要忽略“无数据”值：

np.nanmean(elev1)

18.857142857142858

栅格值统计数据可以通过多种方式可视化。
一种方法是将栅格值“扁平化”为一维数组，然后使用图形函数，例如 plt.hist 或 plt.boxplot（来自 matplotlib.pyplot）。
例如：

x = elev.flatten()
plt.hist(x);

在这里插入图片描述

rgb2gray
关注
  • 0
    点赞
  • 踩
  • 0
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
专栏目录
geocompkg:带有R的地理计算一书中使用的包含小插图的元数据包
05-26
土工复合材料 这是《 》一书的“元包”支持。 安装它会为您提供运行本书中的代码所需的依赖关系。 该资料库还包含最新的解决方案以及补充本书的其他材料。 请参阅以获取源代码。 安装 您可以通过安装软件包来安装...
Geocomputation （3）空间数据操作
May the Force be with you
06-23 498
空间子集化是获取空间对象并返回仅包含与另一个对象在空间中相关的特征的新对象的过程。类似于属性子集（在 @sec-vector-attribute-subsetting 中介绍），可以使用方括号（“[”）运算符使用语法“x[y]”创建“GeoDataFrame”的子集，其中“x”是 一个“GeoDataFrame”，将从中返回行/特征的子集，“y”是“子集对象”。
Geocomputation （4）几何图形操作
May the Force be with you
06-23 301
到目前为止，本书已经解释了地理数据集的结构（第 2 章），以及如何根据其非地理属性（第 3 章）和空间关系（第 4 章）来操作它们。本章重点介绍如何操作地理对象的地理元素，例如简化和转换矢量几何、裁剪栅格数据集以及将矢量对象转换为栅格以及将栅格转换为矢量。阅读完它并尝试最后的练习后，您应该理解并控制 sf 对象中的几何列以及栅格中表示的像素相对于其他地理对象的范围和地理位置。@sec-geo-vec 涵盖了使用“一元”和“二元”运算转换矢量几何形状。
Geocomputation （1）Python中的地理数据
May the Force be with you
06-20 349
正如我们将在本章后面的示例代码中看到的那样，和GeoSeries生态系统提供了一种在 Python 中处理矢量图层的综合方法，许多包都基于它构建。然而，栅格数据的情况并非如此：有几个部分重叠的包用于处理栅格数据，每个包都有自己的优点和缺点。numpydict另一个值得一提的与栅格相关的包是。它是用于处理带标签数组的通用包，因此有利于处理“复杂”栅格格式（例如 NetCDF），使用其自己的本机类表示栅格，即和本章将简要介绍基本的地理数据模型：矢量和栅格。
Geocomputation （6）重新投影地理数据
May the Force be with you
06-29 336
sec-coordinate-reference-systems 引入了坐标参考系统（CRS），重点关注两种主要类型：地理（“经度/纬度”，单位为经度和纬度）和投影（通常以米为单位） 基准）坐标系。本章以这些知识为基础并进行了更深入的探讨。它演示了如何设置地理数据并将其从一种 CRS 转换为另一种 CRS，此外，还强调了您应该注意的由于忽略 CRS 而可能出现的特定问题，特别是当您的数据使用经度/纬度坐标存储时。在许多项目中，无需担心不同的 CRS，更不用说在不同的 CRS 之间进行转换。
geocompr, 关于 R的地理数据，即将出版的书.zip
10-09
geocompr, 关于 R的地理数据，即将出版的书 带 R的 Geocomputation 简介这个库使用 R/Robin /Nowosad和 Jannes Muenchow插件来承载底层的代码，并带有。本书的在线版本是在 https://geocompr.robinlovel
计算笔记本中的城市街道网络分析-研究论文
06-10
它说明了如何为各种研究站点下载城市数据和模拟街道网络、计算网络指标、可视化街道中心性、计算路线以及使用其他空间数据（例如建筑物覆盖区和兴趣点）。 计算笔记本有助于向新用户介绍方法，并帮助研究人员接触对...
qgisprocess:使用“ QGIS”处理算法
05-12
qgisprocess qgisprocess的目标是为流行... 如果您的操作系统无法使用最新版本的QGIS，则可以使用安装了QGIS的的Geocomputation 。 如果自动配置失败（或者，如果你有一个以上的安装，并想选择哪一个是使用qgisproces
GeoComputational Modelling
08-27
Geocomputation may be viewed as the application of a computational science paradigm to study a wide range of problems in geographical systems contexts. This volume presents a clear, comprehensive and ...
【Python】Python中变量或者函数加下划线‘_‘的含义
yangweipeng708的博客
05-14 672
理解这些约定可以帮助Python开发者编写更清晰、更Pythonic的代码，并理解其他Python代码中的命名习惯。）通常用来指示特定的编程约定或用法，这些约定主要关乎变量和函数的可访问性以及预期的使用范围。在Python中，变量或函数前面加一个下划线（
python数据分析numpy基础之intersect1d求数组交集
最新发布
sinat_34735632的博客
05-19 190
python的numpy库的intersect1d(x,y)函数，计算x和y的公共元素，并返回去重后的有序结果。相当于求两个数组中的交集，并且进行去重和排序。numpy.intersect1()查找两个数组的公共元素，并且返回去重后的有序结果。ar1，ar2：必选，列表、元组、数组；如果是多维则自动转一维数组；assume_unique：可选，bool，默认False，若为True则可以加快计算素材，此时ar1和ar2需为元素唯一的数组，否则结果不会去重；
windows11 Django环境安装
云深海阔专栏
05-15 154
3）在这个文件夹中使用manager.py和runserver命令运行开发 web 服务器。2）使用 django-admin 工具创建一个名为“mytestsite”的新框架站点。浏览器上访问http://localhost:800。2、windows10 虚拟环境设置。1、验证python和pip3环境。1）创建文件夹及项目目录。4）web浏览器上验证。3、安装Django。
python数据分析——时间序列
qq_74013365的博客
05-14 1992
时间序列是按照时间顺序排列的一系列随时间变化而变化的数据点或观测值。时间序列可以是离散的，例如每月的销售数据，也可以是连续的，例如气温和股票价格等。时间序列常用于预测和分析未来的趋势，例如经济增长、股票走势、天气变化等。时间序列分析是数据分析中的重要部分，它涉及到对随时间变化的数据进行研究，以揭示其内在规律、趋势和周期性变化。Python作为一种强大的编程语言，拥有众多的数据处理和可视化库，如pandasnumpymatplotlib和seaborn等，这些库在处理时间序列数据时表现出色。
代码随想录算法训练营第二十八天|​216.组合总和III​、17.电话号码的字母组合
qq_52149213的博客
05-14 526
回溯法第二天
Gradio 案例——将 dicom 文件转为 nii文件
蒋含竹的博客
05-15 356
Gradio 的WEB界面案例：利用 SimpleITK 库，将 dicom 文件转为 nii文件
【华为OD机试-C卷D卷-200分】 简易内存池【C++/Java/Python】
goldarmour的博客
05-14 160
请实现一个简易内存池,根据请求命令完成内存分配和释放。 内存池支持两种操作命令，REQUEST和RELEASE，其格式为： REQUEST=请求的内存大小 表示请求分配指定大小内存，如果分配成功，返回分配到的内存首地址；如果内存不足，或指定的大小为0，则输出error。 RELEASE=释放的内存首地址 表示释放掉之前分配的内存，释放成功无需输出，如果释放不存在的首地址则输出error。
python EEL + vue3.js 简单示例
bigcarp的专栏
05-15 185
eel官方示例中暴露的js函数是全局函数，vue中的自定义函数作用域通常都是组件范围内。要让eel.js调用，需要将其升为全局可用。一般方法有 app.config.globalProperties 或 mixin等。
Python slice() 使用方法及示例说明
ym62033的博客
05-15 193
处。
现在地理信息系统方面有些什么前沿的理论和论文
02-16
地理信息系统（Geographic Information System，简称GIS）是一个广泛的领域，涉及许多不同的研究方向。以下是一些当前GIS领域的前沿理论和论文： 1. 地理计算（Geocomputation）：地理计算是GIS领域的一个新兴分支，它涉及到计算方法、技术和工具的发展，旨在推动GIS应用的创新和发展。相关论文：《地理计算与GIS：发展、挑战和机遇》（Geocomputation and GIS: development, challenges and opportunities） 2. 空间数据科学（Spatial Data Science）：空间数据科学是GIS领域的一个新兴研究方向，它关注如何使用大数据、机器学习、人工智能等技术来解决空间数据的挑战。相关论文：《空间数据科学：概念、挑战和机遇》（Spatial Data Science: concepts, challenges and opportunities） 3. 空间语义（Spatial Semantics）：空间语义是一种研究空间信息语义表达的方法，它涉及到语义、知识表示和推理等方面的研究。相关论文：《空间语义：概念、模型和应用》（Spatial Semantics: concepts, models and applications） 4. 空间智能（Spatial Intelligence）：空间智能是一种新兴的研究领域，它旨在开发新的空间数据处理和分析方法，以便更好地理解和利用空间信息。相关论文：《空间智能：理论、方法和应用》（Spatial Intelligence: theory, methods and applications） 5. 3D GIS：3D GIS是GIS领域的一个研究方向，涉及到3D建模、可视化和分析等方面的研究。相关论文：《3D GIS：发展、挑战和未来方向》（3D GIS: development, challenges and future directions） 这只是GIS领域的一部分前沿理论和论文，随着技术的不断发展和应用的不断深入，还会出现更多的研究成果。

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

rgb2gray

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值