熊猫（子）数据框中的最大值和最小值答案

【问题标题】：Max and Min values within pandas (sub)Dataframe熊猫（子）数据框中的最大值和最小值
【发布时间】：2020-09-23 20:58:45
【问题描述】：

我有以下数据框 -df ：

                     crs         Band1 level
lat       lon                               
34.595694 32.929028  b''  4.000000e+00  1000
          32.937361  b''  1.200000e+01  950
          32.945694  b''  2.900000e+01  925
34.604028 32.929028  b''  7.000000e+00  1000
          32.937361  b''  1.300000e+01  950
                 ...           ...   ...
71.179028 25.679028  b''  6.000000e+01  750
71.187361 25.662361  b''  1.000000e+00  725
          25.670694  b''  6.000000e+01  1000
          25.679028  b''  4.000000e+01  800
71.529028 19.387361  b''  1.843913e-38  1000

[17671817 rows x 3 columns]

还有两个数组：

lon1=np.arange(-11,47,0.25)
lat1=np.arange(71.5,34.5,-0.25)

这两个数组（lat1，lon1）产生间隔为 0.25 度的坐标对。

数据框df 包含点（lat，lon），这些点在lon1 和lat1 数组定义的点内密集分布。我想做的是：

在lat1,lon1 定义的点的 0.125 度范围内查找（过滤）来自 df 的所有点
从此子数据框中获取level 的max 和min 值，并将它们存储在与lon1 和lat1 大小相同的单独数组中。

到目前为止我所做的是过滤数据框：

for x1 in lon1:
    for y1 in lat1:
        df3=df[(df.index.get_level_values('lon')>x1-0.125) & (df.index.get_level_values('lon')<x1+0.125)]
        df3=df3[(df3.index.get_level_values('lat')>y1-0.125) & (df3.index.get_level_values('lat')<y1+0.125)]

但这性能很慢。我相信有一个更快的。我也标记了 scikit-learn，因为可能可以用它来完成，但我缺乏这个包的经验。任何帮助表示赞赏。

【问题讨论】：

lon1 和 lat1 是什么关系？是笛卡尔多派系（lat1 中的每个值都具有 lon1 的所有值）还是匹配索引？（lon1[0] 与 lat1[0]，与索引 1 等相同）
lon1 和 lat1 是点的经度和纬度。 lon1[0], lat1[0] 用经纬度描述一个点。
也许this 有助于提高速度；显然对数据进行排序应该已经有所帮助
还有一个问题：你想要一个子数据帧，还是多个子数据帧（以 lon1 数组的大小），每个找到一个最大值和最小值？我的意思是，你想要范围内所有点的最大值和最小值，还是每对的 df 的最大值和最小值？
df1 快照的链接（30k 记录）。：1drv.ms/u/s!Aj5DfuTWg1OMjtskbyjp5bR-Ga6tbg?e=Zfq5oZ

标签： python pandas scikit-learn

【解决方案1】：

在开始之前，让我们将您的 bin 转换为每个 bin 的开头而不是中心：

lon1=np.arange(-11.125,47.125,0.25)
lat1=np.arange(71.625,34.125,-0.25)

为每一行分配纬度和经度箱（注意lat1的倒序，否则你需要将ordered=False传递给pd.cut()）。

df['latcat'] = pd.cut(df.index.get_level_values(0), lat1[::-1])
df['loncat'] = pd.cut(df.index.get_level_values(1), lon1)

对于您的示例数据，我们现在拥有：

                     crs         Band1  level            latcat            loncat
lat       lon                                                                    
34.595694 32.929028  b''  4.000000e+00   1000  (34.375, 34.625]  (32.875, 33.125]
          32.937361  b''  1.200000e+01    950  (34.375, 34.625]  (32.875, 33.125]
          32.945694  b''  2.900000e+01    925  (34.375, 34.625]  (32.875, 33.125]
34.604028 32.929028  b''  7.000000e+00   1000  (34.375, 34.625]  (32.875, 33.125]
          32.937361  b''  1.300000e+01    950  (34.375, 34.625]  (32.875, 33.125]
71.179028 25.679028  b''  6.000000e+01    750  (71.125, 71.375]  (25.625, 25.875]
71.187361 25.662361  b''  1.000000e+00    725  (71.125, 71.375]  (25.625, 25.875]
          25.670694  b''  6.000000e+01   1000  (71.125, 71.375]  (25.625, 25.875]
          25.679028  b''  4.000000e+01    800  (71.125, 71.375]  (25.625, 25.875]
71.529028 19.387361  b''  1.843913e-38   1000  (71.375, 71.625]  (19.375, 19.625]

现在使用 groupby 来获取每个区域的最低和最高等级：

res = df.groupby([df.latcat.cat.codes, df.loncat.cat.codes])['level'].agg(['min', 'max'])

这给了你：

          min   max
0   176   925  1000
147 147   725  1000
148 122  1000  1000

索引的第一级是反向lat1 数组中的位置，-1 表示您的一些示例数据所在的“超出范围”。第二层是lon1数组中的位置。

按要求转换为矩阵：

minlevel = np.full((len(lat1), len(lon1)), np.nan)
maxlevel = np.full((len(lat1), len(lon1)), np.nan)
x = len(lat1) - res.index.get_level_values(0) - 1 # reverse to original order
y = res.index.get_level_values(1)
minlevel[x, y] = res['min']
maxlevel[x, y] = res['max']

【讨论】：

你给你的方法计时了吗？不知怎的，我自己无法让它工作
17671817 行大约需要 10 秒，如您的答案所示。毫不奇怪，这更快，因为它避免了所有 for 循环。

【解决方案2】：

首先让我们回顾一下您的解决方案：对于 lon1 中的每个值和 lat1 的每个值（如果它们的大小为 n，则为 n^2 次迭代），您尝试过滤数据帧，从而扫描整个 df：您的代码在数据帧中运行了 n^2 次，效率低下。

我的解决方案只需要扫描一次数据帧，每次扫描都会执行 n 次操作。它使用 pandas apply 函数，效率不高，但没有它我找不到办法。我很想听听不使用 apply 的过滤解决方案。

我使用了一个可重现的小示例，您可能需要调整索引以匹配您的代码。我相信这个例子更容易理解。

import pandas as pd
import numpy as np

df = pd.DataFrame({"lat":[22.5, 10.76, 7.341, 22.5], "log":[3.64, 7.234, 135, 3.644], "level":[2, 8, 19, 9]})

lat1 = np.array([22.51, 7.33])
lon1 = np.array([3.6, 135.02])

接下来的几行创建一个元组列表，每个元组都包含一个pandas.Interval 对象。这里的元组表示 (lat1[i]+-x, lon1[i]+-x)。请注意，我不必使用 pandas.Interval - 我可以构建另一个 (lat1[i]-x, lat1[i]+x) 元组。但我决定用 pandas 间隔去，其实没关系。结果：对于每一对[lat1, lon1]，我们有一个两个pandas区间的元组，每个都是+-0.125

interval_list = []
const_add = 0.125
for i, item in enumerate(lat1):
    interval_list.append((pd.Interval(left=lat1[i]-const_add, right=lat1[i]+const_add),pd.Interval(left=lon1[i]-const_add, right=lon1[i]+const_add)))

现在我们要过滤数据框。为了使用apply，我创建了一个自定义函数：它检查当前行是否在一个元组中，如果是，则返回 lat1 数组中的索引（稍后您会看到它为什么有用）

def within_range(row, interval_list):
    for i, item in enumerate(interval_list):
        if row[0] in item[0] and row[1] in item[1]:
            return i
    return np.nan

df["point"] = df.apply(lambda x: within_range(x, interval_list), axis=1)

在代码的那个点，我们有一个列名“point”。它的值如下：如果行靠近点 i（其中 i 是 lat1[i] 和 lon1[i] 中的索引），则值为 i。如果没有近点，则值为nan。

现在剩下的就是找到每个点的最大值和最小值，这可以使用groupby轻松实现：

max_series = df.groupby(by="point")["level"].max()
min_series = df.groupby(by="point")["level"].min()

您有两个系列的索引与 lat1 和 lon[1] 中的索引相同。您可以使用Series.array 轻松地将它们转换为数组。值得一提的是，您没有说如何处理缺失值 - 如果 df 中没有任何点接近点 (lat1[50], lon1[50])，那么最大值和最小值数组中的值是多少？这就是我将其保留为系列的原因，我相信在将其更改为数组之前对其进行操作会更容易。

整个代码在一起：

import pandas as pd
import numpy as np

df = pd.DataFrame({"lat":[22.5, 10.76, 7.341, 22.5], "log":[3.64, 7.234, 135, 3.644], "level":[2, 8, 19, 9]})

lat1 = np.array([22.51, 7.33])
lon1 = np.array([3.6, 135.02])

interval_list = []
const_add = 0.125

for i, item in enumerate(lat1):
    interval_list.append((pd.Interval(left=lat1[i]-const_add, right=lat1[i]+const_add),pd.Interval(left=lon1[i]-const_add, right=lon1[i]+const_add)))

def within_range(row, interval_list):
    for i, item in enumerate(interval_list):
        if row[0] in item[0] and row[1] in item[1]:
            return i
    return np.nan

df["point"] = df.apply(lambda x: within_range(x, interval_list), axis=1)
max_arr = df.groupby(by="point")["level"].max()
min_arr = df.groupby(by="point")["level"].min()
# or:
# max_arr = df.groupby(by="point")["level"].max().array

【讨论】：

由于apply 确实是伪装的for 循环，这似乎很慢。
我知道，我云没有找到更好的解决方案。我看到你的没用过，我会尽快过去的

【解决方案3】：

我使用this 答案中描述的技巧来有效地获取对应于一维 bin 的索引，然后遍历 lon 和 lat 的组以获取两者的交集。我在这里使用numpy，不直接申请min / max，而是专注于索引。

import numpy as np
from scipy.sparse import csr_matrix

def digitize_group(x, bins):
    idx_x = np.digitize(x, bins)
    n, m = len(x), len(bins) + 1
    s = csr_matrix((np.arange(n), [idx_x, np.arange(n)]), shape=(m, n))
    return [group for group in np.split(s.data, s.indptr[1:-1])]

# Create dummy data
n = 100000  # 17671817
step = 0.25  # Note the shift by step/2 to transform your arrays to bins
bins_lon = np.arange(-11-step/2, 47+step/2, step) 
bins_lat = np.arange(71.5+step/2, 34.5-step/2, -step)
lon = np.random.uniform(low=bins_lon.min(), high=bins_lon.max(), size=n)
lat = np.random.uniform(low=bins_lat.min(), high=bins_lat.max(), size=n)

# Get the 1D groups
group_lon = digitize_group(lon, bins_lon)
group_lat = digitize_group(lat, bins_lat)

# Combine to 2D groups
group_lonlat = np.zeros((len(group_lon), len(group_lat)), dtype=object)
for i, lo in enumerate(group_lon):
    for j, la in enumerate(group_lat):
        group_lonlat[i, j] = np.intersect1d(lo, la, assume_unique=True)

print(group_lonlat[13, 17])
# array([   15606,   131039,   168479,   171734,   174281,   266717,   ....

通过访问group_lonlat[i, j]，您将获得一个索引列表K，其中每个元素k 都满足：

bins_lon[i] < lon[k] < bins_lon[i+1] & bins_lat[j] < lat[k] < bins_lat[j+1]

通过这些索引，您可以访问您的数据框并执行所有进一步的计算。

我的一台笔记本电脑用180s 计算n=17671817 的索引。

这种方法的一个瓶颈是交集搜索的次优处理。 sortednp 承诺在这里比 numpy 做得更好。而对于较大的n，删除已使用的索引以加快搜索速度会更有效。

import sortednp as snp
for i in range(len(group_lon)):
    for j in range(len(group_lat)):
        group_lonlat[i, j], (ii, jj) = snp.intersect(group_lon[i], group_lat[j], 
                                                     indices=True)
        group_lon[i] = np.delete(group_lon[i], ii)
        group_lat[j] = np.delete(group_lat[j], jj)

这使我们下降到 20s 为 n=17671817 和 0.3s 为 n=30000 。

【讨论】：