Срез динамического массива на основе минимальных / максимальных значений

Question

У меня есть 3-х мерный массив хапа (365, x, y), где 36 соответствует = ежедневные данные. В некоторых случаях все элементы вдоль оси времени axis=0 равны np.nan.

Временные ряды для каждой точки вдоль axis=0 выглядят примерно так:

Мне нужно найти индекс, при котором возникает максимальное значение (данные пика), а затем два минимальных значения на каждой стороне пика.

import numpy as np

a = np.random.random(365, 3, 3) * 10
a[:, 0, 0] = np.nan

peak_mask = np.ma.masked_array(a, np.isnan(a))
peak_indexes = np.nanargmax(peak_mask, axis=0)

Я могу найти минимум перед пиком, используя что-то вроде этого:

early_minimum_indexes = np.full_like(peak_indexes, fill_value=0)

for i in range(peak_indexes.shape[0]):
    for j in range(peak_indexes.shape[1]):
        if peak_indexes[i, j] == 0:
            early_minimum_indexes[i, j] = 0
        else:
            early_mask = np.ma.masked_array(a, np.isnan(a))
            early_loc = np.nanargmin(early_mask[:peak_indexes[i, j], i, j], axis=0)   
            early_minimum_indexes[i, j] = early_loc

С полученным пиком и впадиной, построенным так:

Этот подход очень неоправдан по времени для больших массивов (1 м + элемент). Есть ли лучший способ сделать это с помощью numpy?

Answer 1

Хотя использование замаскированных массивов может быть не самым эффективным решением в этом случае в этом случае, оно позволит вам выполнять замаскированные операции на определенных осях при более или менее сохранении формы, что является большим удобством. Имейте в виду, что во многих случаях замаскированные функции по-прежнему копируют замаскированные данные.

В текущем коде у вас в основном правильная идея, но вы пропустили несколько приемов, например, возможность отрицатьи комбинировать маски. Кроме того, тот факт, что распределение масок как булевых значений впереди более эффективно, и маленькие мелкие мелочи, такие как np.full(..., 0) -> np.zeros(..., dtype=bool).

Давайте рассмотрим это в обратном порядке. Допустим, у вас был хорошо себя ведущий 1-D массив с пиком, скажем a1. Вы можете использовать маскирование, чтобы легко находить максимумы и минимумы (или индексы) следующим образом:

peak_index = np.nanargmax(a1)
mask = np.zeros(a1.size, dtype=np.bool)
mask[peak:] = True
trough_plus = np.nanargmin(np.ma.array(a1, mask=~mask))
trough_minus = np.nanargmin(np.ma.array(a1, mask=mask))

Это учитывает тот факт, что замаскированные массивы отражают смысл маски относительно обычного логического индексирования numpy. Также нормально, что максимальное значение появляется при расчете trough_plus, поскольку оно гарантированно не будет минимальным (если у вас нет ситуации с полностью нан).

Теперь, если a1 был замаскированным массивомуже (но все еще 1D), вы можете сделать то же самое, но временно объединить маски. Например:

a1 = np.ma.array(a1, mask=np.isnan(a1))
peak_index = a1.argmax()
mask = np.zeros(a1.size, dtype=np.bool)
mask[peak:] = True
trough_plus = np.ma.masked_array(a1, mask=a.mask | ~mask).argmin()
trough_minus  (np.ma.masked_array(a1, mask=a.mask | mask).argmin()

Опять же, поскольку маскированные массивы имеют обратные маски, важно комбинировать маски, используя | вместо &, как вы это делали бы для обычных булевых масок. В этом случае нет необходимости вызывать версию nan для argmax и argmin, поскольку все nans уже замаскированы.

Надеемся, что обобщение на несколько измерений станет понятным, еслираспространенность ключевого слова axis в функциях numpy:

a = np.ma.array(a, mask=np.isnan(a))
peak_indices = a.argmax(axis=0).reshape(1, *a.shape[1:])
mask = np.arange(a.shape[0]).reshape(-1, *(1,) * (a.ndim - 1)) >= peak_indices

trough_plus = np.ma.masked_array(a, mask=~mask | a.mask).argmin(axis=0)
trough_minus = np.ma.masked_array(a, mask=mask | a.mask).argmin(axis=0)

Техника N-мерного маскирования исходит из Эффективно заполнить маску на основе стартовых индексов , которые запрашивались только для этой цели.

Answer 2

Вот метод, который

1. копирует данные
2. сохраняет все позиции nan и заменяет все nans на глобальные min-1
3. находит строковый argmax
4. вычитает свое значение из всей строки
   • обратите внимание, что у каждой строки теперь есть только неположительные значения с максимальным значением, равным нулю
5. нули всеНан положения
6. переворачивает знак всех значений справа от максимума
   • это основная идея ;он создает новый глобальный максимум в той позиции, где раньше была правая рука min;в то же время он гарантирует, что левая минута теперь глобальна для строки
7. извлекает строковые argmin и argmax, это позиции левой и правой минут в исходном массиве
8. находит все-нановые строки и перезаписывает индексы max и min в этих позициях с помощью INVALINT

Code:

INVALINT = -9999
t,x,y = a.shape
t,x,y = np.ogrid[:t,:x,:y]
inval = np.isnan(a)
b = np.where(inval,np.nanmin(a)-1,a)
pk = b.argmax(axis=0)
pkval = b[pk,x,y]
b -= pkval
b[inval] = 0
b[t>pk[None]] *= -1
ltr = b.argmin(axis=0)
rtr = b.argmax(axis=0)
del b
inval = inval.all(axis=0)
pk[inval] = INVALINT
ltr[inval] = INVALINT
rtr[inval] = INVALINT

# result is now in ltr ("left trough"), pk ("peak") and rtr