Вопрос по numpy, matlab, python – найти длину последовательностей с одинаковыми значениями в массиве numpy (кодирование длин серий)

40

В программе Pylab (которая, возможно, также может быть программой Matlab) у меня есть массив цифр, представляющих расстояния:d[t] этоdistance вовремяt (и временной интервал моих данныхlen(d) единицы времени).

Интересующие меня события - это когда расстояние ниже определенного порога, и я хочу вычислить продолжительность этих событий. Легко получить массив логических значений сb = d<thresholdи проблема сводится к вычислению последовательности длин истинных слов вb, Но я не знаю, как это сделать эффективно (например, с помощью numpy примитивов), и я прибегнул к обходу массива и обнаружению изменений вручную (т. Е. Инициализировал счетчик, когда значение переходит из False в True, увеличивая счетчик, пока значение равно True и выведите счетчик в последовательность, когда значение вернется к False). Но это очень медленно.

How to efficienly detect that sort of sequences in numpy arrays ?

Ниже приведен код Python, иллюстрирующий мою проблему: четвертая точка занимает очень много времени (если нет, увеличьте размер массива)

from pylab import *

threshold = 7

print '.'
d = 10*rand(10000000)

print '.'

b = d<threshold

print '.'

durations=[]
for i in xrange(len(b)):
    if b[i] and (i==0 or not b[i-1]):
        counter=1
    if  i>0 and b[i-1] and b[i]:
        counter+=1
    if (b[i-1] and not b[i]) or i==len(b)-1:
        durations.append(counter)

print '.'

Ваш Ответ

5   ответов
0
durations = []
counter   = 0

for bool in b:
    if bool:
        counter += 1
    elif counter > 0:
        durations.append(counter)
        counter = 0

if counter > 0:
    durations.append(counter)
уверен, это более разумно, но столь же неэффективно; то, что я хочу сделать, это переместить цикл вниз на уровень C, используя некоторую умную комбинацию пустых вызовов ... Gyom
проверьте мой отредактированный ответ, теперь я предлагаю одну такую «умную комбинацию» (всегда стараясь не быть СЛИШКОМ умной ;-) - но, измерьте скорость этого И И решения на основе itertools.groupby, и позвольте мы знаем, какой из них быстрее (и на сколько) в реалистичных для вас примерах! Alex Martelli
5

ом), вот один из способов решить эту проблему в MATLAB:

threshold = 7;
d = 10*rand(1,100000);  % Sample data
b = diff([false (d < threshold) false]);
durations = find(b == -1)-find(b == 1);

Я не слишком знаком с Python, но, может быть, это поможет тебе дать некоторые идеи. знак рав

diff () тоже существует в numpy, так что это более или менее то, что вы хотите, хотя замените find (foo) на where (foo) [0]. dwf
43

В то время как неnumpy примитивы,itertools функции часто бывают очень быстрыми, поэтому попробуйте еще раз (и, конечно, измерьте время для различных решений, включая это):

def runs_of_ones(bits):
  for bit, group in itertools.groupby(bits):
    if bit: yield sum(group)

Если вам нужны значения в списке, конечно, вы можете использовать list (run_of_ones (bits)); но, возможно, понимание списка может быть еще немного быстрее:

def runs_of_ones_list(bits):
  return [sum(g) for b, g in itertools.groupby(bits) if b]

Перемещение к «тупым» возможностям, как насчет:

def runs_of_ones_array(bits):
  # make sure all runs of ones are well-bounded
  bounded = numpy.hstack(([0], bits, [0]))
  # get 1 at run starts and -1 at run ends
  difs = numpy.diff(bounded)
  run_starts, = numpy.where(difs > 0)
  run_ends, = numpy.where(difs < 0)
  return run_ends - run_starts

Снова: обязательно сравнивайте решения друг с другом на реалистичных для вас примерах!

@ gnovice, я не делаю matlab (как ни странно, моя дочь, теперь кандидат наук в области продвинутой радиотехники, делает это ;-), но теперь, глядя на ваш ответ, я вижу аналогии - получаю конец пробежек минус начало прогонов, найдите их, найдя точки <0 и> 0 в разнице, и добавьте биты к нулям, чтобы убедиться, что все прогоны заканчиваются. Думаю, не так уж много способов избавиться от этой проблемы "длины пробега"! -) Alex Martelli
@ Gyom, милости просим - как подсказывает @gnovice, решение matlab также похоже, или, я думаю, так было бы, если бы кто-то знал matlab - так должно быть, что ни один из них не очень умен; -) ... это больше вопрос того, что раньше мне приходилось делать кодирование по длине прогона (большую часть времени в моем редакторе речь шла о переводе с Numeric, к которому я до сих пор склоняюсь инстинктивно обращаться, к гораздо лучшему NumPy - но где я на самом деле впервые узнал о таких вещах в APL, 30 лет назад, когда я еще был разработчиком аппаратного обеспечения ...! -). Alex Martelli
Большое спасибо ! Решение diff / where - это именно то, что я имел в виду (не говоря уже о том, что оно примерно в 10 раз быстрее, чем другие решения). Назовите это «не слишком умно», если хотите, но я бы хотел, чтобы я был достаточно умен, чтобы придумать это: -) Gyom
36

(работает также со строками, логическими значениями и т. Д.).

Выходит кортеж длин серий, стартовых позиций и значений.

import numpy as np

def rle(inarray):
        """ run length encoding. Partial credit to R rle function. 
            Multi datatype arrays catered for including non Numpy
            returns: tuple (runlengths, startpositions, values) """
        ia = np.asarray(inarray)                  # force numpy
        n = len(ia)
        if n == 0: 
            return (None, None, None)
        else:
            y = np.array(ia[1:] != ia[:-1])     # pairwise unequal (string safe)
            i = np.append(np.where(y), n - 1)   # must include last element posi
            z = np.diff(np.append(-1, i))       # run lengths
            p = np.cumsum(np.append(0, z))[:-1] # positions
            return(z, p, ia[i])

Довольно быстро (i7):

xx = np.random.randint(0, 5, 1000000)
%timeit yy = rle(xx)
100 loops, best of 3: 18.6 ms per loop

Несколько типов данных:

rle([True, True, True, False, True, False, False])
Out[8]: 
(array([3, 1, 1, 2]),
 array([0, 3, 4, 5]),
 array([ True, False,  True, False], dtype=bool))

rle(np.array([5, 4, 4, 4, 4, 0, 0]))
Out[9]: (array([1, 4, 2]), array([0, 1, 5]), array([5, 4, 0]))

rle(["hello", "hello", "my", "friend", "okay", "okay", "bye"])
Out[10]: 
(array([2, 1, 1, 2, 1]),
 array([0, 2, 3, 4, 6]),
 array(['hello', 'my', 'friend', 'okay', 'bye'], 
       dtype='|S6'))

Те же результаты, что и у Алекса Мартелли выше:

xx = np.random.randint(0, 2, 20)

xx
Out[60]: array([1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1])

am = runs_of_ones_array(xx)

tb = rle(xx)

am
Out[63]: array([4, 5, 2, 5])

tb[0][tb[2] == 1]
Out[64]: array([4, 5, 2, 5])

%timeit runs_of_ones_array(xx)
10000 loops, best of 3: 28.5 µs per loop

%timeit rle(xx)
10000 loops, best of 3: 38.2 µs per loop

Немного медленнее, чем Алекс (но все же очень быстро), и гораздо более гибкий.

8

оно принимает массив, содержащий последовательность bools, и подсчитывает длину переходов.

>>> from numpy import array, arange
>>> b = array([0,0,0,1,1,1,0,0,0,1,1,1,1,0,0], dtype=bool)
>>> sw = (b[:-1] ^ b[1:]); print sw
[False False  True False False  True False False  True False False False
  True False]
>>> isw = arange(len(sw))[sw]; print isw
[ 2  5  8 12]
>>> lens = isw[1::2] - isw[::2]; print lens
[3 4]

sw содержит истину, где есть переключатель,isw преобразует их в индексы. Затем элементы isw попарно вычитаются вlens.

Обратите внимание, что если последовательность начинается с 1, она будет считать длину последовательностей 0 с: это можно исправить в индексировании для вычисления линзы. Кроме того, я не проверял угловые случаи таких последовательностей длины 1.

Полная функция, которая возвращает начальные позиции и длины всехTrue -Subarrays.

import numpy as np

def count_adjacent_true(arr):
    assert len(arr.shape) == 1
    assert arr.dtype == np.bool
    if arr.size == 0:
        return np.empty(0, dtype=int), np.empty(0, dtype=int)
    sw = np.insert(arr[1:] ^ arr[:-1], [0, arr.shape[0]-1], values=True)
    swi = np.arange(sw.shape[0])[sw]
    offset = 0 if arr[0] else 1
    lengths = swi[offset+1::2] - swi[offset:-1:2]
    return swi[offset:-1:2], lengths

Испытано для разных bool 1D-массивов (пустой массив; один / несколько элементов; четные / нечетные длины; начинаются сTrue/False; только сTrue/False элементы).

Похожие вопросы