-1

I'm currently trying to compute an exact jacobian for scipy's optimize.minimize function. For the fitting routine, it's a mixture of Gaussian peaks, and I am optimizing the number of Gaussians as well as the peak parameters. For my error function, I am using the residual sum of squares to guide the fit.

My question is on whether I am correctly computing the jacobian to guide the gradient of the minimization routine.

I'm coding in python, so this is very pythonesque because my Latex skills are just a bit past non-existent. If it's too much trouble, I can reformat to Latex.

My error function is:

sum((fitted_data-raw_data)^2)

For my jacobian, i calculated the derivative (using Wolfram Alpha) of each parameter:

amp_jac = []
mu_jac = []
sigma_jac = []
for x,y in zip(xdata, ydata):
   gauss = gaussian(x, params)
   amp_jac.append(2*gauss-2*y*np.exp(-(x-mu)**2/(2*sigma**2)))
   mu_jac.append(2*gauss*(1-mu)*(gauss-y)/(sigma**2))
   sigma_jac.append(2*gauss*((1-mu)**2)*(gauss-y_val)/(sigma**3)))
jacobian = [sum(amp_jac), sum(mu_jac), sum(sigma_jac)]

My jacobian function is quite large compared to the numerical jacobian I can compute with scipy. The reason I am computing my own is for for 1) speed and 2) because the delta in the jacobian function is not always able to fit the peaks I am trying to fit.

So in short -- am I correct in taking the derivative of my residual sum of squares as my jacobian for a minimization routine? Or should I be taking the jacobian of my Gaussian?

Chrismit
  • 99

1 Answers1

0

I'm an idiot. I put in a mock value of 1 for x in Wolfram Alpha and forgot to substitute it back in. I figure this might be helpful to others who want to compute a jacobian with scipy's minimzation function. If the community disagrees, feel free to delete this.

Chrismit
  • 99