Выполнение линейного перемещения к 1D-данным в Python - PullRequest
4 голосов
/ 02 сентября 2011

У меня есть одномерный массив данных, и я хочу извлечь пространственное изменение.Стандартный способ сделать это, который я хочу Pythonize, состоит в том, чтобы выполнить движущуюся линейную регрессию к данным и сохранить градиент ...

def nssl_kdp(phidp, distance, fitlen):
    kdp=zeros(phidp.shape, dtype=float)
    myshape=kdp.shape
    for swn in range(myshape[0]):
        print "Sweep ", swn+1
        for rayn in range(myshape[1]):
            print "ray ", rayn+1
            small=[polyfit(distance[a:a+2*fitlen], phidp[swn, rayn, a:a+2*fitlen],1)[0] for a in xrange(myshape[2]-2*fitlen)]
            kdp[swn, rayn, :]=array((list(itertools.chain(*[fitlen*[small[0]], small, fitlen*[small[-1]]]))))
    return kdp

Это работает хорошо, но МЕДЛЕННО ... Мне нужно сделатьэто 17 * 360 раз ...

Я полагаю, что служебные данные находятся в итераторе в строке [for in arange] ... Есть ли следствие движущейся подгонки в numpy / scipy?

Ответы [ 3 ]

3 голосов
/ 03 сентября 2011

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

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

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

поэтому я предполагаю, что я просто предположу, что это так.

наклон равен (NΣXY - (ΣX) (ΣY)) / (NΣX2 - (ΣX) 2) где "2 "равно" в квадрате "- http://easycalculation.com/statistics/learn-regression.php

для равномерно распределенных данных, знаменатель фиксирован (поскольку вы можете сместить ось x в начало окна без изменения градиента).(ΣX) в числителе также является фиксированным (по той же причине).так что вам нужно иметь дело только с ΣXY и ΣY.последнее тривиально - просто добавьте и вычтите значение.первый уменьшается на ΣY (каждый вес X уменьшается на 1) и увеличивается на (N-1) Y (при условии, что x_0 равен 0, а x_N равен N-1) на каждом шаге.

Я подозреваю, что это не ясно.То, что я говорю, - то, что формулу для наклона не нужно полностью пересчитывать каждый шаг.в частности потому, что на каждом шаге вы можете переименовывать значения X в 0,1, ... N-1 без изменения наклона.так что почти все в формуле одинаково.все эти изменения - два слагаемых, которые зависят от Y, когда Y_0 «выпадает» из окна и Y_N «перемещается».

1 голос
/ 06 сентября 2011

Хорошо ... У меня есть то, что кажется решением ... не личный ответ, но способ сделать многоточечный дифференциал ... Я проверил это, и результат выглядит очень похожим на движущийсярегрессия ... Я использовал 1D фильтр sobel (линейное изменение от -1 до 1, свернутое с данными):

def KDP(phidp, dx, fitlen):
    kdp=np.zeros(phidp.shape, dtype=float)
    myshape=kdp.shape
    for swn in range(myshape[0]):
        #print "Sweep ", swn+1
        for rayn in range(myshape[1]):
            #print "ray ", rayn+1
            kdp[swn, rayn, :]=sobel(phidp[swn, rayn,:], window_len=fitlen)/dx
    return kdp

def sobel(x,window_len=11):
    """Sobel differential filter for calculating KDP
    output:
        differential signal (Unscaled for gate spacing
        example:
    """
    s=np.r_[x[window_len-1:0:-1],x,x[-1:-window_len:-1]]
    #print(len(s))
    w=2.0*np.arange(window_len)/(window_len-1.0) -1.0
    #print w
    w=w/(abs(w).sum())
    y=np.convolve(w,s,mode='valid')
    return -1.0*y[window_len/2:len(x)+window_len/2]/(window_len/3.0)

это работает БЫСТРО!

1 голос
/ 03 сентября 2011

С некоторым успехом я использовал эти движущиеся оконные функции из довольно старого модуля scikits.timeseries.Они реализованы в C, но мне не удалось использовать их в ситуации, когда подвижное окно различается по размеру (не уверен, нужна ли вам эта функциональность).

http://pytseries.sourceforge.net/lib.moving_funcs.html

Отправляйтесь сюда для загрузки (если вы используете Python 2.7+, вам, вероятно, понадобится скомпилировать само расширение - я сделал это для 2.7, и оно отлично работает):

http://sourceforge.net/projects/pytseries/files/scikits.timeseries/0.91.3/

I/ мы могли бы помочь вам больше, если вы немного очистите свой пример кода.Я бы рассмотрел определение некоторых аргументов / объектов в строках 7 и 8 (где вы определяете 'small') как переменные, чтобы вы не заканчивали строку 8 столь многими трудными для следования скобками.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...