title: How to determine the uncertainty of fit parameters with Python?
body: <p>I have the following data for x and y:</p>
  
  <pre><code>x   y
  1.71    0.0
  1.76    5.0
  1.81    10.0
  1.86    15.0
  1.93    20.0
  2.01    25.0
  2.09    30.0
  2.20    35.0
  2.32    40.0
  2.47    45.0
  2.65    50.0
  2.87    55.0
  3.16    60.0
  3.53    65.0
  4.02    70.0
  4.69    75.0
  5.64    80.0
  7.07    85.0
  9.35    90.0
  13.34   95.0
  21.43   100.0
  </code></pre>
  
  <p>For the above data, I am trying to fit the data in the form:</p>
  
  <p><img src="https://chart.googleapis.com/chart?cht=tx&amp;chl=y%20%3D%20a%7Bx%7D%5E%7Bn%7D%20%2B%20b" alt="formula"></p>
  
  <p>However, there are certain uncertainties associated with x and y, where x has uncertainty of 50% of x and y has a fixed uncertainty. I am trying to determine the uncertainty in the fit parameters with this <a href="https://pypi.python.org/pypi/uncertainties" rel="noreferrer">uncertainties package</a>. But, I am having issues with curve fitting with scipy optimize's curve fit function. I get the following error:</p>
  
  <blockquote>
    <p>minpack.error: Result from function call is not a proper array of
    floats.</p>
  </blockquote>
  
  <p><strong><em>How do I fix the following error and determine the uncertainty of the fit parameters (a,b and n)?</em></strong></p>
  
  <p><strong>MWE</strong></p>
  
  <pre><code>from __future__ import division
  import numpy as np
  import re
  from scipy import optimize, interpolate, spatial
  from scipy.interpolate import UnivariateSpline
  from uncertainties import unumpy
  
  
  def linear_fit(x, a, b):
      return a * x + b
  
  
  uncertainty = 0.5
  y_error = 1.2
  x = np.array([1.71, 1.76, 1.81, 1.86, 1.93, 2.01, 2.09, 2.20, 2.32, 2.47, 2.65, 2.87, 3.16, 3.53, 4.02, 4.69, 5.64, 7.07, 9.35, 13.34, 21.43])
  x_uncertainty = x * uncertainty
  x = unumpy.uarray(x, x_uncertainty)
  y = np.array([0.0, 5.0, 10.0, 15.0, 20.0, 25.0, 30.0, 35.0, 40.0, 45.0, 50.0, 55.0, 60.0, 65.0, 70.0, 75.0, 80.0, 85.0, 90.0, 95.0, 100.0])
  y = unumpy.uarray(y, y_error)
  
  
  n = np.arange(0, 5, 0.005)
  coefficient_determination_on = np.empty(shape = (len(n),))
  for j in range(len(n)):
      n_correlation = n[j]
      x_fit = 1 / ((x) ** n_correlation)
      y_fit = y
      fit_a_raw, fit_b_raw = optimize.curve_fit(linear_fit, x_fit, y_fit)[0]
      x_prediction = (fit_a_raw / ((x) ** n_correlation)) + fit_b_raw
      y_residual_squares = np.sum((x_prediction - y) ** 2)
      y_total_squares = np.sum((y - np.mean(y)) ** 2)
      coefficient_determination_on[j] = 1 - (y_residual_squares / y_total_squares)
  </code></pre>