A comparative study of feature selection and feature extraction methods for financial distress identification

Feature selection

Feature selection (FS) is a dimensionality reduction approach focusing on removing redundant, noisy and non-informative features from the original feature set (Al-Tashi et al., 2020). The implementation of this approach is widespread in the financial distress and bankruptcy concepts due to its possibility to achieve knowledge of the features’ importance. Next, examples of feature selection techniques for identifying financial distress are presented. It’s worth noting that the majority of studies focused on 3–5 feature selection methods to reveal the best set of features.

Filter methods analyze each feature individually. One of the simplest examples is the correlation feature selection (CFS). The premise of CFS is that a useful feature subset comprises features that exhibit significant correlations with the class (predictive of it), but are not correlated with one another (not predictive of it). The CFS algorithm analyzes and ranks a set of feature subsets, rather than carrying out these steps for a single attribute. The CFS algorithm was analyzed in the Spanish market study (Faris et al., 2020), which aimed to improve the accuracy of the classification of bankrupt enterprises. The data was balanced using the SMOTE method, and different FS methods (information gain (IG), Relief-f, correlation Attribute Evaluator (CorrAE), Classifier Attribute Evaluator (CAE)) were applied. Results revealed that CFS achieved the best results, i.e., it achieved G-mean value of 0.66 and AUC of 0.89. The initial dataset contained 33 attributes, of which nine were selected using CFS. Another study (Séverin & Veganzones, 2021) attempted to improve the accuracy of bankruptcy prediction by using income information. This study used a similar principle as CFS, i.e., if the attributes had a correlation higher than 0.6, then one of them was removed. Unfortunately, the authors did not analyze the impact of the CFS method on the classification results.

The idea behind the Relief-F method is to assign a weight to each attribute corresponding to its ability to distinguish between different classes. It is assumed that the feature values for the two closest neighbors who belong to the same class have a probability of being the same, while the feature values for the two closest neighbors who belong to different classes have a probability of being different. Researchers (Faris et al., 2020) used the Relief-F method, which achieved a 0.63 G-mean value and a 0.91 AUC score. However, further analysis revealed that the Relief-F method was outperformed by CFS and IG methods when used in combination with ML methods. The Relief-F method was also used in another study (Kou et al., 2021) that attempted to predict SME bankruptcy cases using transaction data. The authors of this study used the ensemble technique, i.e., the attribute importance index was estimated as a combination of the outputs of Relief-f, Chi-squared, IG and gain ratio (GR) methods. Relief’s ability to find features that interact gets worse as more features are added. This happens because the way it calculates neighbors and weights becomes more random as more features are added (Urbanowicz et al., 2018).

In statistics, the Chi-square test (χ²) is used to determine whether two categorical or nominal variables are independent. We can obtain the observed count O and the expected count E from the data of two variables. Chi-square calculates the difference between the expected count E and the observed count O (Gajawada, 2019). If the observed count is close to the expected count when two features are independent, then the Chi-square value will be lower. Accordingly, a high Chi-square score suggests that the independence hypothesis is false. In simple terms, features with higher Chi-square values are associated with the response variable and might be selected for model training. The Chi-squared test was used in the previously mentioned study Kou et al. (2021) as one of the ensemble techniques. In this case, its influence on the modelling performance was not estimated. The Chi-squared test was also used in combination with classification and regression trees (Azayite & Achchab, 2018). However, the authors did not analyze the influence of feature selection on the modelling performance.

Kruskal–Wallis (K–W) and T-test are used to determine if there is a statistically significant difference between different groups. These tests were applied in a bankruptcy prediction study using earnings management information (Séverin & Veganzones, 2021). The research focused on only seven features. T-test and K–W showed that most attributes are statistically significant at the level of 1%. The K–W method identified that the ratio of Financial Expenses and Total Assets is significant only at the 5% level. A study utilized the K–W approach to examine whether the FS methods suggest employing statistically significant features (Papíková & Papík, 2022b). Results indicated that the use of FS approaches in most cases improved the AUC score, and a statistical significance of less than 5% was rarely observed. Ben Jabeur & Serret (2023) focused on the usage of the T-test. Features having statistically significant coefficients at the 1% level were selected for further research. Experiments with different ML methods revealed that the improvement of the AUC score was only observed for the T-test in combination with SVM or PLS-DA ML methods; the achieved improvement was about 3%.

The feature selection based on wrapper methods depends on a particular machine learning algorithm (Verma, 2020). It utilizes a greedy search methodology, examining all possible feature combinations. The suitability of different feature sets is based on a performance metric that changes depending on the problem type. For instance, classification tasks use accuracy, precision, recall, F1-score, etc. The most popular wrapper methods are backward/forward stepwise logistic regression (stepwise LR), stepwise discriminant analysis (stepwise DA) and generic algorithm (GA). An unspecified backward wrapper method was applied to the Italian SME’s data (Perboli & Arabnezhad, 2021). Based on the assumption that important variables should increase the accuracy by at least 1%, the method identified a subset of 15 important features out of 170 possible.

A fairly new Population Diversity Controller-Genetic Algorithm (PDC-GA), which is an improved version of the simple GA FS method, combining GA with k-mean clustering, was used to forecast the bankruptcy of Polish enterprises (Al-Milli, Hudaib & Obeid, 2021). The results indicate that the new FS algorithm improves the results of the GA FS method from 1% to 4%. After analyzing a year’s worth of data, the PDC-GA approach, when paired with the KNN ML approach, yielded a 0.9573 AUC score, whereas the GA FS approach yielded a 0.9391 AUC score. Furthermore, the utilization of the GA method in the Tsai et al. (2021) study was unsuccessful due to its ineffective performance in comparison with other dimensionality reduction methods (PCA, t-test), i.e., GA resulted in an AUC score of 0.641.

Two wrapper methods, i.e., stepwise LR and stepwise DA, were used to forecast French enterprises’ bankruptcy (Ben Jabeur, Stef & Carmona, 2022). Both wrapper methods selected less than 10 features due to the high correlation between features. The results indicated that stepwise LR and stepwise DA achieved similar results, but the XGBoost feature selection technique achieved a higher AUC score of 0.958, while the highest AUC for stepwise LR and stepwise DA was 0.935 and 0.940, respectively. In Papíková & Papík (2022a), the highest AUC value (99.95%) was achieved using the CatBoost classification method for random oversampling with a wrapper (stepwise regression) feature selection method.

Embedded methods combine the advantages of filter and wrapper methods. These methods perform feature selection and training of the algorithm in parallel. One of the embedded techniques used in financial distress or bankruptcy feature selection is the least absolute shrinkage and selection operator (LASSO). This technique eliminates uninformative predictors from the model by shrinking their coefficients to zero (Altman et al., 2022). Research using the LASSO technique can be divided into two categories: (1) studies that apply LASSO as a step in data preparation, and (2) studies that use LASSO for comparing dimensionality reduction methods. Huang, Wang & Kochenberger (2017) and Li et al. (2021) belong to the first category. In this case, the training sample consisted of up to 200 companies described by 90 features. However, the final models included only seven–eight features. Research analyzing 2,040 SMEs from Croatia using the LASSO algorithm selected 21 features out of 87 possible (Altman et al., 2022). Moreover, the LASSO technique could be used as one of the FS voting classifier techniques (Du et al., 2020). Studies belonging to the second category (Jiang et al., 2021; Papíková & Papík, 2022a; Volkov, Benoit & Van den Poel, 2017) show that the LASSO technique efficiency is lower compared with other FS techniques such as stepwise regression or random forest. However, the results obtained in Papíková & Papík (2022a) can be considered controversial because the authors claim that the impact of individual feature selection methods has no statistically significant effect on the models’ performance. Moreover, in this case, the total feature space consisted only of 27 features. An interesting analysis perspective is the stability of the LASSO method over the years (Zizi et al., 2021); the first year, LASSO identified seven important features, the next year –9, and only four of them remained the same.

Another embedded technique often used for financial distress identification is XGBoost. Ben Jabeur, Stef & Carmona (2022) and Du et al. (2020), using XGBoost method, identified nine–10 features, which cover 10–20% of the entire feature set. In comparison, tree-based, L1-based and f–classif methods selected 20–30% of the entire feature set (Du et al., 2020). However, the impact on modelling accuracy was not compared. The final model was developed using a single combined set of features. The CatBoost boosting technique was used for feature selection (Jabeur et al., 2021). CatBoost automatically ranks the most important features and uses them for forecasting. The results show that CatBoost, with its feature importance ranking, performs well as it achieves an AUC score of 0.994. However, this study did not consider other FS techniques.

Logistic regression analysis (LRA) is an embedded feature selection technique used in the financial distress prediction study (Doğan, Koçak & Atan, 2022). After applying the feature selection method, the cross-validation rate increased from 87.21% to 90.06%. However, the results of the study should be interpreted with caution because the data set is small, consisting of 172 stock market enterprises, of which 71 experienced financial difficulties. Another small data set (266 enterprises with 17 financial ratios) was used for French enterprises’ bankruptcy prediction (Ben Jabeur & Serret, 2023). This study used a Partial Least Squares Discriminant Analysis (PLS-DA) FS approach. The results were assessed over a period of one year by employing diverse classification ML techniques. The proposed feature selection method improved the AUC score results of SVM and PLS-DA ML models about 2%. After the year, the authors expanded their research by including additional FS methods, namely XGBoost, stepwise LR, stepwise DA and PLS-DA (Ben Jabeur, Stef & Carmona, 2022). In this study, the PLS-DA method achieved similar results to the wrapper methods but was inferior to the XGBoost FS method, scoring 0.3%–2.7% lower AUC.

The Filter and wrapper FS approach and their hybrid methods were examined in a study focusing on the financial distress in Taiwanese electronics industry businesses (Lin & Hsu, 2017). The hybrid FS method achieved 87.77% accuracy, the filter FS method –75.10%, and the wrapper –82.33%. The authors claim that this hybrid FS method keeps wrapper approaches’ advantages while saving computational resources.

In conclusion, the authors use feature selection techniques for the reduction of the main dimensionality problems, such as data sparsity, multiple testing, and overfitting (Kuiziniene et al., 2022). Moreover, based on the analysis of feature selection, three main directions of research can be distinguished: (1) research that uses one feature selection method as a data preprocessing step; (2) research that uses several feature selection methods as a data preprocessing step; (3) research that compares and analyzes several feature selection methods. The third direction proved the effectiveness of CFS (Faris et al., 2020), stepwise LR for a one-year period, t-test for a two-year period (Ben Jabeur & Serret, 2023), XGBoost (Ben Jabeur, Stef & Carmona, 2022), and FS voting classifier (Du et al., 2020) methods.

Feature extraction

The feature extraction (FE) approach is a process of transforming high-dimensional data into a lower-dimensional feature space (de Freitas & de Freitas, 2013). This approach is frequently used for visual information representation (Kuiziniene et al., 2022; Ye, Ji & Sun, 2013). For instance, multidimensional scaling MDS (Khoja, Chipulu & Jayasekera, 2019; Mokrišová & Horváthová, 2020) and tSNE (Zoričák et al., 2020) approaches were used for data visualization in the context of financial distress and bankruptcy. A comparative study of different feature extraction methods, including linear discriminant analysis (LDA), principal component analysis (PCA), locality preserving projections (LPP), neighborhood preserving embedding (NPE), and robust subspace learning (RSL) were presented in Wang, Liu & Pu (2019). Feature extraction approaches include linear and non–linear methods.

Linear methods reduce dimensionality by using linear functions (Ayesha, Hanif & Talib, 2020). For example, principal component analysis (PCA) is a projection-based method in which the main assumption is that the behavior of data could be represented by a significant variance (Sulistiani, Widodo & Nugraheni, 2022). It is noted that authors tend to select reduced dimensions in three manners: (1) experimentally(by selecting several reduced dimensions) (Adisa et al., 2019); (2) based on a pre-defined covariance ratio threshold (usually in the range of 60%–90% of the original information) (Jiang et al., 2021); (3) hybrid (numerical determination by experimentally choosing a percentage of the explained variance) (Sulistiani, Widodo & Nugraheni, 2022). Moreover, PCA is usually followed by a balancing technique. A more interesting application of PCA is discussed in Jiang et al. (2021). In this case, PCA was used to reduce the dimension of word vectors generated by BERT (Bidirectional Encoder Representations from Transformers). Another interesting study (Acharjya & Rathi, 2021) conducted through the use of PCA identified 15 attributes as the main features out of the total of 64 conditional attributes.

The main advantage of the non–linear methods is a better representation of the real-world data, which is more often non–linear than linear (van der Maaten, Postma & Herik, 2007). Kohonen’s self-organizing map (SOM) is an example of such methods being used in the financial distress domain (Mora García et al., 2008). However, in this research, the SOM technique was applied before splitting the data into training and testing sets and thus violating the requirement of unseen data. Another non–linear method often used in bankruptcy prediction is the autoencoder. This technique reduced the dimensionality of the data from 64 to 16 features (Soui et al., 2020). Moreover, the stacked autoencoder with softmax classification function achieved the highest AUC (87.9%) score, compared to other methods (SVM, RF, XGBoost, etc.).

The comparison of feature selection and extraction methods provided in Table 1 shows that comprehensive and high-performance feature extraction is inseparable from high computational costs.

Class

In this study, we analyze a binary classification problem, where 0 indicates “good” enterprise conditions and 1 stands for “bad” enterprise conditions. We define a “bad” situation as financial distress that is detected at least in one Lithuanian government register and meets the following conditions:

1. a bankruptcy case is filed against the enterprise;

2. the enterprise’s status changed to going bankrupt, bankrupt, under restructuring, under liquidation, liquidation being initiated, removed, liquidated, liquidation due to bankruptcy;

3. the enterprise has made announcements to the register center about bankruptcy, liquidation, restructuring, insolvency, etc.;

4. the enterprise is included in the State tax inspectorate’s lists of (1) companies temporarily exempted from submitting declarations to the STI, (2) companies that have declared temporary inactivity of the STI; (3) companies for which the STI has submitted a proposal to the State Register for deregistration in accordance with Article 2.70 of the Civil Code.

5. enterprise had no socially insured employees for the past six months.

The final dataset consists of 274,105 unique records, of which only 2.88% cases represent financial distress.

Features

The data can be divided into ten different categories, depending on the provider of data or its held information (see Table 2). For example, three data source providers are combined in the Sector’s category: (1) sectors type identified by NACE category; (2) information on sector profitability, competitiveness, etc. from the Lithuanian Statistics Department; (3) sectoral indicators calculated by combining financial statement data and NACE types. The “No.” column indicates how many features fall under each category. The data type indicates whether the feature is a primary data feature or an additional feature derived from primary data. For instance, company assets or the number of employees in December, are assigned to the primary data type, whereas ROA, gross profit, change of employees per year, etc., belong to the additional data. The data frequency is divided into three categories:

1. stable: information is constant, e.g., legal status, types of sectors;

2. depending on an event: changes when the event occurs, the number of courts, the number of changes of directors, the time elapsed since the last event, etc.;

3. periodic data (annual, quarterly, monthly): information that is updated at the indicated periodicity, e.g., financial reports, macro indicators, the number of employees.

All categorical features have been transformed into binary features by expanding the feature space. In the final data set, each enterprise is described by 972 features for each analyzed year. The list of all features is presented in different tables by data category (see Table 2) and shown in Appendix A: Table A.1, A.2, A.3, A.4, A.5, A.6, A.7, A.8, A.9, A.10, A.11, A.12, A.13, A.14, A.15¹. However, not all features are included in the analysis due to not having the variance condition. These features in Appendix A are crossed out and eliminated from the study.

The data contains 972 features, of which 951 are continuous and 21 binary. Categorical features describe the enterprise type and define belonging to one of the 17 possible sectors (based on the Nace classifier). Most enterprises (∼87% ) are private limited liability enterprises, ∼12% of companies are small communities, and only ∼1% of data consists of public limited liability and individual enterprises. The distribution of Nace values shows that most of the enterprises belong to the Nace C–manufacturing (10.2%), Nace F–construction (11.35%), Nace G– wholesale and retail trade; repair of motor vehicles and motorcycles (28.61%), Nace H–transportation and storage (10.54%), Nace M–professional, scientific and technical activities (12.1%).

Figure 1 represents a distribution of feature values according to classes, i.e., financially distressed and healthy enterprises. Each figure shows one of the most significant features in each data category (see Table 2). In this case, significance is determined using a vote-based feature selection covering all feature selection methods (see ‘Dimensionality reduction methods’). Moreover, feature names contain information regarding the time perspective, i.e., Y indicates the number of years before the outcome; Q represents a quarter of the analyzed year; M shows a month of the analyzed year. If there is no time indication, the feature is affected by the occurrence of certain events. For example, a director change would affect variables representing the frequency of this event and the time since the last change. The distribution of the most significant features justifies the presence of an invisible “optimal cut-off point” among financial distress and healthy enterprises, discussed in Beaver (1966). Hence, more advanced artificial intelligence (AI) techniques should be used to identify the signs of financial distress.

Density plots are rescaled histograms, which use kernel density estimation to reveal the probability density function (PDF) (Y-axis) of the specific feature (Wilke, 2023). A better understanding of the data comes from comparing density plots for different data sets. For example, density plots of features like the history of bad events (9_Num_RgFD_stat_{_}docs) for training and testing data sets are almost the same. This tendency is advantageous due to the ML patterns learning process in the training set, which is also reflected in the test sample. The same tendency could be seen in features 8_Tax_payment_Y1, 7_SD15_Delay_Y1_M12, 6_SEC_R_Retention_ratio_Y1 and 2_ROA_Y1. Nevertheless, upon examination of the 1_Director_change_Num and 4_LawS_def_Num variables, it is observed that the training component is divided into two distinct segments of the trend, of which only one remains in the testing sample. Furthermore, the values of the trend in the training and testing samples for the 5_GDP_perc_Y1_Q3 variable differ substantially.

Dimensionality reduction methods

A large number of features can cause problems related to data sparsity, multiple testing, multicollinearity and overfitting (Kuiziniene et al., 2022). Different dimensionality reduction techniques have been developed to overcome these problems (see Fig. 3). In this study, several dimensionality reduction techniques have been selected based on their relevance concerning different data types (continuous and categorical features), running time, overfitting, and suitability to unseen data.

These dimension reduction techniques can be classified into several categories:

1. Feature selection approach (FS) is used to determine a small subset of informative features from the original large range of data (Al-Tashi et al., 2020) by removing irrelevant, redundant or noisy features. Let Y = {y₁, y₂, …, y_n} denote the class labels, where n is the number of enterprises and y_i ∈ {0, 1}, i =1 , …, n belongs to one of the two classes c₁ – non-financial distress or c₂ - financial distress. Each enterprise is defined by a number of features x_ij, j = 1, …, d. We assume that x_ij are realizations of features X_j that can be ordered according to their importance (2)${\displaystyle \mathit{FI}=\left(X_{\pi \left(1\right)},X_{\pi \left(2\right)},\dots ,X_{\pi \left(d\right)}\right)}$ where π(j), j = 1, 2, …, d is the new descending index of feature X_j.

   1. Group—Filter methods assume complete independence between the model and the data. In this case, features are selected based on their scores in various statistical tests, determining their correlation with the outcome variable.

      1. Cohen’s D (ChD) coefficient is a standardized effect size for comparing the means of Almasri (2022); Imik Tanyildizi & Tanyildizi (2022) two groups (financial distress enterprises and non-financial distressed). Cohen’s D is defined as (3)${\displaystyle \mathit{ChD}=\left|\frac{{\overline{X}}_{x\in c_{1}j}-{\overline{X}}_{x\in c_{2}j}}{\sqrt{\frac{{\sigma }_{x\in c_{1}j}^2+{\sigma }_{x\in c_{2}j}^2}{2}}}\right|}$ where ${\overline{X}}_{\ast j}$ and σ_∗j are means and standard deviations of corresponding groups. The order of the feature’s importance is determined by the effect size, i.e., the greater the ChD, the larger the effect size.

      2. Correlation (Corr) coefficient analyzes possible linear association between variables (Cheng et al., 2020; Jiarpakdee, Tantithamthavorn & Treude, 2020; Séverin & Veganzones, 2021). Based on the type of feature, we use Pearson and point-biserial correlation: (4)${\displaystyle \mathit{Pearson correlation}=r_{X_j{X}_k}=\left|\frac{\sum _{i=1}^n\left(x_{ij}-{\overline{X}}_j\right)\left(x_{ik}-{\overline{X}}_k\right)}{{\sigma }_{X_j}{\sigma }_{X_k}}\right|}$ (5)${\displaystyle \mathit{Point\; -\; biserial correlation}=r_{pb}=\left|\frac{{\overline{X}}_{x\in c_{1}j}-{\overline{X}}_{x\in c_{2}j}}{{\sigma }_{X_j}}\sqrt{pb}\right|}$ where p and b denotes the proportion of cases in the classes of â financial distressâù and “non-financial distress”. σ_Xj is a standard deviation of feature j for all enterprises. Correlation with the outcome variable determines the order of the feature’s relevance (higher is preferable). A high correlation (> 0.75) between features indicates the problem of multicollinearity. To deal with this problem, from the group of highly correlated features, only one is left, having the highest correlation with the outcome, while others are removed from the sample. Finally, all the remaining features are ranked by the strongest relationship with the outcome.

      3. The Kruskal–Wallis (K–W) test is used to compare medians among two or more groups. The Kruskal–Wallis test is chosen because it does not require any prior assumptions about data distribution and is good for all types of features (Devi Priya et al., 2022). The selected significance level is 0.05, i.e., smaller p values suggest that there are substantial differences between the two groups (Linja et al., 2023). The test statistic for the Kruskal–Wallis test is defined as follows: (6)${\displaystyle \mathit{Kruskal\; -\; Wallis}{\chi }^2=H=\left(n-1\right)\frac{\sum _{l=1}^c{n}_l{\left({\overline{R}}_l-\overline{R}\right)}^2}{\sum _{l=1}^c\sum _{j=1}^{n_l}{\left(R_{lj}-\overline{R}\right)}^2}}$ where n – the sample size, n_l is the number of observations in class l, R_lj is the rank of feature j from class l, ${\overline{R}}_l$ the mean rank of all features in class l, ${\overline{R}}_l$ is the mean of all the R_lj, $\overline{R}$ is the mean of all the R_lj. The order of feature importance can be defined by the lowest p − value or the highest |H| value.

      4. Group—Embedded methods include the feature selection in the model fitting process. Classification methods that perform embedded feature selection include:

      1. The least absolute shrinkage and selection operator (LASSO) is a method that combines feature selection and regularization (Altman et al., 2022; Zizi et al., 2021; Huang, Wang & Kochenberger, 2017). The LASSO has the advantage of maintaining the stability of ridge regression and outperforming stepwise regression models (Li et al., 2021). The LASSO function aims to solve (Tibshirani, 1996): (7)${\displaystyle \underset{{\beta }_0,{\beta }_i}{min}=\left\{\sum _{i=1}^n{\left(y_i-{\beta }_0-x_{ij}^T\beta \right)}^2\right\}\text{s.t.}\sum _{j=1}^d\left|\right.{\beta }_j\left|\right.\le \lambda .}$ where β: = (β₁, β₂, …, β_d) is a coefficient vector and λ > 0 adjusts the sparsity of the estimator by eliminating features from the model forcing β_j → 0. The importance of a feature is the absolute value of the LASSO regression coefficient.

      2. Random Forest (RF) uses a permutation importance measure, which is shown to be an effective tool for feature selection (Gregorutti, Michel & Saint-Pierre, 2015). RF offers two techniques for calculating the relevance of all features, resulting in a feature importance rank. These techniques create a rank by considering both the feature cost and differentiating ability (Zhou, Zhou & Li, 2016):

         1. Mean Decreasing accuracy (MDA) (Han, Guo & Yu, 2016) is a feature importance measure based on OOB (out of bag) error. Assume h_t(X_i) and h_t(x_ij) refer to the predicted label for OOB instance X_i before and after feature permutation (Wang, Yang & Luo, 2016). The relevance of a feature X_j is assessed by MDA, which estimates the mean drop in OOB accuracy before and after feature X_j permutation: (8)${\displaystyle \mathit{MDA}\left(X_j\right)=\frac{1}{n_{tree}}\sum _{t=1}^{n_{tree}}\frac{\sum _{i\in \mathit{OOB}}I\left(y_i=h_t\left(X_i\right)\right)-\sum _{i\in \mathit{OOB}}I\left(y_i=h_t\left(x_{ij}\right)\right)}{\left|\right.\mathit{OOB}\left|\right.}}$

         2. Mean Decreasing Gini (MDG) calculates the overall reduction in node impurity (e.g., Gini index) as a result of splitting on the feature and averaging it across all trees (Wang, Yang & Luo, 2016). The MDG is defined as follows: (9)${\displaystyle \mathit{MDG}\left(X_j\right)=\frac{1}{n_{dot}}\left[\right.1-\sum _{j=1}^{n_{dot}}\mathit{Gini}{\left(l\right)}_j\left]\right.}$ where Gini(l)_j is the jth Gini index of feature X_j among the n_dot tree nodes. The Gini index is defined as follows: (10)${\displaystyle Gini=\sum _{l=1}^c{l}_j\left(1-l_j\right)=1-\sum _{jl=1}^c{l}_k^2}$ Where l_j is the probability of the feature X_j classification in the distinct class (Zhang et al., 2019).

      3. Extreme Gradient Boosting Machine (XGBoost or XGB) method assigns a relevance value to each feature based on its involvement in the outcome decision-making using boosted decision trees by Gain metrics (Zheng, Yuan & Chen, 2017) ((11)). (11)${\displaystyle Gain=\frac{1}{2}\left[\right.\frac{{\left(\sum _{i\in I_L}{g}_i\right)}^2}{\sum _{i\in I_L}{h}_i+{\gamma }_1}+\frac{{\left(\sum _{i\in I_R}{g}_i\right)}^2}{\sum _{i\in I_R}{h}_i+{\gamma }_1}+\frac{{\left(\sum _{i\in I}{g}_i\right)}^2}{\sum _{i\in I}{h}_i+{\gamma }_1}\left]\right.-{\gamma }_2}$ where g_i, h_i indicate the first and second-order gradients, γ₁ and γ₂ are regularization parameters, I = I_L∪I_R, while I_L and I_R represent all the left and all the right nodes after each split, respectively (Jiang et al., 2023).The greater the XGBoost gain, the more useful and significant the feature (Chen et al., 2021).

      4. Group—Hybrid combines different methods together:

      1. Overlapping Features (Over_feat)–features having a rating in all feature selection methods, e.g., if LASSO excludes a feature it will not be further analyzed.

      2. Voted importance (V oted_imp)–a joint rank combining results from all feature selection methods. Each feature selection method ranks the features from the most important to the least. Based on these ranks, we can calculate the cumulative weight for each feature: (12)${\displaystyle \mathit{V\; oted\_imp}\left(X_j\right)=\sum _{j=1}^m{\mathit{FI}}_j}$ where FI_j is the rank of the feature and m is the number of used feature selection approaches (in I and II group).

2. Feature extraction approach (FE) maps high-dimensional data X into a new, lower-dimensional space (Z = [z₁, z₂, …, z_n] ∈ ℝ^n×k, where k < d) retaining as much information as possible (Ayesha, Hanif & Talib, 2020):

   1. Group—Linear methods a transformation of data, which projects data linearly:

      1. Linear Discriminant Analysis (LDA) projects the data into a new feature space by minimizing the scatter within classes (S_w) and maximizing the scatter between classes(S_b) (Djoufack Nkengfack et al., 2021). This method only creates as many linear discriminants as there are number of classes, minus one (Anowar, Sadaoui & Selim, 2021). (13)${\displaystyle S_w=\frac{1}{n}\sum _{l=1}^c\sum _{i=1}^{n_l}{\left(X_{il}-{\overline{X}}_g\right)}^T\left(X_{il}-{\overline{X}}_l\right)}$ (14)${\displaystyle S_b=\frac{1}{n}\sum _{l=1}^c{n}_l{\left({\overline{X}}_l-\overline{X}\right)}^T\left({\overline{X}}_l-\overline{X}\right)}$ where c is the number of classes, $\overline{X}$ – the global mean of all classes (Anowar, Sadaoui & Selim, 2021; Djoufack Nkengfack et al., 2021), and eigenvectors are defined as: (15)${\displaystyle W=eig\left(S_w^{-1}{S}_b\right)}$ where the largest eigenvalues correspond to k-eigenvectors of $\left(S_w^{-1}{S}_b\right)$, and of linear discriminants forms W = (w₁, w₂, …, w_k) (Djoufack Nkengfack et al., 2021; Ayesha, Hanif & Talib, 2020) (16)${\displaystyle Z^{\ast }=X{W}^T}$

      2. Principal component analysis (PCA) is the most extensively used dimensionality reduction approach summarizing a large set of features in low dimension with minimum loss of information (Ma & Park, 2022). The PCA method transforms a set of features into uncorrelated Principal Components (PCs), which are ranked and ordered by their importance, e.g., the first PC explains the most variance, the second PC explains the most variance in what is left once the effect of the first PC is removed, etc. (Ayesha, Hanif & Talib, 2020). The general form of PCA transformation can be estimated as a linear weighted combination of features: (17)${\displaystyle Z^{\ast }=XW}$ where Z = [z₁, z₂, …, z_n] ∈ ℝ^n×d are principal components and W = [w₁, w₂, …, w_d] ∈ ℝ^d×d are eigenvectors by Ayesha, Hanif & Talib (2020). The first eigenvector is determined as follows: (18)${\displaystyle w_1=\underset{\left|\right|w\left|\right|=1}{argmax}\left\{\left|\right|Xw|{|}^2\right\}}$

      3. The factorial analysis of mixed data (FAMD) was proposed for handling data containing continuous and categorical features (Ran, 2019). FAMD is a combination of PCA and MCA (multiple correspondence analysis) methods (Ma & Park, 2022). Based on the types of features, relationships between them are measured using Pearson’s correlation coefficient, the squared correlation ratio or the chi-squared test (Momenzadeh et al., 2021; Josse, 2016). (19)${\displaystyle F_s=\underset{F_s\in {\mathbb{R}}^n}{argmax}\sum _{j_{cont}=1}^{p_{cont}}{r}^2\left(F_1,X_{j_{cont}}\right)+\sum _{j_{cat}=1}^{p_{cat}}{\eta }^2\left(F_1,X_{j_{cat}}\right)}$ where F_s is the orthogonal condition to $F_s^{\prime }$ for all s′ < s, p_cont and p_cat are the number of continues and categorical features, respectively (Josse, 2016). The general transformation of FAMD can be expressed as: (20)${\displaystyle Z^{\prime }=X{F}_s}$

      4. Union of separate PCA & FAMD models (Union_PCA_FAMD) select the appropriate data transformation depending on a data type in the analyzed category (see Table 2). In total, there were 12 different categories analyzed, and for three of them, subcategories were created: (1) Financial statements: records and ratios; (2) Sector: class and ratios; (3) Social insurance: employees and delay of taxes (including analyzes of time and amount). The PCA or FAMD method is chosen depending on the data type (continuous or categorical). (21)${\displaystyle f\left(x\right)=\left\{\begin{array}{cc}{Z}^{\ast },\hfill & \text{if all features belonging to the category are continuous}\hfill \\ Z^{\prime },\hfill & \text{if category contains mixed type features}\hfill \end{array}\right.}$

      5. Group – Non–linear methods reduce the dimensionality by reflecting intrinsic non–linear relations between variables.

      1. Autoencoder (AE) is a neural network which is very efficient in generating abstract features of high–dimensional data by minimizing reconstruction loss (Meng et al., 2017). This loss function sets the target value equal to the input by applying a backpropagation algorithm. The main idea of the process is the usage of the encoder–decoder. Firstly, the encoder process transforms input data to lower-dimensional space (to a bottleneck or middle layer); secondly, the decoder converts it back into high-dimensional space (Yan & Han, 2018). The minimization of reconstruction loss assures that the middle layer incorporates most of the information from the original input space (Phadikar, Sinha & Ghosh, 2023). Thus, the middle layer holds the reduced representation of the input data (Kunang et al., 2018). (22)${\displaystyle X^{\prime }=g_{\theta }\left(Z\right)=s\left(W^{\prime }Z+b_Z\right)}$ where f(X) is an encoder function, g(Z) is a decoder, s is an activation function, W represents weights and b is the bias for data x. The training process of the autoencoder finds parameters θ = (W, b_X, b_Z) by minimizing the reconstruction loss (Kunang et al., 2018): (23)${\displaystyle \theta =minL\left(X,X^{\prime }\right)=minL\left(X,g\left(\left(f\left(X\right)\right)\right)\right)}$

      2. Union of separate autoencoders (Union_AE) trains several autoencoders depending on the analyzed data category. Analyzed data categories are the same as in Union_PCA_FAMD.

Number of features

The scientific literature lacks information on how much reduction can be performed on the original data without negatively impacting the model accuracy, i.e., how to identify the minimum data set necessary for efficient modelling. In order to fill this gap, this study analyzes different subsets of features and their influence on the model’s performance. For methods based on certain thresholds, a fixed set of effective features is determined, i.e., LASSO algorithm k = 24, Over_feat k = 24, and for ChD k = 14. Nonetheless, there is no established methodology for determining the threshold, thus an established set of effective features has been determined experimentally by utilizing a diverse number of features from a lower-dimensional space, i.e., k ∈ {15, 30, 50, 100}.

In the case of feature extraction, the number of features is predefined by a specific threshold, e.g., for LDA k_LDA = number of classes − 1 (Anowar, Sadaoui & Selim, 2021), for PCA we analyze features whose cumulative variance exceeds a predefined threshold (70–90% of the total variance) (Aggarwal, 2015), for autoencoder the number of features depends on the minimization of reconstruction loss of the training sample (Phadikar, Sinha & Ghosh, 2023). However, for the comparison with feature selection approaches, the selection of k is the same for PCA and FAMD methods, which explains from 89.8% (100 PC’s) to 67.5% (15 PC’s) of the total variance. In the analysis of Union_PCA_FAMD we make two data sets, which explain ≥70% and ≥80% of the total variance in each data category, and the number of features is k_{Union_PCA_FAMD} ∈ {64, 93}, respectively. To analyze the efficiency of non–linear methods, we develop several autoencoder structures. The different number near autoencoder – I, II, III identifies a number of hidden dense layers between the input layer and the middle layer; the number of units in the middle layer is k_AE ∈ {32, 64, 100}. Finally, the Union_AE uses k_{Union_AE} ∈ {42, 74}. It should be noted that in the case of k_{Union_AE} = 74, the obtained number of features is high. Therefore, the quantity of features is reduced until all of them satisfy th e < 0.7 correlation requirement, resulting in a feature set comprising of 42 features.

Machine learning methods

Machine learning (ML) techniques cover a wide range of algorithms addressing various types of problems. There are two main groups of ML algorithms: supervised and unsupervised learners (Algren, Fisher & Landis, 2021).

1. Supervised machine learning techniques are used to construct predictive models. These algorithms try to model relationships and dependencies between the dependent variable and independent variables. In machine learning, dependent variables are frequently referred to as labels, targets or classes and independent variables are known as features (Algren, Fisher & Landis, 2021). This study employs boosting and neural network techniques:

   1. The Boosting technique seeks to increase the prediction accuracy by combining results from multiple classifiers (Nettleton, 2014). In this study, we focus on AdaBoost, Categorical Boosting (CatBoost) and Extreme Gradient Boosting Machine (XGBoost) techniques.

      1. AdaBoost—an ensemble approach that trains and deploys trees one after the other. It implements the boosting technique by combining multiple weak classifiers to build one strong classifier. Each weak classifier attempts to improve the classification of samples that were misclassified by the previous weak classifier (Misra & Li, 2020).

      2. Categorical Boosting (CatBoost)—a machine learning technique that is based on gradient boosting decision trees. A powerful machine learning technique called gradient boosting can handle challenges with diverse features, noisy data, and complex dependencies. Moreover, this approach can easily handle categorical features by substituting the original categorical variables with one or more numerical values. Ordered boosting is a new technique that replaces the traditional gradient estimation methods (Zhang, Zhao & Zheng, 2020).

      3. Extreme Gradient Boosting Machine (XGBoost)—a tree ensemble model that can be applied for both classification and regression problems. Its main idea is to make the target function as minimal as possible while employing the gradient descent method to produce new trees based on all previous trees. While using XGBoost to solve regression issues, incrementally new regression trees are added, and the residuals of the prior model are subsequently fitted using the newly created CART tree. The final predicted value is the sum of the outcomes of each tree (Pan et al., 2022). Selecting the loss function is a crucial step in setting XGBoost models. Loss functions are used for classification and regression problems; specifically, we use the “binary” loss function for binary classification (Brownlee, 2021).

   2. Discriminant analysis (DA)—a classification and dimension reduction technique described in ‘Dimensionality reduction methods’ (Fisher, 1936). The discriminant analysis creates linear combinations separating objects into classes, assuring that the variance within the class would be minimized and the variance between classes—maximized. These directions are known as discriminant functions, and their number is equal to the number of classes minus one (Canizo et al., 2019). In this study, we use two types of DA: Linear DA and Quadratic DA, which enables the non–linear separation of data (Fisher, 1936).

   3. Decision trees (DT)—a powerful data mining tool that’s often used for feature selection and classification tasks. A decision tree may be used in data mining as a primary model-building method, or for automated feature selection (Bunge & Judson, 2005). Decision trees assist in identifying the most important characteristics that support precise classification by assessing the relevance of each feature at each node split (Sugumaran, Muralidharan & Ramachandran, 2007).

   4. K-nearest neighbor (KNN)—probably the most well-know non-parametric classification technique (Chanal et al., 2022). This method classifies each unlabeled observation by the majority class among its k-nearest neighbors in the training set (Chanal et al., 2022). In this study, we perform experiments with k = 3, 5, 7 nearest neighbors.

   5. Logistic regression (LR)—contrary to its name, LR is a classification model rather than a regression model (Subasi, 2020). For binary and linear classification issues, logistic regression provides a faster and more efficient solution. This classification approach performs very well with linearly separable classes and is relatively simple to implement. For the LR method, the assumption of multicollinearity is fulfilled by removing features, highly correlated with other features. Furthermore, a significance level of p-value is utilized, and features with a p-value greater than 0.05 are also eliminated. For example, using the K_W feature selection approach, we have identified 100 significant features that have been used for LG model development. However, to avoid the problem of multicollinearity, we have performed additional analysis and have removed highly correlated features. Thus, the final LG model has been developed using 20 statistically significant features.

   6. Naive Bayes (NB)—based on the Bayes theorem, the naive Bayes classifier is a member of the family of probabilistic classifiers. It is based on the idea that a feature’s presence in a class is independent of any other features that may also be present in that class. For example, if an enterprise has a debt of 100,000 euros to a social insurance institution, suffers financial losses, and the turnover is consistently decreasing, it can be categorized as bankrupt. All three of these independent contributing factors are considered by a naive Bayes algorithm (Krishnan, 2021).

   7. Neural networks are artificial intelligence models that are designed to mimic the human brain’s functions (Casas, 2020). We use three different neural network methods:

      1. Artificial neural networks (ANN (I-III))—the model that uses calculations and mathematics to imitate human psychology. The unique architecture format used by ANN models mimics a biological nervous system. The ANN models are made up of neurons that interact in a complex, non–linear manner, just like the human brain. Weighted links connect the neurons to one another (Malekian & Chitsaz, 2021). The hidden structure of the neuron network has been marked I-III, which indicates hidden layers between input and dense layers. After each layer (except a dense one) is implemented, a drop-out layer is excluded.

      2. Convolutional neural network (CNN (II, III, V))—a specific kind of feed-forward neural network used in AI. The input data to CNN is displayed as multidimensional arrays. It works well for a large number of labelled data. Based on the important function played by the receptive field, it gives weights to each neuron, such that it can distinguish between the relative importance of different neurons. Three different layer types create CNN’s architecture: convolution, pooling, and fully linked (Shajun Nisha & Nagoor Meeral, 2021). The hidden structure is indicated the same as in ANN. However, this is only an indication for used conv_1d and flatten layers; input, drop out, max-pooling and dense layer are not included in the hidden layers structure calculation.

      3. Extreme learning machine (ELM (100,150,200,300))—is essentially a single feed-forward neural network. The weights between inputs and hidden nodes are distributed at random throughout its single layer of hidden nodes. As a result, the parameters of the model can be calculated without the need for a learning process, and they remain constant throughout the training and prediction stages. On the contrary, the weights that connect hidden nodes to outputs can be learned exceedingly fast (Garza-Ulloa, 2022).

      4. Self-organizing map (SOM)—takes a set of input data and maps it onto neurons of a (usually) two-dimensional grid. Each neuron in the 2D grid is assigned a weight vector with the same dimensionality as that of the input vector. The weights represent the associated memory. Finding a winning neuron and modifying the weights of the winner neuron and its nearby neurons are the two fundamental steps in the SOM computational process (Köküer et al., 2007).

   8. Random Forest (RF)—an algorithm that uses classification and decision regression trees as the foundation for its ensemble machine learning approach (Xia, 2020). Using the values of the predictor variables as inputs, decision tree learning attempts to build a statistical prediction model to predict the values of the response variable(s). The model is created by establishing a value for the response variable within each of the recursively partitioned predictor variable spaces, which generates a decision tree.

   9. Support vector machine (SVM)—an approach for supervised machine learning that can be applied to both classification and regression problems. The SVM algorithm seeks a hyperplane in the data space that creates the largest minimum distance (referred to as a “margin”) between the samples belonging to different classes (Satapathy et al., 2019). In this research, SVM kernel polynomial and radial basis functions are used.

2. Data labels are not used in unsupervised learning methods. Instead, discovering a structure within a dataset, “cleaning up” the data, or grouping data points into groups are the tasks that the algorithm is expected to complete (Algren, Fisher & Landis, 2021). In this study, three different unsupervised methods are employed:

   1. Isolation Forest (IF)—an unsupervised, non-parametric method for anomaly detection in multivariate data. This outlier detection method seeks to isolate anomalies from the rest of the data using an ensemble of decision trees. It constructs a forest of isolation trees, or iTrees. Partitioning data is required for creating iTrees. Furthermore, iTrees are binary trees that use recursive partitioning to isolate abnormalities. These anomalies, which are “rare and diverse”, are therefore “more sensitive to isolation and hence have short path lengths” (Chater et al., 2022).

   2. K-means clustering (Kmeans) divides the set of samples into K clusters, each described by the centroid, i.e., the mean of the samples in a cluster. The k-means algorithm uses an iterative approach to find the optimal cluster assignments by minimizing the sum of squared distances between observations and assigned cluster centroid (Ashour, Guo & Hawas, 2019).

   3. One-class support vector machine (One class SVM)—an approach that identifies unusual observations compared to the instances of the known class. The one-class SVM attempts to classify one class of objects and separate it from all other possible objects. The classifier can correctly classify certain objects, but the others will be classified as outliers. Thus, a one-class classifier is trained to label these observations as outliers (Seo, 2007).

Evaluation metrics

Different performance metrics can be obtained based on the confusion matrix (see Table 3). In Table 3, TP denotes the number of true positives, TN is the number of true negatives, FP is the number of false positives, FN is the number of false negatives, POS is the number of actual positives, PPOS is the number of predicted positives, NEG is the number of actual negatives, PNEG is the number of predicted negatives, N is the number of all instances (Fernández et al., 2018). In this study, classes are defined as follows: positive class—non-financial distress enterprises and negative class– financial distress enterprises.

The most commonly used evaluation metrics are provided in Eqs. (24)–(30) (Kuiziniene et al., 2022). Different accuracy measures were used to get insights into different aspects of model performance. The interpretation of the results is common for all metrics – higher values indicate better performance.

Precision—the ratio of true positives (TP) to predicted positives (PPOS). (24)${\displaystyle \mathit{Precision}=\frac{TP}{TP+FP}.}$

Recall—the ratio of the true positives (TP) to actual positives (POS), also known as sensitivity or TPR (true positive ratio). (25)${\displaystyle \mathit{Recall}=\frac{TP}{TP+FN}.}$

Specificity—the ratio of the true negatives (TN) to the actual negatives (NEG). Also known as TNR (true negative ratio). (26)${\displaystyle \mathit{Specificity}=\frac{TN}{TN+FP}.}$

The area under the ROC curve (AUC)—a measure of how well a model can distinguish between two classes and is expressed as follows: (27)${\displaystyle \mathit{AUC}={\int }_0^1\left(\mathit{TPR}\right)d\left(\mathit{FPR}\right)}$

where false positive ratio FPR = 1–specificity.

The Gini index is an alternative to AUC, which is more popular in the context of bankruptcy prediction. Gini is a metric that indicates the model’s discriminatory power. Moreover, the simple expression of Gini is: (28)${\displaystyle \mathit{Gini}=2\mathit{AUC}-1.}$

The proof of formula transformation from Eqs. (10) to (28) can be found in the articles (Kaymak, Ben-David & Potharst, 2012; Liu et al., 2018).

Accuracy (ACC)—the proportion of correctly classified instances. However, under highly imbalanced problems, it is simple to achieve high accuracy (Fernández et al., 2018). (29)${\displaystyle \mathit{ACC}=\frac{TP+TN}{TP+TN+FP+FN}}$

F − score is the harmonic mean of precision and recall, where the most common value of β is 1. Therefore, the estimate is often called F1 or F1 score. However, this measurement focuses on the analysis of the positive class (Fernández et al., 2018). (30)${\displaystyle \mathit{F\; -\; score}=\frac{\left(1+{\beta }^2\right)Precision\cdot Recall}{{\beta }^2\cdot Precision+Recall}.}$

Different from F1, the G–mean focuses on the balance between instances of majority and minority classes (Fernández et al., 2018); it is a geometric mean of the true positive rate and the true negative rate. (31)${\displaystyle \mathit{G\; -\; mean}=\sqrt{TPR\cdot TNR}.}$

A benchmark model

A benchmark financial distress model has been added to the research to improve comparability with the proposed methodology. The logistic regression approach employing financial ratios data serves as a benchmark model, since it is the most commonly employed approach in research on financial distress or bankruptcies (Kuiziniene et al., 2022).

Foremost, the benchmark model has been implemented in correlation analysis. If two features are highly correlated (>0.7), then the feature with a lower correlation to the target is removed. In the logistic regression, the Wald test is commonly used to identify statistically significant variables. The selected significance level for this research is 0.05. Features used in the benchmark model are in bold and underlined (see Appendix A, Table A.3).

The benchmark model’s results are presented in Tables 5 and 4. The confusion matrix of the benchmark model shows that the effectiveness criteria are fulfilled – the model correctly classified 1,180 Financial distress cases and 56,994 Non-financial distress cases (see Table 5). However, the AUC score was found to be only 0.6991. The benchmark model is outperformed by the combinations of methods and application techniques proposed in the methodology. Moreover, when ranked by the AUC metric, the benchmark model is ranked 984 out of 1071 experiments (see Table 4).

Feature selection

In this section, only the results of the feature selection approach and the analysis of two research questions are provided, i.e., RQ2–which data categories are the most important in financial distress prediction? and RQ3–does the composition of the final feature set depend on the feature selection technique?

Figure 6A shows the AUC results for the best machine learning models for each feature selection technique with a certain number of features, i.e., constant parameters are: feature selection technique and the number of features. To give an example, 26 experiments were done using the RF_MDG feature selection technique with 30 significant features. The one with the highest AUC score was found to be 0.8687. The best AUC results for ML models with different numbers of features are distributed between 0.8012 and 0.8687 for all used FS methods with different features. Embedded and hybrid feature selection methods achieve a relatively high accuracy, i.e., AUC score ≥ 0.85. Therefore, future experiments should be conducted with these groups on the topic of financial distress. Wrapper methods in this research have not been analyzed due to the long estimation time and limited computation resources. In further research, this method could be implemented after the embedded method feature ranking and the selection of an initial number of features. This study indicates that the best results are achieved when using sets of up to 50 features or fewer. Also, using simple (not deep) neural network methods to recognize financial distress in enterprises is found to give the best results.

Figure 6B shows how much each data category (see Table 2) is used in respect to the other categories. It can be observed that all FS methods use features from Financial statements, Social Insurance and Other category. Contrary to this, Macro, Enterprise type and Sector’s data categories are rarely selected and only by a few FS methods. A comparison of filter and embedded FS techniques shows that the filter methods avoid using features from the Taxes data category.

Different FS methods tend to select features from the same data category as the number of features increases, e.g., Social insurance could be split into employees and delay of taxes (including delay in days and amount information). In that case, it can be seen that if FS methods choose Social insurance → employees subcategory, increasing the number of the features would include more features from this particular subcategory.

Feature extraction

This section presents the results from feature extraction experiments as well as the analysis of the research question—RQ4: which feature extraction strategy is more efficient: analyzing all features at once vs analyzing features by categories?

The idea represented in Fig. 7A is the same as for Fig. 6A, but it is intended for the analysis of feature extraction experiments. Upon comparing the results obtained with FE and FS (see Figs. 6A and 7A), it is observed that FE methods’ superior AUC results are distributed lower than those obtained from FS, ranging from 0.6513 to 0.8381 (see Fig. 7A). In addition, the best results for the AE feature extraction technique, with 32 and 100 relevant features, were found to be lower than the benchmark model. The highest AUC ≥ 0.83 from the FE methods category has been achieved with the implementation of the union strategy for PCA_FAMD and AE methods. This also resulted in the highest number of used features in PCA and FAMD methods.

In this study, the union of separate autoencoders (see the ‘Dimensionality reduction methods’ section; Union_PCA_FAMD and Union_AE) has been implemented for both linear and non–linear data transformations separately, with respect to the analyzed data category. This approach involves the creation of data categories with reduced dimensionality and their connection to a single feature space. The usage of the Union strategy for linear and non–linear data transformation has given the best AUC score when compared to all ML methods separately. According to the distributions presented in Fig. 4A, only 53.11% of the experiments were completed using FE methods because the minimum classification quotas were not reached. Figure 7B shows the AUC score distribution of ML models using FE methods. Moreover, Union_PCA_FAMD and Union_AE achieved efficient results—73.21% and 77.19%, respectively (see Fig. 4C). This confirms the union of separate autoencoders strategy’s viability and proves that its development in further research could lead to better results. In addition, better Union results are achieved when a smaller number of features is present, rather than a larger one. Moreover, a linear transformation provides better results than a non–linear one.

Machine learning models

This part presents the analysis of the last research question—RQ5: which machine learning model performs better in financial distress identification?

In this study, experiments with 17 different machine learning models and their 28 modifications have been provided. However, not all experiments have succeeded due to the time limit for a single method (24 h). Usually, the cycle of all ML models, for the one-dimensionality reduction method, takes from 2 to 4 days. A comparison of models and their modifications is presented in Fig. 8. Blue color indicates the better models. Comparing dimensionality reduction categories, i.e., FE and FS, the efficiency trends remain the same between models. However, higher AUC score is achieved in the FS approach. What concerns the comparison of tree-based algorithms (Fig. 8C), it is difficult to select the most effective methods as their AUC scores are similar. Moreover, in the FS approach the achieved AUC score variance is small, despite the used FS method, which shows the stability of the ML method and equally important feature selection.

As the FS approach has a similar tendency, a larger sample size with more significant variables, and their number of included features in the analysis, the comparison of neural network methods is given. The main conclusion remains the same – simpler (not deep) neural network methods give better results in recognition of enterprises that have financial distress. However, the more features are included in the analysis, the less stable the neural network results become if all conditions remain the same, except when the FS method is used. This implies that the FS method is important. In further research, it is necessary to expand the neural network structure searches for each dimensionality reduction method separately in order to achieve a deeper comparison of the results.