-
Notifications
You must be signed in to change notification settings - Fork 0
/
least_linear_squares_no_uncertainties.py
70 lines (50 loc) · 1.54 KB
/
least_linear_squares_no_uncertainties.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
57
58
59
60
61
62
63
64
65
66
67
68
69
import matplotlib.pyplot as plt
import numpy as np
#let's generate data for the fit
x = np.linspace(0,10,10)
y = np.linspace(0,10,10)
#plotting this data
plt.scatter(x,y, label = "testing", color = 'blue')
# Add axis labels
plt.xlabel('x-axis')
plt.ylabel('y-axis')
# Add a title
plt.title('Graph of test')
# Add a legend
plt.legend()
# Show the plot
plt.grid(True)
plt.show()
# ok, now to fit to this data set using linear least squares regression method......
x_sum = sum(x)
y_sum = sum(y)
x_squared_sum = sum(x**2)
x_y_sum = sum(x*y)
N = len(x)
slope = (N*x_y_sum - x_sum*y_sum)/(N*x_squared_sum -(x_sum)**2)
y_intercept = (y_sum - slope*x_sum)/(N)
#Finding uncertainties in fit constants
y_predicted = slope*x + y_intercept
sum_residuals = np.sum((y - y_predicted)**2)
sum_residuals = sum_residuals/(N - 2)
cool_x = np.sum((x - np.mean(x))**2)
slope_uncertainty = np.sqrt(sum_residuals/cool_x)
weirdo = np.sum((y - (slope*x + y_intercept))**2)
sigma_y = np.sqrt(weirdo/(N - 2))
y_intercept_uncertainty = weirdo*np.sqrt(x_squared_sum/(N*x_squared_sum - (x_sum)**2))
print("Slope:", slope, "+/-", slope_uncertainty)
print("Y-Intercept:", y_intercept, "+/-", y_intercept_uncertainty)
# Plotting the fitted line
y_fit = slope * x + y_intercept
plt.scatter(x, y, label='testing', color='blue')
plt.plot(x, y_fit, label=f'Fit: y = {slope:.2f}x + {y_intercept:.2f}', color='red')
# Add axis labels
plt.xlabel('x-axis')
plt.ylabel('y-axis')
# Add a title
plt.title('Graph of Linear Fit')
# Add a legend
plt.legend()
# Show the plot
plt.grid(True)
plt.show()