So, suppose we have an equation of the form \(y = m x + b\), and we observe a bunch of data points \((x_i, y_i)\) that we think should occur on this line, except for errors in \(y _ i\). \[\begin{align} m u_1 + c &= v_1 \\ m u_2 + c &= v_1 \\ \vdots \\ m u_i + c &= v_i \\ \vdots \end{align}\] The parameters \(m\) and \(b\) are our unknowns that we'd like to solve for, given our data points, but we have way more than 2 equations for our 2 unknowns, so we'll use linear leastsquares. If we rearrange our equations into matrix form, \[\begin{gather} \begin{bmatrix} u_1 & 1 \\ u_2 & 1 \\ \vdots & \vdots \\ u_i & 1 \\ \vdots & \vdots \end{bmatrix} \begin{bmatrix} m \\ c \end{bmatrix} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_i \\ \vdots \end{bmatrix} \end{gather}\] So, then \[\begin{gather} A = \begin{bmatrix} u_1 & 1 \\ u_2 & 1 \\ \vdots & \vdots \\ u_i & 1 \\ \vdots & \vdots \end{bmatrix}, \quad x = \begin{bmatrix} m \\ c \end{bmatrix}, \quad b = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_i \\ \vdots \end{bmatrix}. \end{gather}\] We can now calculate \(A^T A\) and \(A^T b\) and use
x = scipy.linalg.solve(A.T.dot(A), A.T.dot(b))
to solve the normal equation \(A^T A x = A^T b\). (Expert note: It is often recommended to use a QR algorithm or SVD decomposition to find linear leastsquares solutions rather than solving the normal equations directly by elimination with partial pivoting.)
Interestingly, we can use the same methods to fit exponential and power formulas if we first apply a logtransformation.
\[y=kx^a \quad \text{or} \quad \log y = \log k + a \log x \quad \text{or} \quad \frac{dy}{dx} = a \frac{y}{x}\] \[y=e^{ax+c} \quad \text{or} \quad \log y = ax + c \quad\text{or} \quad \frac{dy}{dx} = a y\]
mass of adult  mass of egg 

113380  1700 
4536  165.4 
3629  94.5 
1020  34 
283  14 
3.6  0.6 
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 

In Dialogues Concerning Two New Sciences, Galileo proposed the "squarecube" law for how bones must change as the size of an animal changes. The naive theory is that if you double the mass of a person or animal, and you double mass of the skeleton correspondingly, everything will work out just fine (As suggested in Gulliver's travels, for instance.) But Galileo held that it isn't so. The argument was simple. The strength of a bone in a body will be proportional to its crosssectional area (\(A\)) of the bone. Since the bones have to support the full mass of the anmial (\(M_{body}\)), \(A \propto M_{body}\). Now, the crosssectional area is proportional to the skeletal mass \(M_{skel}\) raised to the 2/3rds power (\(A \propto M_{skel}^{2/3}\)), so with a little algebra, \(M_{skel} \propto M_{body}^{3/2}\). Galileo included a sketch in his book to illustrate  interestingly, the sketch is actually inaccurate, making the bone much thicker than the theory predicts.
While Galileo's general argument may be right  skeleton mass can not scale isometrically with body mass  the specific formula he proposed does not seem to match the data.