I am trying to derive the coefficients $a$ and $b$ for a least squares curve of the form $y=\frac{a}{x^2+b}$
To start I have written $$S=\sum_i(y_i-y)^2=\sum_i\left(y_i-\frac{a}{x_i^2+b}\right)^2=\sum_i\left(y_i^2-\frac{2a y_i}{x_i^2+b}+\frac{a^2}{(x_i^2+b)^2}\right) \\S=\sum_iy_i^2-2a\sum_i\frac{y_i}{x_i^2+b}+a^2\sum_i\frac{1}{(x_i^2+b)^2}$$
Now the way I have seen it done for the line $y=a+bx$ is they minimize $S$ with respect to both $a$ and $b$, applying that for this curve gives the two constraint equations:
$$0=a\sum_i\frac{1}{(x_i^2+b)^2}-\sum_i\frac{y_i}{x_i^2+b}\\ 0=-a\sum_i\frac{1}{(x_i^2+b)^3}+\sum_i\frac{y_i}{x_i^2+b}$$
Summing these together we find $$0=a\left(\sum_i\frac{1}{(x_i^2+b)^2}-\sum_i\frac{1}{(x_i^2+b)^3}\right)$$ Since $a\neq 0$ we must have $$\sum_i\frac{1}{(x_i^2+b)^2}=\sum_i\frac{1}{(x_i^2+b)^3}$$ Which seems to imply $$(x_i^2+b)^2=(x_i^2+b)^3$$ and unfortunately there doesnt appear to be a constant $b$ that solves this equation. Did I make a mistake somewhere or is this the result?
Starting with a rather silly sentence : since the model is nonlinear with respect to $b$, by the end, you will need a nonlinear regression and, as usual, you will need reasonable starting guesses for the parameters.
You can get these estimates in a preliminary step $$y=\frac{a}{x^2+b} \implies \frac 1y=\frac {x^2} a+\frac ba$$ So, defining for each $(x_i,y_i)$ data points, new variables $t_i=x_i^2$ and $z_i=\frac 1{y_i}$, you have $$z=\alpha\, t+\beta$$ and, from the corresponding linear regression, you will get, as estimates $a=\frac 1 \alpha$ and $b=\frac \beta \alpha$.
If you do not want to use nonlinear regression, you can reduce the problem to a nasty equation in $b$.
Using Rohan's expressions, let us define $$S_1=\sum \frac{y_i}{x_i^2+b} \qquad S_2=\sum \frac{1}{(x_i^2+b)^2}$$ $$S_3=\sum \frac{y_i}{(x_i^2+b)^2} \qquad S_4= \sum \frac{1}{(x_i^2+b)^3}$$ Eliminating $a$, you then need to to find the zero of function $$f(b)=S_2\, S_3-S_1 S_4$$ Since you have the estimate, you could use it as the starting point of Newton method with analytical derivatives since $$f'(b)=S_2S'_3+S'_2S_3-S_1S'_4-S'_1S_4$$ where will appear other summations (I let you the task of writing them - they are not difficult). The iterative process will converge very fast. So $b$ and then $a$.
Edit
For illustration purposes, consider the following data set $$\left( \begin{array}{cc} x & y \\ 0 & 46 \\ 1 & 46 \\ 2 & 45 \\ 3 & 43 \\ 4 & 41 \\ 5 & 38 \\ 6 & 36 \\ 7 & 33 \\ 8 & 30 \\ 9 & 28 \end{array} \right)$$ The preliminary step leads to $z=0.000176905 t+0.0216373$ from which the estimates $a=5652.76$ and $b=122.310$.
Starting with this guess, the iterates of Newton method will be $$\left( \begin{array}{cc} n & b_n \\ 0 & 122.31000 \\ 1 & 121.84307 \\ 2 & 121.85174 \end{array} \right)$$ from which $a=5635.22$. A nonlinear regression would have be providing exactly the same values for the parameters.