diff --git a/scripts/find_grade_gain.py b/scripts/find_grade_gain.py
new file mode 100644
index 0000000..f4d4a47
--- /dev/null
+++ b/scripts/find_grade_gain.py
@@ -0,0 +1,29 @@
+import matplotlib.pyplot as plt
+from scipy import stats
+import numpy as np
+
+
+def main():
+    """
+    [1] H. Schuman, E. Walsh, C. Olson, and B. Etheridge, “Effort and Reward:
+        The Assumption that College Grades Are Affected by Quantity of Study*,”
+        Social Forces, vol. 63, no. 4, pp. 945–966, June 1985.
+    """
+    # [1, p. 950]
+    hours_studied = np.array([1, 2.5, 3.5, 4.5, 5.5, 6.5])
+    gpa = np.array([2.94, 2.91, 2.97, 2.86, 3.25, 3.18])
+
+    slope, intercept, r, p, std_err = stats.linregress(hours_studied, gpa)
+
+    print(f"GPA/hour (slope) of best fit line: {slope}")
+
+    plt.plot(hours_studied, gpa, label="Plot from publication")
+    plt.plot(hours_studied, slope * hours_studied + intercept, label="Best fit")
+    plt.xlabel("Hours studied")
+    plt.ylabel("GPA")
+    plt.legend()
+    plt.show()
+
+
+if __name__ == "__main__":
+    main()