Импортировать большой файл .tiff в виде разреженной матрицы
У меня есть большой файл .tiff (4,4 ГБ, значения 79530 x 54980) с 1 полосой. Поскольку действительны только 16% значений, я подумал, что лучше импортировать файл как разреженную матрицу для экономии оперативной памяти. Когда я впервые открываю его как np.array, а затем преобразовываю в разреженную матрицу с помощью csr_matrix(), мое ядро уже падает. См. код ниже.
from osgeo import gdal
import numpy as np
from scipy.sparse import csr_matrix
ds = gdal.Open("file.tif")
band = ds.GetRasterBand(1)
array = np.array(band.ReadAsArray())
csr_matrix(array)
Есть ли лучший способ работы с этим файлом? В конце концов, мне приходится производить расчеты на основе значений в растре. (К сожалению, из соображений конфиденциальности я не могу прикрепить соответствующий файл.)
Насколько я знаю, в Scipy нет разреженных матриц, что приводит к значительно меньшему объему используемой памяти. Вы можете найти список здесь: docs.scipy.org/doc/scipy/reference/sparse.html . Единственные решения, которые я вижу, это: работа в сопоставленная память или работа по частям. Первый проще, но медленнее, чем второй.
Ответы 2
Предполагая, что вы знаете размер своей матрицы, вы можете создать пустую разреженную матрицу, а затем установить только допустимые значения одно за другим.
from osgeo import gdal
import numpy as np
from scipy.sparse import csr_matrix
ds = gdal.Open("file.tif")
band = ds.GetRasterBand(1)
matrix_size = (1000, 1000) # set you size
matrix = csr_matrix(matrix_size)
# for each valid value
matrix[i, j] = your_value
Изменить 1
Если вы не знаете размер своей матрицы, вы сможете проверить это следующим образом:
from osgeo import gdal
ds = gdal.Open("file.tif")
width = ds.GetRasterXSize()
height = ds.GetRasterYSize()
matrix_size = (width, height)
Редактировать 2
Измерил метрики, предложенные в комментариях (заполнил полностью). Вот как я измерял использование памяти.
размер 500x500
| матрица | пустой размер | полный размер | время заполнения |
|---|---|---|---|
| csr_matrix | 2856 | 2992 | 477,67 с |
| doc_matrix | 726 | 35807578 | 3,15 с |
| lil_matrix | 8840 | 8840 | 0,54 с |
размер 1000x1000
| матрица | пустой размер | полный размер | время заполнения |
|---|---|---|---|
| csr_matrix | 4856 | 4992 | 7164,94 с |
| doc_matrix | 726 | 150578858 | 12,81 с |
| lil_matrix | 16840 | 16840 | 2,19 с |
Вероятно, лучшим решением было бы использовать lil_matrix.
Спасибо, но я не знаю размер моей разреженной матрицы
@justinfiedler, вы должны иметь возможность проверить размер своих данных. Проверить мои изменения
Редактирование разреженной матрицы CSR должно быть безумно медленный, поскольку оно не предназначено для этого использования (особенно из-за низкой разреженности). Я думаю, что это должно занять больше оперативной памяти из-за того, что Scipy изменяет размер вектора nnz для каждого элемента (требуется двойное использование памяти CSR) и типа данных по умолчанию. Заполнение разреженной матрицы CSR 1000x1000 заняло на моей машине 3м30. Я ожидаю, что это займет несколько десятков часов для ввода OP. Матрицы LIL подходят для этого лучше, но использование памяти выше, так что в итоге это займет больше места и будет намного медленнее, чем плотная матрица...
Каждый раз, когда вы добавляете значение в матрицу CSR, вся существующая матрица копируется и все указатели вычисляются заново. Вы взяли O(N) проблему и решили ее O(2^N).
Только @CJR O(N^2) который уже большой ;) .
Чтобы было ясно, lil_matrix может быть быстрым по сравнению с итеративным построением csr_matrix, но будет использовать больше памяти, чем массив numpy (при такой разреженности), и медленнее, чем просто использование этого массива numpy. Вы только что неправильно измерили размер объектов вложенного списка (что должно быть очевидно из того факта, что добавление к нему значений не меняет использование памяти).
Как ты делаешь это? # for each valid value? csr_matrix использует np.nonzero для поиска ненулевых значений. Но если у вас уже есть массивы i, j, v для всех ненулевых значений, почему бы просто не создать массив, используя входные данные стиля coo, np.csr_matrix((v, (i,j))). Нет необходимости повторять, устанавливая одно значение за раз.
Можете ли вы сказать, где происходит сбой?
band = ds.GetRasterBand(1)
temp = band.ReadAsArray()
array = np.array(temp) # if temp is already an array, you don't need this
csr_matrix(array)
Если array равно 4,4 ГБ, (79530, 54980)
In [62]: (79530 * 54980) / 1e9
Out[62]: 4.3725594 # 4.4gB makes sense for 1 byte/element
In [63]: (79530 * 54980) * 0.16 # 16% density
Out[63]: 699609504.0 # number of nonzero values
создание csr требует выполнения np.nonzero(array) для получения индексов. Это создаст 2 массива этого размера 0,7 * 8 ГБ (индексы - 8-байтовые целые числа). coo на самом деле требуются эти 2 массива плюс 0,7 для ненулевых значений — около 12 Гб. Преобразовав в csr, атрибут row уменьшится до 79530 элементов — примерно 7 Гб. (с поправкой на 8 байт/элемент)
Таким образом, при плотности 16% разреженный формат в лучшем случае все еще больше, чем плотная версия.
Ошибка памяти при преобразовании матрицы в разреженную матрицу, указанный dtype недействителен
это недавний случай ошибки памяти, которая произошла на шаге nonzero.
Ни 0,7Гб, ни 0,7Гб на индексы массива. Вы забыли умножить размер на 4 или 8, так как индексы закодированы в int32 (по умолчанию в Windows) или int64 (по умолчанию в Linux). Так что это не 2,1 Гб, а 6,3 Гб, что слишком много для ОП и вызывает сбой из-за нехватки памяти.
Честно говоря, разреженность в 16%, вероятно, того не стоит — матрица CSR должна хранить два значения для каждого ненулевого значения. Вы экономите примерно 65 % памяти всего массива, но теряете большую часть удобства, связанного с плотным массивом, непрерывным в памяти. Если все значения 0-255, вы можете преобразовать плотный массив в
np.uint8(что составляет ~ 16% от размера стандартногоnp.int64)