title: Linear fit including all errors with NumPy/SciPy
body: <p>I have a lot of x-y data points with errors on y that I need to fit non-linear functions to. Those functions can be linear in some cases, but are more usually exponential decay, gauss curves and so on. SciPy supports this kind of fitting with <code>scipy.optimize.curve_fit</code>, and I can also specify the weight of each point. This gives me weighted non-linear fitting which is great. From the results, I can extract the parameters and their respective errors.</p>
  
  <p>There is just one caveat: The errors are only used as weights, but not included in the error. If I double the errors on all of my data points, I would expect that the uncertainty of the result increases as well. So I built a test case (<a href="https://github.com/martin-ueding/test-fit-error/blob/master/fit-error">source code</a>) to test this.</p>
  
  <p>Fit with <code>scipy.optimize.curve_fit</code> gives me:</p>
  
  <pre><code>Parameters: [ 1.99900756  2.99695535]
  Errors:     [ 0.00424833  0.00943236]
  </code></pre>
  
  <p>Same but with <code>2 * y_err</code>:</p>
  
  <pre><code>Parameters: [ 1.99900756  2.99695535]
  Errors:     [ 0.00424833  0.00943236]
  </code></pre>
  
  <p>Same but with 2 * y_err:</p>
  
  <p>So you can see that the values are identical. This tells me that the algorithm does not take those into account, but I think the values should be different.</p>
  
  <p>I read about another fit method here as well, so I tried to fit with <code>scipy.odr</code> as well:</p>
  
  <pre><code>Beta: [ 2.00538124  2.95000413]
  Beta Std Error: [ 0.00652719  0.03870884]
  </code></pre>
  
  <p>Same but with <code>20 * y_err</code>:</p>
  
  <pre><code>Beta: [ 2.00517894  2.9489472 ]
  Beta Std Error: [ 0.00642428  0.03647149]
  </code></pre>
  
  <p>The values are slightly different, but I do think that this accounts for the increase in the error at all. I think that this is just rounding errors or a little different weighting.</p>
  
  <p>Is there some package that allows me to fit the data and get the actual errors? I have the formulas here in a book, but I do not want to implement this myself if I do not have to.</p>
  
  <hr>
  
  <p>I have now read about <code>linfit.py</code> in another question. This handles what I have in mind quite well. It supports both modes, and the first one is what I need.</p>
  
  <pre><code>Fit with linfit:
  Parameters: [ 2.02600849  2.91759066]
  Errors:     [ 0.00772283  0.04449971]
  
  Same but with 20 * y_err:
  Parameters: [ 2.02600849  2.91759066]
  Errors:     [ 0.15445662  0.88999413]
  
  Fit with linfit(relsigma=True):
  Parameters: [ 2.02600849  2.91759066]
  Errors:     [ 0.00622595  0.03587451]
  
  Same but with 20 * y_err:
  Parameters: [ 2.02600849  2.91759066]
  Errors:     [ 0.00622595  0.03587451]
  </code></pre>
  
  <hr>
  
  <p>Should I answer my question or just close/delete it now?</p>