Использование метода наименьших квадратов Scipy для наилучшего сопоставления данных со смещенной кривой мощности

Я новичок в Python/Scipy/Numpy.

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

Похоже, что смещенная кривая мощности подходит лучше, но у меня проблемы с синтаксисом функции a(x-b)^c, и я не уверен, что это будет работать с методом log10.

Текущий код:

from scipy import *
from scipy import optimize
#                                         # READ DATA FROM CSV FILE
DIA, FlowRate =  genfromtxt('#ThisDataFile#', delimiter=',', unpack='true')
#                                         # CONVERT DATA TO LOG ITS MORE ACCURATE
logx = log10(DIA)
logy = log10(FlowRate)                    #crlf#
#                                         # define our (line) fitting function
fitfunc = lambda p, x: p[0] + p[1] * x
errfunc = lambda p, x, y: (y - fitfunc(p, x))
pinit = [5.0, -5.0]
out = optimize.leastsq(errfunc, pinit,
                       args=(logx, logy), full_output=1)
pfinal1 = out[0]
amp1 = 10.0**pfinal1[0]

print(amp1, pfinal1[1])

Может кто-то указать мне верное направление.

Спасибо Питер

Обновлять:

Пример данных (по измеренному) давлению:

DIA  = [ 1, 2, 3, 4, 5, 6, 7, 8 ]
flow = [ 58.33254,
         30.11954,
         16.02723,
          9.47614,
          5.75362,
          3.63373,
          2.37532,
          1.58426
         ]

numpy least-squares
person PeterC    schedule 25.03.2017    source источник
comment
Я подгоняю данные, которые вы опубликовали, к двойной экспоненте, журналы не нужны, с отличными результатами: y = a * exp (bx) + c * exp (dx). 2,1159998407921812E+01 b = -3,3158309848989914E-01 c = 9,7577758407802179E+01 d = -8,1563989368844170E-01   -  person James Phillips    schedule 26.03.2017
comment
Спасибо, Джеймс. Я попробовал метод двойной экспоненты в Excel, но не смог приблизиться.   -  person PeterC    schedule 26.03.2017
comment
Пробовал: из импорта scipy * из импорта scipy оптимизировать импорт numpy as np ##ЧИТАТЬ ДАННЫЕ ИЗ ФАЙЛА CSV x, y = genfromtxt('#ThisDataFile#', delimiter=',', unpack='true') def f(x, a, b, c, d): вернуть a * np.exp(b * x) + c * np.exp(d * x) def остаток(p, x, y): вернуть y - f(x, *p ) p0 = [1., 1., 1.,1.] popt, pcov = optim.leastsq(residual, p0, args=(x, y)) print (popt) ...... но получил предупреждение слишком много итераций, и результат был далеко ............ как я говорю, я новичок в этом методе, поэтому я уверен, что это мой подход больше, чем ваш метод.   -  person PeterC    schedule 26.03.2017
comment
Просто чтобы уточнить, данные представляют собой изменения давления через отверстия разного размера для фиксированного расхода. Я упростил x, так как каждый из 8 четных шагов представляет 1/32. Фактические размеры в дюймах были 0,15625,0,1875,0,21875,0,25,0,28125,0,3125,0,34375,0,375. Используя простой топор ^ b, я нахожу результаты близкими, но недостаточно близкими для наименьшего и наибольшего диаметров. Если я наберу их в Curve Expert Pro, кривая Shifted Power будет очень точной, но я не могу воспроизвести ее в Scimpy. У меня есть около 5000 наборов данных для запуска, и я автоматизировал их с помощью процедуры командной строки, вызванной из ColdFusion.   -  person PeterC    schedule 26.03.2017
comment
Решатель Excel не очень известен - я использовал решатель Python с интерфейсом генетического алгоритма для решателя Левенберга-Марквардта scipy. Генетический алгоритм искал в пространстве ошибок начальные оценки параметров для решателя LM, решатель Excel этого не делает. Если вы проверите значения уравнений и параметров из моего комментария, вы увидите отличные плавные результаты.   -  person James Phillips    schedule 27.03.2017
comment
Я только что понял, что решатель Excel должен работать, если вы запустите его с начальными значениями, близкими к тем, которые я дал для двойной экспоненты. Это стоит быстрого теста. Вам не нужно будет вводить все значение параметра, просто используйте его до первого десятичного знака.   -  person James Phillips    schedule 27.03.2017

Ответы (1)


arrow_upward
1
arrow_downward

Пример данных давления = "1,2,3,4,5,6,7,8" поток="58.33254,30.11954,16.02723,9.47614,5.75362,3.63373,2.37532,1.58426"

person PeterC    schedule 25.03.2017
comment
Извините, новичок на форумах - я ответил вместо того, чтобы добавить комментарий! - person PeterC; 25.03.2017