-
Notifications
You must be signed in to change notification settings - Fork 0
/
least_square.py
56 lines (42 loc) · 1.35 KB
/
least_square.py
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
#!/usr/bin/env python3
import numpy as np
import scipy.linalg
def least_square(data, model):
"""
Least square fit of a 2D array
args:
data: nxn array
model: nxn array
returns:
[LP0 LP1]: LP0 + LP1 * model
None: if input is not correct
"""
"""
if len(data.shape) != 2 or data.shape != model.shape or data.shape[0] != data.shape[1]:
print("The matrix cannot be used for fitting!")
return None
"""
# Make sure the input are arrays without single dimension
data = np.squeeze(np.array(data))
model = np.squeeze(np.array(model))
n = data.size
aa = np.array([[n, np.sum(model)], [np.sum(model), np.sum(model**2)]])
b = np.array([np.sum(data), np.sum(data*model)])
x = scipy.linalg.inv(aa)@b
return x
def least_square_1d(data, model):
return least_square(data, model)
if __name__ == "__main__":
import scipy.io
tt = scipy.io.loadmat('data/least_square_in.mat')
data = tt['data']
model = tt['model']
paras = least_square(data, model) # MATLAB output: [0.8361, 0.0449]
print(paras)
print('-----------')
data1 = np.copy(data)
model1 = np.copy(model)
data1 = data1.reshape(data1.size)
model1 = model1.reshape(model1.size)
paras = least_square_1d(data1, model1) # MATLAB output: [0.8361, 0.0449]
print(paras)