Python подходит две кривые одновременно с одной и той же переменной и параметрами - PullRequest
0 голосов
/ 15 ноября 2018

Я пытаюсь изогнуть фитинг с двумя функциями одновременно.У меня есть xdata , и оба из mdata и ndata являются функциями xdata с такими же параметрами (r1, r2, r3, l), как показано в коде.Теперь я могу построить график mdata и ndata относительно xdata, а также рисунок ndata относительно t mdata.но результаты не были хорошими.код ниже

две функции:

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import leastsq
from scipy.optimize import curve_fit
from scipy.optimize import differential_evolution
# from scipy.signal import medfilt
import warnings

def func(x, r1, r2, r3,l):
    # "brick wall" ensuring all parameters are positive
    if r1 < 0.0 or r2 < 0.0 or r3 < 0.0 or l < 0.0 :
        return 1.0E10 # large value gives large error, curve_fit hits a brick wall
    c = 47e-8
    w=2*np.pi*x
    m=r1+(r2*l*w)/(r2**2+l**2*w**2)+r3/(1+r3*c**2*w**2)
    return m

def gunc(x, r1, r2, r3,l):
    c = 47e-8
    w = 2 * np.pi * x
    n = (r2 ** 2 * l * w) / (r2 ** 2 + l ** 2 * w ** 2) - r3 ** 3 * c * w / (1 + r3 * c ** 2 * w ** 2)
    return n

, а затем была определена остаточная функция для минимизации

def residual_two_functions(pars, x, mdata, ndata):
    r1 = pars[0]
    r2 = pars[1]
    r3 = pars[2]
    l = pars[3]
    # c = pars[4]
    diff1 = mdata - func(x, r1, r2, r3, l)
    diff2 = ndata - gunc(x, r1, r2, r3, l)
    return np.concatenate((diff1, diff2))

остальной код

def readdata(filename):
    x = filename.readlines()
    x = list(map(lambda s: s.strip(), x))
    x = list(map(float, x))
    return x

 # test data
f_x= open(r'C:\Users\Desktop\ageingmodel\simpletry\fre.txt')
xdata = readdata(f_x)

f_m= open(r'C:\Users\Desktop\ageingmodel\simpletry\real.txt')
mdata = readdata(f_m)

f_n= open(r'C:\Users\Desktop\ageingmodel\simpletry\imag.txt')
ndata = readdata(f_n)

xdata = np.array(xdata)
mdata = np.array(mdata)
ndata = np.array(ndata)

# mdata.flatten(order='C')
# medfilt(mdata)

x = xdata

# function for genetic algorithm to minimize (sum of squared error)
def sumOfSquaredError(parameterTuple):
    warnings.filterwarnings("ignore") # do not print warnings by genetic algorithm
    val = func(xdata, *parameterTuple)
    return np.sum((mdata - val) ** 2.0)


def generate_Initial_Parameters():
    # min and max used for bounds
    maxX = max(xdata)
    minX = min(xdata)
    maxY = max(mdata)
    minY = min(mdata)
    minBound = min(minX, minY)
    maxBound = max(maxX, maxY)
    parameterBounds = []
    parameterBounds.append([minBound, maxBound]) # search bounds for r1
    parameterBounds.append([minBound, maxBound]) # search bounds for r2
    parameterBounds.append([minBound, maxBound]) # search bounds for r3
    parameterBounds.append([minBound, maxBound]) # search bounds for l
    # parameterBounds.append([minBound, maxBound]) # search bounds for c

    # "seed" the numpy random number generator for repeatable results
    result = differential_evolution(sumOfSquaredError, parameterBounds, seed=3)
    return result.x

# by default, differential_evolution completes by calling curve_fit() using parameter bounds
par_init = generate_Initial_Parameters()
print (par_init)

# initial values
# par_init = np.array([0.1,10,1,10])

best, cov, info, message, ier = leastsq(residual_two_functions,
                                        par_init, args=(x, mdata, ndata),
                                        full_output=True)

print(" Best-Fit Parameters: ",  best)
# print(info)
m_fit = func(xdata, *best)
n_fit = gunc(xdata, *best)


plt.figure()
plt.plot(mdata, ndata, 'bD', label='data')
plt.plot(m_fit,n_fit, 'r-',label='fitted curve')
plt.xlabel('mdata')
plt.ylabel('ndata')
plt.legend()

plt.figure()
plt.plot(xdata, mdata, 'bD', label='mdata')
plt.plot(xdata, m_fit, 'r-',label='fitted curve')
plt.xlabel('xdata')
plt.ylabel('mdata')
plt.legend()

plt.figure()
plt.plot(xdata, ndata, 'bD', label='ndata')
plt.plot(xdata, n_fit, 'r-',label='fitted curve')
plt.xlabel('xdata')
plt.ylabel('ndata')
plt.legend()
plt.show()

результаты подгонки далеки от совершенства.как оптимизировать функцию?спасибо заранее!

рисунок: данные против mdata данные против xdata данные против xdata

канал данных: testdata_link

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