title: fitting two lorentz by lmfit Model minimize in python
body: <p>maybe someone can help me. I spent days on it, but I couldn't get through the problem. Thanks in advance.</p>
  
  <p>I want to fit 2 lorentzians to my experimental data. I broke my equations down to simple form of two lorentzians <code>lorentz1</code> and <code>lorentz2</code> functions. Then I defined two other functions <code>L1</code> and <code>L2</code> to only multiply a constant <code>cnst</code> to them. I have all 4 parameters to fit: <code>cnst1</code>,<code>cnst2</code>,<code>tau1</code>,<code>tau2</code>.</p>
  
  <p>I use <code>lmfit</code>: Model and minimize (probably both use same method).</p>
  
  <p>The initial fit parameters are set in a way that visually it is closer to a fine fit. But minimizing using lmfit gets lost (first image below):</p>
  
  <hr>
  
  <p><img src="https://i.stack.imgur.com/eHZ1V.png" alt="unconstrained parameters"></p>
  
  <p>using these parameters:</p>
  
  <pre><code>params.add('cnst1', value=1e3 , min=1e2, max=1e5)
  params.add('cnst2', value=3e5, min=1e2, max=1e6)
  params.add('tau1', value=2e0, min=0, max=1e2)
  params.add('tau2', value=5e-3, min=0, max=10)
  </code></pre>
  
  <p>but error percentages are low:</p>
  
  <pre><code>cnst1:   117.459806 +/- 14.67188 (12.49%) (init= 1000)
  cnst2:   413.959032 +/- 44.21042 (10.68%) (init= 300000)
  tau1:    11.0343531 +/- 1.065570 (9.66%) (init= 2)
  tau2:    1.55259664 +/- 0.125853 (8.11%) (init= 0.005)
  </code></pre>
  
  <hr>
  
  <p>On the other hand, constaining the parameters to very close to initial (force to be like initial):</p>
  
  <p><img src="https://i.stack.imgur.com/DTwcR.png" alt="enter image description here"></p>
  
  <p>using parameters:</p>
  
  <pre><code>#params.add('cnst1', value=1e3 , min=0.1e3, max=1e3)
  #params.add('cnst2', value=3e5, min=1e3, max=1e6)
  #params.add('tau1', value=2e0, min=0, max=2)
  #params.add('tau2', value=5e-3, min=0, max=10)
  </code></pre>
  
  <p>fit is visually better but error values are huge:</p>
  
  <pre><code>[[Variables]]
  cnst1:   752.988629 +/- 221.3098 (29.39%) (init= 1000)
  cnst2:   3.0159e+05 +/- 3.05e+07 (10113.40%) (init= 300000)
  tau1:    1.99684317 +/- 0.600748 (30.08%) (init= 2)
  tau2:    0.00497806 +/- 0.289651 (5818.56%) (init= 0.005)
  </code></pre>
  
  <p>here is the total code:</p>
  
  <pre><code>import numpy as np
  from lmfit import Model, minimize, Parameters, report_fit
  import matplotlib.pyplot as plt
  
  x = np.array([0.02988, 0.07057,0.19365,0.4137,0.91078,1.85075,3.44353,6.39428,\
          11.99302,24.37024,52.58804,121.71927,221.53799,358.27392,464.70405])
  
  y = 1.0 / np.array([4.60362E-4,5.63559E-4,8.44538E-4,0.00138,0.00287,0.00657,0.01506,\
              0.03119,0.0584,0.09153,0.12538,0.19389,0.34391,0.68869,1.0])
  
  def lorentz1(x, tau):
      L =  tau  / ( 1 + (x*tau)**2 )  
      return(L)
  
  def lorentz2(x, tau):
      L =  tau**2  / ( 1 + (x*tau)**2 )  
      return(L)
  
  def L1(x,cnst1,tau1):
      L1 =  cnst1 * lorentz1(x,tau1)
      return (L1)
  
  def L2(x, cnst2, tau2):
      L2 =  cnst2 * lorentz2(x,tau2)
      return (L2)    
  
  def L_min(params, x, y):
      cnst1 = params['cnst1'].value
      cnst2 = params['cnst2'].value
      tau1 = params['tau1'].value
      tau2 = params['tau2'].value
  
      L_total = L1(x, cnst1, tau1) + L2(x, cnst2, tau2)
      resids = L_total - y
      return resids
  
  #params  = mod.make_params( cnst1=10e2, cnst2=3e5, tau1=2e0, tau2=0.5e-2)
  params = Parameters()
  #params.add('cnst1', value=1e3 , min=0.1e3, max=1e3)
  #params.add('cnst2', value=3e5, min=1e3, max=1e6)
  #params.add('tau1', value=2e0, min=0, max=2)
  #params.add('tau2', value=5e-3, min=0, max=10)
  
  params.add('cnst1', value=1e3 , min=1e2, max=1e5)
  params.add('cnst2', value=3e5, min=1e2, max=1e6)
  params.add('tau1', value=2e0, min=0, max=1e2)
  params.add('tau2', value=5e-3, min=0, max=10)
  
  
  #1-----Model--------------------
  mod = Model(L1) + Model(L2)
  result_mod = mod.fit(y, params, x=x)
  print('---results from lmfit.Model----')
  print(result_mod.fit_report())
  
  #2---minimize-----------
  result_min = minimize(L_min, params, args=(x,y))
  final_min = y + result_min.residual
  print('---results from lmfit.minimize----')
  report_fit(params)
  
  #-------Plot------
  plt.close('all')
  plt.loglog(x, y,'bo' , label='experimental data')
  plt.loglog(x, result_mod.init_fit, 'k--', label='initial')
  plt.loglog(x, result_mod.best_fit, 'r-', label='final')
  plt.legend()
  plt.show()
  </code></pre>