Wie scipy Kurvenanpassung mit Fehlerbalken und Standardfehler bei der Anpassung Parameter erhalten?

Question

Ich versuche meine Datenpunkte anzupassen. Es sieht so aus, als ob die Anpassung ohne Fehler nicht so optimistisch ist, daher versuche ich jetzt, die Daten anzupassen, die die Fehler an jedem Punkt implementieren. Meine Fit-Funktion ist unter:

def fit_func(x,a,b,c): 
    return np.log10(a*x**b + c) 

dann meine Datenpunkte sind unter:

r = [ 0.00528039,0.00721161,0.00873037,0.01108928,0.01413011,0.01790143,0.02263833, 0.02886089,0.03663713,0.04659512,0.05921978,0.07540126,0.09593949, 0.12190075,0.15501736,0.19713563,0.25041524,0.3185025,0.40514023,0.51507869, 0.65489938,0.83278859,1.05865016,1.34624082] 
logf = [-1.1020581079659384, -1.3966927245616112, -1.4571368537041418, -1.5032694247562564, -1.8534775558300272, -2.2715812166948304, -2.2627690390113862, -2.5275290780299331, -3.3798813619309365, -6.0, -2.6270989211307034, -2.6549656159564918, -2.9366845162570079, -3.0955026428779604, -3.2649261507250289, -3.2837123017838366, -3.0493752067042856, -3.3133647996463229, -3.0865051494299243, -3.1347499415910169, -3.1433062918466632, -3.1747394718538979, -3.1797597345585245, -3.1913094832146616] 

Weil meine Daten in Log-Skala ist, logf, dann ist die Fehlerbalken für jeden Datenpunkt nicht symmetrisch ist. Die oberen Fehlerbalken und untere Fehlerbalken sind unter:

upper = [0.070648916083227764, 0.44346256268274886, 0.11928131794776076, 0.094260899008089094, 0.14357124858039971, 0.27236750587684311, 0.18877122991380402, 0.28707938182603066, 0.72011863806906318, 0, 0.16813325716948757, 0.13624929595316049, 0.21847915642008875, 0.25456116079315372, 0.31078368240910148, 0.23178227464741452, 0.09158189214515966, 0.14020538489677881, 0.059482730164901909, 0.051786777740678414, 0.041126467609954531, 0.034394612910981337, 0.027206248503368613, 0.021847333685597548] 
lower = [0.06074797748043137, 0.21479225959441428, 0.093479845697059583, 0.077406149968278104, 0.1077175009766278, 0.16610073183912188, 0.13114254113054535, 0.17133966123838595, 0.57498950902908286, 2.9786837094190934, 0.12090437578535695, 0.10355760401838676, 0.14467588244034646, 0.15942693835964539, 0.17929440903034921, 0.15031667827534712, 0.075592499975030591, 0.10581886912443572, 0.05230849287772843, 0.04626422871423852, 0.03756658820680725, 0.03186944137872727, 0.025601929615431285, 0.02080073540367966] 

Ich habe die Armatur als:

popt, pcov = optimize.curve_fit(fit_func, r, logf,sigma=[lower,upper]) 
logf_fit = fit_func(r,*popt) 

Aber das ist falsch, wie kann ich die Kurve implementieren aus scipy Einpassen der oberen aufzunehmen und geringere Fehler? Wie könnte ich die Anpassungsfehler der Anpassungsparameter a, b, c bekommen?

Answer 1

Sie können scipy.optimize.leastsq mit benutzerdefinierten Gewicht:

import scipy.optimize as optimize 
import numpy as np 
# redefine lists as array 
x=np.array(r) 
y=np.array(logf) 
errup=np.array(upper) 
errlow=np.array(lower) 
# error function 
def fit_func(x,a,b,c): 
    return np.log10(a*x**b + c) 
def my_error(V): 
    a,b,c=V 
    yfit=fit_func(x,a,b,c) 
    weight=np.ones_like(yfit) 
    weight[yfit>y]=errup[yfit>y] # if the fit point is above the measure, use upper weight 
    weight[yfit<=y]=errlow[yfit<=y] # else use lower weight 
    return (yfit-y)**2/weight**2 
answer=optimize.leastsq(my_error,x0=[0.0001,-1,0.0006]) 
a,b,c=answer[0] 
print(a,b,c) 

Es funktioniert, aber sehr empfindlich auf den Anfangswert, da ein Protokoll ist, die in falscher Domain (negative Zahlen) gehen und dann scheitert es. Hier finde ich a=9.14464745425e-06 b=-1.75179880756 c=0.00066720486385, die ziemlich nah an Daten ist.