This study presents a detailed analysis of shear and peel stresses in adhesively bonded single lap joints using the Goland and Reissner analytical model. The investigation evaluates the effects of key parameters, including adhesive thickness, adhesive material, adherend material, and overlap length on stress distribution. A General Linear Model (GLM) and Analysis of Variance (ANOVA) are used to assess the significance of each factor. Results show that adhesive thickness contributes 36.55% to shear stress variation, followed by adhesive material (31.08%) and adherend material (25.83%). For peel stress, adhesive thickness accounts for 38.01% of the variation. A second-order polynomial regression model is employed to capture non-linear relationships between the input parameters and stress outcomes. The predicted shear stress of 8.676 MPa closely matches the actual value of 8.64 MPa, with a relative error of 0.42%, while the predicted peel stress of 10.9901 MPa aligns with the actual value of 11.04 MPa, with a relative error of 0.45%. The analysis highlights that thinner adhesive layers lead to higher stress concentrations, while thicker layers distribute stress more effectively. The choice of adhesive material and adherend material also significantly impacts stress levels. The study concludes that optimizing adhesive thickness, material selection, and overlap length is essential for improving the performance and reliability of adhesively bonded joints. The polynomial regression model successfully captures the non-linear stress behavior, offering a robust tool for predicting joint performance.