Я подумал о настраиваемом формате разреженных данных. Он предназначен для эффективного хранения и загрузки, а не для вычислений. Суть в том, чтобы хранить индексы и значения ненулевых записей. Мне было интересно, есть ли какие-то настройки, которые могут улучшить производительность.
Потребность возникла из-за обработки таких данных: N «изображений» (32x32) по четыре канала каждое. Изображения содержат в среднем ~ 5% ненулевых значений. Поскольку N становится очень большим, хранение всех изображений в ОЗУ неэффективно. Таким образом, сохраняется только количество ненулевых записей, их индексы и значения, а также исходная форма.
Вот пример того, как это можно реализовать:
import numpy as np
def disassemble_data(data):
# get some dense data and make it sparse
lengths = np.count_nonzero(data, axis=(1, 2, 3))
idxs = np.flatnonzero(data)
vals = data.ravel()[idxs]
return lengths, idxs, vals, data.shape
def assemble_data(lengths, idxs, vals, shape):
# get some sparse data and make it dense
data = np.zeros(shape)
lower_idx = 0
for length in lengths:
upper_idx = lower_idx + length
data.ravel()[idxs[lower_idx:upper_idx]] = vals[lower_idx:upper_idx]
lower_idx = upper_idx
return data
# Create some dummy data
my_data = np.random.uniform(0, 1, (10, 4, 32, 32))
my_data[my_data > 0.05] = 0
# make data sparse and then dense again
my_new_data = assemble_data(*disassemble_data(my_data))
# assert that this actually works
assert np.allclose(my_data, my_new_data)
Теперь мы можем сразу увидеть преимущество: данные уплотняются изображение за изображением. Это позволяет нам загружать весь набор данных в ОЗУ и генерировать плотные изображения по запросу с помощью генератора:
def image_generator(lengths, idxs, vals, shape):
idxs %= np.prod(shape[1:])
lower_idx = 0
for length in lengths:
upper_idx = lower_idx + length
data = np.zeros(shape[1:])
data.ravel()[idxs[lower_idx:upper_idx]] = vals[lower_idx:upper_idx]
lower_idx = upper_idx
yield data
Кроме того, также можно создавать пакеты изображений:
def image_batch_generator(lengths, idxs, vals, shape, batch_size):
idxs %= np.prod((batch_size, *shape[1:]))
lengths = np.sum(lengths.reshape(-1, batch_size), axis=1)
lower_idx = 0
for length in lengths:
upper_idx = lower_idx + length
data = np.zeros((batch_size, *shape[1:]))
data.ravel()[idxs[lower_idx:upper_idx]] = vals[lower_idx:upper_idx]
lower_idx = upper_idx
yield data
Это довольно удобный подход для моих нужд. И все же мне было интересно, можно ли это ускорить.
Например. Я видел, что набор элементов numpys работает быстрее, чем прямое назначение (согласно документы). Но это работает только для одного элемента, а не для индексных массивов.
Есть ли другие подходы? Я совсем не знаком с cython и им подобным, поэтому буду рад получить несколько советов!
Это немного зависит от вашего варианта использования. Если вы не хотите использовать разреженные матрицы только в качестве «алгоритма сжатия», но работать только с плотными матрицами, используя набор данных HDF5, сжатие также может быть хорошим и более эффективным способом сделать это. (Это должно обеспечить пропускную способность выше 600 МБ / с) stackoverflow.com/a/48997927/4045774
Я знаю, что есть scipy.sparse, но представленные там форматы не соответствуют тому, который я использую. scipy.sparse предназначен для матриц, а этот метод обобщается на тензоры любой формы. Также большим отличием этого метода является возможность получения срезов разреженного тензора, поэтому нет необходимости сначала раздувать весь тензор. Есть ли лучший способ добиться этого с помощью scipy.sparse? Я использовал HDF5 со сжатием (lz4). При использовании только половины пространства декомпрессия эмпирически занимает больше времени и занимает все ядро процессора. Даже распараллеливание более 8 ядер было медленнее.






Я немного протестировал, как это можно сделать более эффективно, и пришел к выводу, что ваш подход совсем неплох для сильно некоррелированных данных, созданных np.random.uniform. На реальных данных это может быть немного иначе.
Я немного увеличил скорость ваших функций, дав около 1,4 ГБ / с для сжатия и 1,2 ГБ / с для распаковки, что совсем неплохо. С h5py (blosclz) я смог достичь только около 450 МБ / с, но записал данные также на диск.
Улучшенный разреженный алгоритм
import numpy as np
import numba as nb
#We can use uint16 on (4,32,32), since max. idx<2**16
@nb.jit()
def to_sparse_data_uint16(data):
data_flat=data.reshape(-1)
idx=np.empty(data.size,dtype=np.uint16)
data_out=np.empty(data.size,dtype=data.dtype)
ii=0
for i in range(data_flat.shape[0]):
if (data_flat[i]!=0):
idx[ii]=i
data_out[ii]=data_flat[i]
ii+=1
return idx[0:ii], data_out[0:ii], data.shape
def to_dense_data(idx,data,shape):
length=np.prod(shape)
data_out=np.zeros(length,dtype=data.dtype)
data_out[idx]=data
return data_out.reshape(shape)
########################
#do you really need float64 here?
images = np.array(np.random.uniform(0, 1, (100000, 4, 32, 32)),dtype=np.float32)
images[images > 0.05] = 0.
res=[]
t1=time.time()
for i in range(100000):
res.append(to_sparse_data_uint16(images[i,:,:,:]))
print(time.time()-t1)
t1=time.time()
for i in range(100000):
data=to_dense_data(res[i][0],res[i][1],res[i][2])
print(time.time()-t1)
Пример HDF5
import numpy as np
import tables #register blosc
import h5py as h5
import h5py_cache as h5c
import time
# I assume here that you don't need float64 for images..
# 1650MB Testdata
images = np.array(np.random.uniform(0, 1, (100000, 4, 32, 32)),dtype=np.float32)
images[images > 0.05] = 0.
#Write data (32,7 GB uncompressed)
hdf5_path='Test.h5'
f = h5c.File(hdf5_path, 'w',chunk_cache_mem_size=1024**2*100) #200 MB cache size
dset_images = f.create_dataset("images", shape=(20*100000, 4, 32, 32),dtype=np.float32,chunks=(1000, 4, 32, 32),compression=32001,compression_opts=(0, 0, 0, 0, 9, 1, 1), shuffle=False)
t1=time.time()
#Don't call h5py to often, this will lead to bad performance
for i in range(20):
dset_images[i*100000:(i+1)*100000,:,:,:]=images
f.close()
print(time.time()-t1)
print("MB/s: " + str(32700/(time.time()-t1)))
Интересно увидеть numba. Можно ли это применить и к декомпрессии? Предположение о том, что индексы малы, неверно, поскольку 2 ** 16 = 65536, к сожалению, слишком мало для набора данных из нескольких миллионов изображений. Хороший улов на float32!
В примере с декомпрессией я обращаюсь к массиву только один раз (с ограничением памяти). Также обратите внимание, что накладные расходы на вызов функции и добавление списка сильно ухудшают производительность. Сжатие занимает всего (0,65 с, список добавляет еще 0,4 с). Ad uint16: количество изображений здесь не проблема, поскольку idx начинается с нуля на каждом изображении. Если изображения не увеличиваются, этого должно быть достаточно.
Кажется, я что-то упускаю, разве цикл не увеличивает i до data.size? Итак, i = data.size-1 на последней итерации. Затем i присваивается idx [ii]. Так что это приведет к переполнению (или усечению?), Насколько я могу судить!
Тип данных рекламного изображения: обычно используются типы данных для изображений uint8, а для некоторых изображений - uint16, для которых датчик не может превышать 10–14 бит на канал. Обычно 16 бит, потому что это легче обрабатывать, но это также можно оптимизировать ... Проверьте, что здесь действительно достаточно. Пожалуйста, также очистите шаблон доступа и, возможно, не очень случайные ненулевые места и значения между ними (это также может привести к распространенным алгоритмам сжатия без потерь на переднем плане). Ваши изображения также могут быть (или нет) коррелированы с другими. Не существует самого быстрого \ лучшего алгоритма сжатия без знания точных данных.
Изображения принимаются несколькими матрицами фотоумножителей и должны быть откалиброваны, поэтому остаются дробные интенсивности пикселей (то есть количество фотоэлектронов) с высоким динамическим диапазоном. Вы совершенно правы насчет корреляции. Однако я уже пробовал использовать сжатое хранилище фрагментов HDF5 с LZ4 (из пакета 'hdf5plugin'), это было примерно в 20 раз медленнее, чем разреженный формат. Я попробую установить blosc и размер блока.
В вашем примере каждое изображение (4x32x32), поэтому idx не может выходить за рамки 4 * 32 * 32-1 = 4095 (и ii не за 4096). Как сказано выше, uint16 подходит только в том случае, если np.prod (dims_of_single_image) <2 ** 16. Пожалуйста, извините, у вас, очевидно, нет 3 каналов ... Использование python-blosc для прямого сжатия фрагментов изображений также возможно и может достигать скорости более 2-3 ГБ / с на четырехъядерном процессоре и намного выше, чем на мощной рабочей станции. Я могу завтра добавить пример. github.com/Blosc/python-blosc
Ааа, я вижу. Я полностью неправильно прочитал ваш код, извините! Ваш код предназначен для сжатия каждого изображения по отдельности и сохранения их в списке. Мой код берет все изображения сразу и сохраняет их все в одном большом массиве. Затем я просто сохраняю length, idxs, vals, data.shape с numpy.savez в одном файле и получаю невероятно высокую скорость чтения. Я добавлю тест позже.
Вы только что заново изобрели scipy.sparse? Единственная разница в том, что он фактически поддерживает вычисления, преобразованные в плотные по мере необходимости.