Catastrophe theory is an area of non-linear science and a mathematical theory that studies the classification of equilibrium behavior of dynamical systems in the neighborhood of singularities (Thom, 1975; Arnold, 1988). From the earlier works of its founders, it was proved that at these critical points, a dynamical system could be described by seven elementary catastrophe models (Castrigiano and Hayes, 2004), from which cusp catastrophe appears the most popular, especially in social and behavioral science applications. The cusp model assumes a dissipating or potential-minimizing system that obeys Eq. 1, which expresses the first derivative of a potential function, U, with respect to the outcome, y:

Setting the first derivative equal to zero resulted in an equilibrium function that could be represented by the three-dimensional surface (Figure 1). The two independent variables that describe the state variable (dependent measure) are called control factors and are named asymmetry (a) and bifurcation (b). The assumed potential-optimization process is theoretically compatible with human systems since the psychological construct involved in decision-making and behaviors always takes part in the self-regulation process that aims to optimize some functions related to cognitive conflict and adaptation.

Catastrophe theory adheres to epistemology complexity theory, a meta-theoretical framework known as Complex Adaptive system (CAS), which ponders the dynamical processes evolving via self-organization mechanisms (Prigogine, 1961; Nicolis and Prigogine, 1977).

Contemporary stochastic catastrophe theory enables researchers to model empirical data and detect potential non-linearities (Cobb, 1978; Guastello, 2002). These are described on the cusp response surface (Figure 1) and include bimodality, hysteresis, inaccessibility area, divergence, bifurcation, and sudden jumps (Gilmore, 1993). The three-dimensional cusp response surface suggests that at the back area of the surface, changes are smooth and linear relationships hold between the state variable and the asymmetry, (a) while in the front area, the surface folds and two distinct regions appear. These might correspond to qualitatively different states or behavioral modes, called attractors in the language of complexity theory. Here, the probability density function of the state variable becomes bimodal, and in the area between the two attractors, the area of inaccessibility, the behavior is unlikely to befall. This defines a bifurcation effect, whereas changes can occur merely as transitions between the two behavioral attractors. From the front of the surface, a sigma-like line denotes that a hysteresis effect is observed, which also explicates the emergent discontinuity. These features are dynamic in nature, occurring when the bifurcation variable, (b) drives beyond a threshold value.

It is important to accentuate that detecting bifurcation effects in empirical data justifies fostering the CAS framework since these phenomena can only occur and characterize complex non-linear systems (Nicolis and Nicolis, 2007).

The present study is exploratory research investigating teachers’ attitudes toward STEM education, which are decisive for exhibiting positive behaviors in the school context and further the successful implementation of STEM in contemporary educational settings.

There are a number of latent constructs, such as cognitive and affective conditions, that, as descendant variables, could influence and shape attitudes. These, along with self-efficacy and commitment, are teachers’ readiness dimensions (Abdullah et al., 2017; Papagiannopoulou et al., 2023). These independent variables are linearly related to attitudes, and high self-efficacy, cognitive preparedness, and commitment most likely determine highly positive attitudes. Similarly, affective conditions should be in line with positive dispositions, which analogously affect attitudes. However, emotions have an idiosyncratic character and a more crucial role (Cacioppo and Gardner, 1999; Sy et al., 2005; van Kleef, 2009). Emotions are the most sensitive factor influenced by internal and external fluctuations and are very involved in the human self-regulation process, a fact that renders them a decisive determinant of behavior (Folkman et al., 1986; Keltner and Haidt, 1999; van Kleef et al., 2004, 2010). This distinct and rather anticipated role is associated with non-linear phenomena hypotheses that have been supported by empirical evidence (Stamovlasis and Sideridis, 2014; Antoniou et al., 2022).

For instance, if the cognitive preparedness of teachers is high, it is expected to influence attitudes positively, but concurrent high levels of negative affective state acting as a strong moderator factor might cause unfavorable results. Emotions, which are involved with the course of self-regulation, might provide unfavorable feedback, conflicts, and clashes to teachers’ cognitive system, and putting stress into the decision-making process might cause interruptions and contribute to inhibitory processes. Teachers might possess high levels of readiness in terms of cognitive conditions, self-efficacy, and commitment; however, personal conceptions, worries, or concerns about feasibility or efficacy issues can be reflected in negative affective states. This can occur independently of the level of readiness, and, from a measurement perspective, the cases where two identical readiness conditions (scores in the independent variables) are associated with different scores or even polarly opposite results in attitudes might come true. This is seemingly a paradox in traditional thought and cannot be captured by linear models. Theoretically, this implies that a threshold should be present in control factors, beyond which rapid changes occur in the dependent variable, which is observed as discontinuity. The proper approach here is to apply catastrophe theory models (Guastello, 2002) and foster complexity and non-linear dynamics as the theoretical framework to explain such potential effects (Nicolis and Nicolis, 2007).

Thus, it is reasonable to hypothesize that the affective state in the present context will also be predictive of attitudinal factors and be better understood within the complexity and non-linear framework. The following hypotheses were posited:

The above hypotheses are interdependent and tested via the empirical data’s exploratory analyses.

The participants were 494 educators, 352 working in primary and 142 in secondary education [mean age 44.9 ± 9.5 years; 21.5% men and 78.5% women; mean years of teaching experience 17.2 ± 10.5]. The data collection was carried out via a web-based form, where the teachers completed the questionnaires voluntarily and anonymously, while the procedure is regarded as an opportunity sampling. The subjects were approached through social media, and an accompanying letter explained the purpose of this survey and its scientific drive. The study was conducted in accordance with the Declaration of Helsinki. It has been approved following the Ethics and Deontology Committee of the Aristotle University of Thessaloniki, Greece protocol (no. 36728/2023).

In Cusp 1 (Table 2), the difference (Cogn – Affe) acts as asymmetry [b = 0.9441972, p < 0.001], and the sum (Cogn + Affe) acts as the bifurcation factor [b = 1.1502093, p < 0.001]. The chi-square test of the linear vs. cusp model gives χ² = 413.2, df = 2, p < 0.001, and the model fit statistics AIC, AICc, and BIC favor the cusp catastrophe model. The values for the cusp [AIC = 963.9392, AICc = 964.2361, and BIC = 997.5595] are minimum compared to the linear [AIC = 1369.4596, AICc = 1369.5414, and BIC = 1386.2698] and logistic [AIC = 1358.4593, AICc = 1358.6897, and BIC = 1387.8770], respectively. The value of pseudo-R² is larger for the cusp; however, this index is not reliable and not interpreted as the % variance explained.

Figure 2 depicts a visual display of the lower part of the cusp surface, including the bifurcation area as a shaded region. Since at least 10% of the observations fall within this area, it is evident that the analysis supports the cusp model (van der Maas et al., 2003). The size of the dots that are empirical data is a function of the observed bivariate density of the bifurcation factor’s values at that point’s location, and the color denotes their position relative to the distance between the two attractors; that is, the darker color denotes observations that are on or closer to the upper attractor and the lighter color denotes observations that are on or closer to the lower attractor. Figure 3 portrays the three-dimensional cusp surface of attitudes toward STEM as a function of the two control factors and encompasses additional visual support for the bifurcation effect, demonstrating that the observations close to the folded surface are not located within the inaccessibility area but on the upper and the lower parts of the surface, that is the attractors’ areas.

In Cusp 2 (Table 3), the difference (Com – Affe) acts as asymmetry [b = 0.97431, p < 0.001] and the sum (Com + Affe) acts as the bifurcation factor [b = 1.00964, p < 0.001]. The chi-square test of the linear vs. cusp model gives χ² = 373.1, df = 2, p < 0.001, and the model fit statistics AIC, AICc, and BIC favor the cusp catastrophe model. The values for the cusp [AIC = 991.6047, AICc = 991.7772, and BIC = 1016.820] are minimum compared to the linear [AIC = 1360.7172, AICc = 1360.7990, and BIC = 1377.527] and logistic [AIC = 1343.2242, AICc = 1343.3471, and BIC = 1364.237], respectively. The value of pseudo-R² is larger for the cusp; however, this index is not reliable and not interpreted as the % variance explained.

Figure 4 depicts a visual display of the lower part of the cusp surface, including the bifurcation area (shaded region). The dots shown on the surface are allocated between the two attractors, providing support for the bifurcation effect. Moreover, Figure 5 which portrays the three-dimensional cusp surface of attitudes toward STEM as a function of the two control factors, encompasses additional visual support, demonstrating that the observations are not located within the area of inaccessibility but on the upper and the lower surface (attractors’ areas).

In Cusp 3 (Table 4), the difference (SEff – Affe) acts as asymmetry [b = 1.19772, p < 0.001] and the sum (SEff + Affe) acts as the bifurcation factor [b = 1.12337, p < 0.001]. The chi-square test of the linear vs. cusp model gives χ² = 386.6, df = 2, p < 0.001, and the model fit statistics AIC, AICc, and BIC favor the cusp catastrophe model. The values for the cusp [AIC = 950.6089, AICc = 950.9058, and BIC = 984.2292] are minimum compared to the linear [AIC = 1329.5382, AICc = 1329.6200, and BIC = 1346.3484] and logistic [AIC = 1317.4811, AICc = 1317.7116, and BIC = 1346.8989], respectively. The value of pseudo-R² is larger for the cusp; however, this index is not reliable and not interpreted as the % variance explained.

The corresponding visual displays for Cusp 3 are Figures 6, 7. The first depicts the lower part of the cusp surface, including the bifurcation area, and the observation is located on the upper and lower part of the surface (dark and lighter colors, respectively). Figure 7 portrays the three-dimensional cusp surface of attitudes toward STEM as a function of the two control factors and demonstrates the bifurcation effect with the observational data displayed close to or on the two attractors’ areas.

The general interpretation of a cusp model suggests that at the lower values, the bifurcation changes are smooth and linear, and beyond a threshold value, the system enters the area of unpredictability. In this area, cases/individuals with the same values (scores) on control factors can be found either at the lower or upper attractor regions. This brings unpredictability to the system and implies that potential changes occur merely as sudden jumps between the two behavioral attractors. In other words, small random fluctuations in the independent variables can cause large changes in the dependent measure, that is, transitions between positive (high) and negative (low) attitudinal values.

Returning to the specific cusp models, employing, e.g., Cusp 1 with (Cogn – Affe) as asymmetry and (Cogn + Affe) as bifurcation, given that the roles of Cogn and Affe are visualized as two opposing forces when the difference (Cogn – Affe) is small, changes are expected to be smooth and linear compared to the effect of their sum (Cogn + Affe). When the net moderating effect of combined high Cogn and high Affe suddenly increases, in that scale, a threshold value is likely to exist, beyond which small fluctuations can induce abrupt changes, resulting in a bifurcation process. The cusp structure revealed in the empirical data is interpreted merely if complexity theory is fostered, and it is due to the systems’ dynamics and the underlying self-organization processes. Thus, the present analysis supports the proposed control variables’ effective role and establishes catastrophe phenomena on teachers’ attitudes regarding the implementation of STEM education. The primary finding is that emotions and affective conditions, which are involved in self-regulation, are not just ordinary linear predictors of attitudes, but their influence can induce discontinuous changes.

Teachers’ attitudes are a crucial determining factor for the successful implementation of STEM education, which has attracted researchers’ attention. The individual differences related to the readiness for STEM have been associated with attitudes, which finally influence the educational process and teachers’ responsibilities in the classroom for effective teaching (Ibrahim et al., 2017). Linear statistical models have provided such evidence, but the present research reveals the unique role of affective conditions, which is captured merely within the complexity and non-linear dynamics framework. Emotions could interact with other variables that influence attitudes, such as cognitive preparedness, self-efficacy, and commitment, in an opposing/moderating way, and their net effect induces non-linear phenomena. The three cusp models presented in the previous sections capture these interactions that are taking place in the relevant self-regulation processes, providing an enhanced theoretical understanding. As mentioned in the rationale section, the idiosyncratic character of emotions justifies the findings, which confirm the association of affective conditions with discontinued changes and unpredictability (Stamovlasis and Sideridis, 2014; Antoniou et al., 2022).

This catastrophe theory application is the first reported in the domain of STEM educational research, and as such, the discoveries have a series of implications for STEM research, theory, and practice. Predominately, the results enrich the theoretical premise being developed on interdisciplinary terrain, where basic psychological research elucidates the decisive role of affective conditions for teachers’ STEM implementation. It is apparent that the findings challenge the traditional epistemological framework based on linear and mechanistic assumptions. Given that bifurcations merely characterize complex adaptive systems (Nicolis and Nicolis, 2007), their detection designates that the underlying system is, ontologically, a non-linear dynamical system, and it should be investigated as such (Stamovlasis, 2010). This epistemological implication, closely related to theory development, suggests that research should not be restricted to traditional modeling but should consider co-examining alternative non-linear models with better explanatory power.

On the other hand, in practice, stakeholders, curriculum developers, and educational leaders are informed that the affective dimension in STEM implementation is crucial to the extent that it could change attitudes and behavior radically, so any designed intervention or educational reform should not ignore it.

Despite its important findings, this study also has some limitations. One is due to cross-sectional data collected via a self-reported survey and opportunity sampling. In addition, it is rather correlational, and thus, conclusions related to causality should be drawn cautiously. Additional limitations might be imposed by validity issues originated, as in usual psychometric studies, from a small number of items measuring latent constructs. Since it is the first study reported on attitudes regarding STEM implementation, the findings should be replicated with additional samples to establish firm conclusions. Of course, the study is restricted to particular variables, while other factors might be involved in the strengthening or reduction of teachers’ attitudes, e.g., environmental support, parental involvement, or burnout (Sideridis and Alghamdi, 2023). The non-linear approach and the findings are encouraging to extend the research with other individual differences, and further work is necessary to understand the role of the affective condition under different circumstances, which might be interesting to teaching effectiveness and educational policies regarding curricula and teacher training for STEM. Moreover, the affective dimension also concerns students, the satisfaction of which is sought after in teaching research, and well-documented investigations (e.g., Chang et al., 2022; Zhao et al., 2022) might be extended by investigating analogous non-linear effects.

In summary, the conclusions to be highlighted are as follows: first, the new knowledge emerged from cusp analysis, i.e., that the affective conditions are influential factors of teachers’ attitudes toward STEM but are associated with nonlinear changes, is a property that should be taken into consideration by researchers and by policies makers when design educational reforms on this matter. In addition, the complexity and non-linear dynamics proved to be an applicable meta-theoretical framework, and new avenues of investigation are open in the area of STEM education.