diff --git a/DIRECTORY.md b/DIRECTORY.md
index c5bac23ad35b..3fd519b72f9a 100644
--- a/DIRECTORY.md
+++ b/DIRECTORY.md
@@ -868,6 +868,7 @@
   * [Sum Of Digits](maths/sum_of_digits.py)
   * [Sum Of Geometric Progression](maths/sum_of_geometric_progression.py)
   * [Sum Of Harmonic Series](maths/sum_of_harmonic_series.py)
+  * [Sum Of Squares](maths/sum_of_squares.py)
   * [Sumset](maths/sumset.py)
   * [Sylvester Sequence](maths/sylvester_sequence.py)
   * [Tanh](maths/tanh.py)
diff --git a/maths/sum_of_squares.py b/maths/sum_of_squares.py
new file mode 100644
index 000000000000..8fec79ae30ca
--- /dev/null
+++ b/maths/sum_of_squares.py
@@ -0,0 +1,38 @@
+"""
+This script demonstrates the implementation of the
+sum of squares of the first n natural numbers.
+
+The function takes an integer n as input and returns the sum of squares
+from 1 to n using the formula n(n + 1)(2n + 1) / 6.
+
+This formula computes the sum efficiently
+without the need for iteration.
+
+https://www.cuemath.com/algebra/sum-of-squares/
+"""
+
+
+def sum_of_squares(num_of_terms: int) -> int:
+    """
+    Implements the sum of squares formulafor the first n natural numbers.
+
+    Parameters:
+        num_of_terms (int): A positive integer representing the limit of the series
+
+    Returns:
+        sum_squares (int): The sum of squares of the first n natural numbers.
+
+    Examples:
+    >>> sum_of_squares(5)
+    55
+
+    >>> sum_of_squares(10)
+    385
+    """
+    return num_of_terms * (num_of_terms + 1) * (2 * num_of_terms + 1) // 6
+
+
+if __name__ == "__main__":
+    import doctest
+
+    doctest.testmod()