Least Squares Data Fitting
Import the libraries
import numpy as np
from scipy.optimize import least_squares
Get the data points: $(u_i, p_i)$
u_i = np.array([0, np.pi /4.0, np.pi / 2.0])
p_i = np.array([5.0, 4.0, 6.0])
Define a function, $p(u)$, parameterised by $a$ and $b$, which fits the data points $u_i$. Let
$$ p(u) = a \cos ( u ) + b u $$def p(beta, u):
return beta[0] * np.cos(u) + beta[1] * u
Compute the function to minimize, $\beta = (a,b)^T$
def fun(beta):
return p(beta, u_i) - p_i
The scipy least square routine solves the problem iteratively from an initial guess, so let
beta0 = np.array([1.0, 1.0])
Compute the least squares solution. Note that the output of the cost function provide by the scipy function is half the sum of the square differences, see here, so
res = least_squares(fun, beta0)
print("beta:", res.x)
print("cost:", 2.0 * res.cost)
print("fun:", res.fun)
print("message:", res.message)
print("success:", res.success)
The analytical results from solving the linear system
a = (25.0 + 2.0 * np.sqrt(2.0)) / 7.0
b = (88.0 - 10.0 * np.sqrt(2.0)) / (7.0*np.pi)
print("a:", (25.0 + 2.0*np.sqrt(2.0)) / 7.0)
print("b:", (88.0 - 10.0 * np.sqrt(2.0)) / (7.0*np.pi) )
That is $a = \dfrac{25 + 2 \sqrt{2} }{7}$ and $b= \dfrac{88 - 10 \sqrt{2}} {7 \pi}$
The error, $E$, which is minimized is:
E = (p_i[0] - p([a,b], u_i)[0])**2 + (p_i[1] - p([a,b], u_i)[1])**2 + (p_i[2] - p([a,b], u_i)[2])**2
print(E)
that is:
$$ \begin{align*} E = \left( \dfrac{10 - 2\sqrt{2}}{7} \right)^2 + \left( \dfrac{4-10\sqrt{2}}{7} \right)^2 + \left( \dfrac{5 \sqrt{2} - 2}{7} \right)^2 \\ & = \dfrac{2}{7} \left( 27 - 10 \sqrt{2} \right) & = 3.673675536076871 $$