Накладные расходы функции Scipy наименьший квадрат ()

Я работаю над программой анализа изображений, и я сузил свое узкое место до попыток много раз подогнать 2D-гауссова к маленькому окну (20x20) пикселей. 90% времени выполнения тратится на этот код.

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

   def gaussian(height, center_x, center_y, width_x, width_y):
           """Returns a gaussian function with the given parameters"""
        width_x = float(width_x)
        width_y = float(width_y)
        return lambda x,y: height*exp(
                    -(((center_x-x)/width_x)**2+((center_y-y)/width_y)**2)/2)


def moments(data):
    """Returns (height, x, y, width_x, width_y)
    the gaussian parameters of a 2D distribution by calculating its
    moments """
    total = data.sum()
    X, Y = indices(data.shape)
    x = (X*data).sum()/total
    y = (Y*data).sum()/total
    col = data[:, int(y)]
    width_x = sqrt(abs((arange(col.size)-y)**2*col).sum()/col.sum())
    row = data[int(x), :]
    width_y = sqrt(abs((arange(row.size)-x)**2*row).sum()/row.sum())
    height = data.max()
    return height, x, y, width_x, width_y


def fitgaussian(data):
    """Returns (height, x, y, width_x, width_y)
    the gaussian parameters of a 2D distribution found by a fit"""
    params = moments(data)
    errorfunction = lambda p: ravel(gaussian(*p)(*indices(data.shape)) -
                                 data)
    p, success = optimize.leastsq(errorfunction, params, maxfev=50, ftol=1.49012e-05)
    return p

Мне удалось вдвое сократить время выполнения, объединив функции errorfunction() и gaussian(), поэтому каждый раз, когда наименьший квадрат() вызывает errorfunction(), происходит один вызов функции вместо двух.

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

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

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


person Craig    schedule 27.01.2013    source источник
comment
если то, что вам подходит, действительно является гауссовым, вам действительно нужен вызов leastsq --- не можете ли вы просто рассчитать моменты, которые вы уже делаете.   -  person ev-br    schedule 27.01.2013
comment
Не могу поверить, что я не проверил это, спасибо за предложение! Беглый взгляд показывает, что в большинстве случаев расчетные моменты и окончательная подобранная гауссовская зависимость отличаются менее чем на 10%. Я собираюсь попробовать использовать их в своем алгоритме, чтобы увидеть, насколько хорошо он работает, и я вернусь.   -  person Craig    schedule 27.01.2013
comment
Ну, я вернулся. Кажется, что использование моментов достаточно хорошо. Спасибо Женя!   -  person Craig    schedule 27.01.2013
comment
Женя, вероятно, должна опубликовать свою рекомендацию в качестве ответа, поскольку, похоже, она ответила на вопрос...   -  person Namey    schedule 13.02.2013


Ответы (1)


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

def log_gaussian(height, center_x, center_y, width_x, width_y):
    """Returns a gaussian function with the given parameters"""
    width_x = float(width_x)
    width_y = float(width_y)
    log_height = log(height)
    return lambda x,y: (log_height - 
             (((center_x-x)/width_x)**2 - ((center_y-y)/width_y)**2)/2)

Это приведет к 1 вызову log на итерацию, а не 1 вызову exp на строку набора данных на каждой итерации.

person mtadd    schedule 01.03.2013