ABSTRACT In the pursuit of a sustainable future with limited resources and growing environmental concerns, adhesive joining stands out as a pivotal technology for the advancement of lightweight structures. This study aims to investigate the influence of various variables on the critical strain energy release rate in mode II (GIIc) for End-Notched Flexure (ENF) bonded composite joints, employing a new methodology based on Design of Experiments (DoE) approach. Two-dimensional finite element models of the ENF bonded composite joints are developed using commercial software, with all numerical models generated via Python™ scripts linked with AbaqusⓇ. Eleven design parameters related to the specimen geometry and material properties are systematically evaluated to determine the influence of each variable on GIIc. Numerical force–displacement curves are obtained, followed by the application of the Compliance-Based Beam Method (CBBM) to estimate the critical fracture energy in mode II, represented by a numerical envelope. The influence of each parameter is assessed using the main effect (ME) metric. The top five most influential variables affecting GIIc and CBBM are identified as the adhesive’s critical strain energy release rate in mode II (GIIc), effective laminate longitudinal elastic modulus (), adhesive thickness (tA), adherent thickness (h), and pre-crack length (). Furthermore, the study highlights the epistemic uncertainty associated with the geometric variables. Despite the limitations inherent in the computational model, the Plackett-Burman method consistently proves effective in conducting sensitivity analysis of the variables in the ENF test. These findings demonstrate promising prospects for the application of this procedure in the design of composite joint structures.