This study presents a simplified analytical expression to estimate the vertical normal stress \(\sigma _z\) beneath the center of an elliptically loaded area, derived by extending classical Boussinesq theory. The formulation employs a polar coordinate transformation and a binomial series expansion in terms of eccentricity e, truncated at the sixth-order term. It reduces exactly to the classical circular case when \(e = 0\) , ensuring theoretical consistency. Validation against reference solutions shows that the expression maintains high accuracy for eccentricities up to \(e = 0.9\) and depth ratios up to \(z/b = 2\) , with deviations remaining below 4%. For \(z/b = 2.00\) , the sixth-order approximation yields only 3.80% error, while including terms up to \(e^{10}\) reduces errors to below 0.5% for \(z/b = 2.55\) . A detailed parametric study highlights the nonlinear influence of eccentricity, footing width, and depth on vertical stress behavior. For instance, increasing e from 0.2 to 0.9 increases stress from approximately 43kPa to 59kPa at fixed depth and semi-axis, while stress sharply decreases with increasing depth under constant geometry. The proposed expression offers a computationally efficient and accurate alternative to elliptic-integral-based approaches, with practical relevance in geotechnical and structural engineering applications involving non-circular surface loads.