iso*_*eel 14 python memory sorting performance numpy
我需要使用numpy对非常大的基因组数据集进行排序.我有一个26亿浮点数的数组,维度= (868940742, 3)一旦加载并且只是坐在那里,我的机器上占用了大约20GB的内存.我有一台2015年初的13英寸MacBook Pro,配备16GB内存,500GB固态高清和3.1 GHz intel i7处理器.只是加载数组溢出到虚拟内存,但没有到我的机器遭受的程度,或者我必须停止我正在做的其他事情.
我从22个较小的(N, 2)子阵列逐步构建了这个非常大的数组.
函数使用我调用的22个子数组中的每一个FUN_1生成2个新(N, 1)数组sub_arr.
第一个输出FUN_1是通过插入来自sub_arr[:,0]数组的值生成的b = array([X, F(X)]),第二个输出是通过sub_arr[:, 0]使用数组放入二进制数来生成的r = array([X, BIN(X)]).我分别称这些输出b_arr和rate_arr.该函数返回一个3元组的(N, 1)数组:
import numpy as np
def FUN_1(sub_arr):
"""interpolate b values and rates based on position in sub_arr"""
b = np.load(bfile)
r = np.load(rfile)
b_arr = np.interp(sub_arr[:,0], b[:,0], b[:,1])
rate_arr = np.searchsorted(r[:,0], sub_arr[:,0]) # HUGE efficiency gain over np.digitize...
return r[rate_r, 1], b_arr, sub_arr[:,1]
我在for循环中调用该函数22次,并full_arr = numpy.zeros([868940742, 3])使用以下值填充预先分配的零数组:
full_arr[:,0], full_arr[:,1], full_arr[:,2] = FUN_1
在这一步节省内存方面,我认为这是我能做的最好的,但我愿意接受建议.无论哪种方式,我都没有遇到问题直到这一点,它只需要大约2分钟.
这是排序例程(有两个连续的排序)
for idx in range(2):
sort_idx = numpy.argsort(full_arr[:,idx])
full_arr = full_arr[sort_idx]
# ...
# <additional processing, return small (1000, 3) array of stats>
现在这种情况一直在起作用,虽然很慢(大约需要10分钟).然而,我最近开始使用的一个更大,更精细分辨率表[X, F(X)]上面插值步长值FUN_1返回b_arr现在的排序确实放缓,但一切仍然是相同的.
有趣的是,我甚至没有在排序现在滞后的步骤中对插值进行排序.以下是不同插值文件的一些片段 - 较小的一个在每种情况下小约30%,在第二列中的值更加均匀; 较慢的一个具有更高的分辨率和更多的唯一值,因此插值的结果可能更独特,但我不确定这是否应该有任何影响......？
更大,更慢的文件:
17399307 99.4
17493652 98.8
17570460 98.2
17575180 97.6
17577127 97
17578255 96.4
17580576 95.8
17583028 95.2
17583699 94.6
17584172 94
更小,更统一的常规文件:
1 24
1001 24
2001 24
3001 24
4001 24
5001 24
6001 24
7001 24
我不确定是什么导致这个问题,我会对任何建议感兴趣,或者只是对这种类型的内存限制情况下的排序进行一般性输入!
ali*_*i_m 11
目前,每次调用np.argsort都会生成一个(868940742, 1)int64索引数组,这将自己占用大约7 GB.此外,当您使用这些索引对列进行排序时,full_arr您正在生成另一个(868940742, 1)浮点数组,因为花式索引始终返回副本而不是视图.
一个相当明显的改进是full_arr使用其.sort()方法进行排序.不幸的是,.sort()不允许您直接指定要排序的行或列.但是,您可以指定要对结构化数组进行排序的字段.因此,您可以通过将view数组作为具有三个浮点字段的结构化数组添加到三个列之一来强制进行就地排序,然后按以下字段之一进行排序:
full_arr.view('f8, f8, f8').sort(order=['f0'], axis=0)
在这种情况下,我正在full_arr通过第0个字段进行排序,该字段对应于第一列.请注意,我假设有三个float64列('f8') - 如果你的dtype不同,你应该相应地改变它.这也要求您的数组是连续的并且采用行主格式,即full_arr.flags.C_CONTIGUOUS == True.
信用此方法应该去乔金顿他的回答在这里.
虽然它需要更少的内存,但是与使用np.argsort生成索引数组相比,按字段排序结构化数组的速度要慢得多,正如您在下面的评论中所提到的(参见前一个问题).如果您使用np.argsort获取一组索引进行排序,则可能会通过使用np.take而不是直接索引来获取排序数组,从而获得适度的性能提升:
%%timeit -n 1 -r 100 x = np.random.randn(10000, 2); idx = x[:, 0].argsort()
x[idx]
# 1 loops, best of 100: 148 µs per loop
%%timeit -n 1 -r 100 x = np.random.randn(10000, 2); idx = x[:, 0].argsort()
np.take(x, idx, axis=0)
# 1 loops, best of 100: 42.9 µs per loop
但是我不希望在内存使用方面看到任何差异,因为两种方法都会生成副本.
关于为什么对第二个数组进行排序更快的问题 - 是的,当数组中的唯一值较少时,您应该期望任何合理的排序算法更快,因为平均来说它的工作量较少.假设我有一个1到10之间的随机数字序列:
5 1 4 8 10 2 6 9 7 3
有10个!= 3628800可能的方法来排列这些数字,但只有一个数字按升序排列.现在假设只有5个唯一数字:
4 4 3 2 3 1 2 5 1 5
现在有2⁵= 32种方式按升序排列这些数字,因为我可以在排序的向量中交换任何一对相同的数字而不会破坏排序.
默认情况下,np.ndarray.sort()使用Quicksort.该qsort算法的变体通过递归地选择数组中的'pivot'元素,然后重新排序数组,使得小于枢轴值的所有元素都放在它之前,并且所有大于枢轴值的元素放在之后它.已经对已经等于枢轴的值进行了排序.具有较少的唯一值意味着平均而言,在任何给定扫描中,更多值将等于枢轴值,因此完全排序阵列需要更少的扫描.
例如:
%%timeit -n 1 -r 100 x = np.random.random_integers(0, 10, 100000)
x.sort()
# 1 loops, best of 100: 2.3 ms per loop
%%timeit -n 1 -r 100 x = np.random.random_integers(0, 1000, 100000)
x.sort()
# 1 loops, best of 100: 4.62 ms per loop
在此示例中,两个数组的dtypes相同.如果较小的阵列与较大的阵列相比具有较小的项目大小,那么由于花式索引而复制它的成本也会更小.
numpy遇到这篇文章,我想指出考虑np.dtype你正在使用的重要性.在我的情况下,我实际上能够使用半精度浮点数,np.float16即将内存中的20GB对象减少到5GB并使排序更易于管理.numpyis 使用的默认值np.float64,这是您可能不需要的很多精度.查看此处的doc,其中描述了不同数据类型的容量.感谢@ali_m在评论中指出了这一点.我做了一个糟糕的工作来解释这个问题,但我发现了一些有用的解决方法,我认为这对于任何需要对真正庞大的numpy数组进行排序的人来说都是有用的.
我正在numpy从包含这些元素的人类基因组数据的22个"子阵列" 构建一个非常大的阵列[position, value].最终,最终数组必须根据特定列中的值进行"就地"数字排序,而不必对行内的值进行混洗.
子数组维度遵循以下形式:
arr1.shape = (N1, 2)
...
arr22.shape = (N22, 2)
sum([N1..N2]) = 868940742 即有近1BN的位置要排序.
首先,我使用函数处理22个子数组,该函数process_sub_arrs返回与输入长度相同的3元组1D数组.我将1D数组堆叠到一个新(N, 3)数组中,并将它们插入到np.zeros为完整数据集初始化的数组中:
full_arr = np.zeros([868940742, 3])
i, j = 0, 0
for arr in list(arr1..arr22):
# indices (i, j) incremented at each loop based on sub-array size
j += len(arr)
full_arr[i:j, :] = np.column_stack( process_sub_arrs(arr) )
i = j
return full_arr
full_arr如下:full_arr = np.zeros([868940742, 3], dtype=np.float16),它只是大小的1/4,更容易排序.结果是一个巨大的20GB阵列:
full_arr.nbytes = 20854577808
正如@ali_m在他的详细帖子中指出的那样,我早期的例行程序效率低下:
sort_idx = np.argsort(full_arr[:,idx])
full_arr = full_arr[sort_idx]
数组sort_idx的大小为33%full_arr,在排序后会挂起并浪费内存full_arr.据推测,这种类型会产生full_arr"花式"索引的副本,可能会将内存使用量提高到已经用于保存大规模阵列的233%!这是一个缓慢的步骤,持续大约十分钟,并严重依赖虚拟内存.
我不确定"花式"排序是否会产生持久性副本.观察我的机器上的内存使用情况,似乎full_arr = full_arr[sort_idx]删除了对未分类原件的引用,因为在大约1秒后剩下的就是分类数组和索引使用的内存,即使存在瞬态副本.
这个argsort()用于节省内存的更紧凑的用法是:
full_arr = full_arr[full_arr[:,idx].argsort()]
这仍然会导致分配时出现峰值,其中瞬态索引数组和瞬态副本都已生成,但内存几乎立即被释放.
@ali_m指出了一个很好的技巧(记入Joe Kington)用于生成带有on的事实上的结构化数组.好处是这些可以"就地"排序,保持稳定的行顺序:viewfull_arr
full_arr.view('f8, f8, f8').sort(order=['f0'], axis=0)
视图对于执行数学数组操作非常有用,但是对于排序来说,即使是来自我的数据集的单个子数组也是如此.一般来说,结构化数组似乎不能很好地扩展,即使它们具有非常有用的属性.如果有人知道为什么这是我有兴趣知道.
使用非常大的数组最小化内存消耗和提高性能的一个好方法是构建一个小而简单的函数管道.函数一旦完成就会清除局部变量,所以如果中间数据结构正在构建并且使内存损坏,这可能是一个很好的解决方案.
这是我用来加速大规模数组排序的管道草图:
def process_sub_arrs(arr):
"""process a sub-array and return a 3-tuple of 1D values arrays"""
return values1, values2, values3
def build_arr():
"""build the initial array by joining processed sub-arrays"""
full_arr = np.zeros([868940742, 3])
i, j = 0, 0
for arr in list(arr1..arr22):
# indices (i, j) incremented at each loop based on sub-array size
j += len(arr)
full_arr[i:j, :] = np.column_stack( process_sub_arrs(arr) )
i = j
return full_arr
def sort_arr():
"""return full_arr and sort_idx"""
full_arr = build_arr()
sort_idx = np.argsort(full_arr[:, index])
return full_arr[sort_idx]
def get_sorted_arr():
"""call through nested functions to return the sorted array"""
sorted_arr = sort_arr()
<process sorted_arr>
return statistics
调用堆栈:get_sorted_arr - > sort_arr - > build_arr - > process_sub_arrs
一旦完成每个内部函数,get_sorted_arr()最后只需保存已排序的数组,然后返回一小组统计信息.
dtype的代表你的巨大数组,你将需要使用更高的精度进行汇总计算.例如,因为full_arr.dtype = np.float16该命令np.mean(full_arr[:,idx])试图计算半精度浮点的平均值,但是当对一个大型数组进行求和时,这会快速溢出.使用np.mean(full_arr[:,idx], dtype=np.float64)将防止溢出.我最初发布这个问题是因为我感到困惑的是,相同大小的数据集突然开始窒息我的系统内存,尽管新"慢"集中的唯一值的比例存在很大差异.@ali_m指出,实际上,具有较少唯一值的更均匀数据更容易排序:
Quicksort的qsort变体通过递归地选择数组中的'pivot'元素,然后重新排序数组,使得小于pivot值的所有元素都放在它之前,并且所有大于pivot值的元素都放在之后它.与枢轴相等的值已经排序,因此直观地说,阵列中的唯一值越少,需要进行的交换次数就越少.
在这一点上,我最终试图解决这个问题的最后一个改变是提前舍入新的数据集,因为插值步骤中存在不必要的高级小数精度.这最终比其他节省内存的步骤具有更大的影响,表明排序算法本身就是这种情况下的限制因素.
期待任何人可能对此主题提出的其他意见或建议,我几乎肯定会对某些技术问题感到错误,所以我很乐意回复:-)