Установка ступенчатой ​​функции

Я пытаюсь подобрать ступенчатую функцию, используя scipy.optimize.leastsq. Рассмотрим следующий пример:

import numpy as np
from scipy.optimize import leastsq

def fitfunc(p, x):
    y = np.zeros(x.shape)
    y[x < p[0]] = p[1]
    y[p[0] < x] = p[2]
    return y

errfunc = lambda p, x, y: fitfunc(p, x) - y # Distance to the target function

x = np.arange(1000)
y = np.random.random(1000)

y[x < 250.] -= 10

p0 = [500.,0.,0.]
p1, success = leastsq(errfunc, p0, args=(x, y))

print p1

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

[  5.00000000e+02  -4.49410173e+00   4.88624449e-01]

когда первый параметр будет оптимален при установке на 250, а второй на -10.

Кто-нибудь знает, почему это может не работать и как заставить его работать?

Если я побегу

print np.sum(errfunc(p1, x, y)**2.)
print np.sum(errfunc([250.,-10.,0.], x, y)**2.)

Я нахожу:

12547.1054663
320.679545235

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


python scipy curve-fitting best-fit-curve
person astrofrog    schedule 03.10.2009    source источник
comment
Вы уверены, что это было бы оптимальным? Каково значение errfunc для возвращенного p1 и какое оно для [250, whatever, -10], которое вы бы предпочли получить?   -  person Alex Martelli    schedule 03.10.2009
comment
добавил эту информацию в исходный пост (а не сюда, потому что код путается в комментарии)   -  person astrofrog    schedule 03.10.2009
comment
Альтернативой подгонке кривой является использование методов вейвлета для выделения признаков.   -  person Peter Mortensen    schedule 05.10.2009

Ответы (5)


arrow_upward
2
arrow_downward

Оказывается, подгонка будет намного лучше, если я добавлю аргумент epsfcn= к наименьшему квадрату:

p1, success = leastsq(errfunc, p0, args=(x, y), epsfcn=10.)

и результат

[ 248.00000146   -8.8273455     0.40818216]

Мое основное понимание состоит в том, что первый свободный параметр должен быть перемещен больше, чем расстояние между соседними точками, чтобы повлиять на квадрат остатков, а epsfcn имеет какое-то отношение к тому, насколько большие шаги использовать для поиска градиента или что-то подобное.

person astrofrog    schedule 03.10.2009
comment
Этот подход сработал и для меня; Я рекомендую принять это как ответ. - person keflavich; 31.03.2012

arrow_upward
1
arrow_downward

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

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

Для чего именно вы собираетесь это использовать? Некоторый контекст может помочь.

person duffymo    schedule 03.10.2009
comment
У меня есть данные о линейном дрейфе в зависимости от времени. В какой-то момент времени t0 дрейф внезапно подскакивает и имеет другой наклон, причем это происходит и во второй раз. Итак, что я действительно хочу, так это три линии в трех разных диапазонах. Проблема в том, что время прыжка заранее неизвестно, и мне нужно сделать это для тысяч наборов данных, поэтому я хочу, чтобы время прыжка, наклоны и пересечения линий были свободными параметрами. Я просто подумал, что начну с более простого случая. - person astrofrog; 04.10.2009
comment
Если это данные, зависящие от времени, тем больше причин использовать преобразования Фурье. Возможно, БПФ был бы более полезным. - person duffymo; 04.10.2009
comment
SciPy имеет возможности быстрого преобразования Фурье: docs.scipy.org/doc /scipy/reference/generated/ - person duffymo; 04.10.2009

arrow_upward
1
arrow_downward

Я предлагаю аппроксимировать ступенчатую функцию. Вместо бесконечного наклона в «точке изменения» сделайте его линейным на расстоянии x (1,0 в примере). Например. если параметр x, xp, для функции определен как средняя точка на этой линии, то значение при xp-0,5 является более низким значением y, а значение при xp+0,5 является более высоким значением y и промежуточными значениями функции в интервал [хр-0,5; xp+0,5] — линейная интерполяция между этими двумя точками.

Если можно предположить, что ступенчатая функция (или ее аппроксимация) переходит от более низкого значения к более высокому значению, то я думаю, что начальным предположением для последних двух параметров должно быть наименьшее значение y и наибольшее значение y соответственно вместо 0,0 и 0.0.

У меня есть 2 исправления:

1) np.random.random() возвращает случайные числа в диапазоне от 0,0 до 1,0. Таким образом, среднее значение равно +0,5, а также является значением третьего параметра (вместо 0,0). И тогда второй параметр равен -9,5 (+0,5 - 10,0) вместо -10,0.

Таким образом

print np.sum(errfunc([250.,-10.,0.], x, y)**2.)

должно быть

print np.sum(errfunc([250.,-9.5,0.5], x, y)**2.)

2) В оригинальной функции fitfunc() одно значение y становится равным 0,0, если x точно равно p[0]. Таким образом, в этом случае это не ступенчатая функция (больше похоже на сумму двух ступенчатых функций). Например. это происходит, когда начальное значение первого параметра равно 500.

person Peter Mortensen    schedule 04.10.2009

arrow_upward
0
arrow_downward

Скорее всего, ваша оптимизация застряла в локальных минимумах. Я не знаю, как на самом деле работает метод наименьшего квадрата, но если вы дадите ему начальную оценку (0, 0, 0), он тоже застрянет там.

Вы можете проверить градиент при начальной оценке численно (оценить +/- эпсилон для очень маленького эпсилон и разделить на 2 * эпсилон, взять разницу), и я уверен, что он будет около 0.

person bayer    schedule 03.10.2009

arrow_upward
0
arrow_downward

используйте statsmodel ols. ols использует метод наименьших квадратов для подбора кривой

person Golden Lion    schedule 13.06.2021