Fuzzy evaluation and explainable machine learning for diagnosis of rheumatic and autoimmune diseases

Dataset collection

The problem of diagnosing rheumatic and autoimmune diseases to avoid sample bias is essential for obtaining an extensive and well-organized dataset to ensure that the model is optimal for all types of patients and diseases. Our work data were collected between 2019 and 2024 from various medical factors in Iraq, including one hospital and three laboratories of a wide range of autoimmune conditions, confirming variety in the subjects of the work in terms of familial and disease presentation (Mahdi, Jahani & Abd, 2025; Abd, Mahdi & Jahani, 2025). Special care was taken to ensure that the dataset was diverse and anonymous to facilitate the privacy of the patients while allowing for successful training and model evaluation. This study was approved by the Saint Raphael Hospital and Al-Hawraa, Al-Sibtain, and Taqadum health centers under approval numbers 25, 66, 98, and 122, respectively. Table 2 lists 14 features corresponding to the rheumatic and autoimmune diseases.

Table 3 shows the patient distributions of certain rheumatic and autoimmune diseases and the normal healthy population. Among these, RA is the most abundant with 2,848 cases, making this work even better with minimum cases, as it is proven to be easier to obtain in an information system. AS and SS also have a reasonable number of cases (2,127 and 1,852, respectively), which helps to model the relevant conditions. In total, 1,783 and 1,604 cases of PA and N, respectively, have been reported, which still provides sufficient cases for classification even with the intricacy of the disease. The SLE class (1,355 cases) is necessary for helper states since it helps to distinguish between healthy and diseased states, reducing the chances of false positives. Although the smallest ReA class with only 516 cases is faced with many challenges regarding the model’s ability to identify its characteristics. This class is imbalance, particularly the small size of the ReA class compared with the RA class, might have a negative effect on the performance of the model in terms of discrimination, especially where these diseases are more common and thus expected. To overcome this problem, methods such as resampling, class sparsity, and data enhancement for class imbalance is applied so that the model performs equally well in all diseases, even those that are rare, such as ReA.

Preprocessing

One of the most important stages in the process of analyzing any data is the preprocessing, as it enables the model to perform well (Zhu, Jiang & Alonso, 2024). It encompasses several steps that aim to cleansing, converting, and structuring the data in a manner that increases the efficiency of the ML algorithms. Three preprocessing steps were used in this work:

• –

  Missing values are handled via the mean and mode (Palanivinayagam & Damaševičius, 2023). The mean imputation for numerical variables is shown in Eq. (1). This method is useful only if the data are missing in some cases owing to randomness and if the population is approximately normal. Additionally, for categorical variables, the missing values are imputed with the mode, that is, the maximum occurring category, as shown in Eq. (2) (1) $$m_i=\frac{{\sum }_{j=1}^n{m}_j}{n}\mathrm{f}\mathrm{o}\mathrm{r}{m}_j\in \{m_1,m_2,m_3,\dots ,m_j\}$$ (2) $${\beta }_{mode}=\underset{{\beta }_k}{\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{a}\mathrm{x}\mathrm{f}\mathrm{r}\mathrm{e}\mathrm{q}\mathrm{u}\mathrm{e}\mathrm{n}\mathrm{c}\mathrm{y}({\beta }_{\mathrm{k}})}$$

• where $m_i$ is the missing value, ${\sum }_{j=1}^n{m}_j$ is the sum of all nonmissing values, $n$ is the number of nonmissing values, ${\beta }_{mode}$ is the most frequent category, and ${\beta }_k$ represents each distance category in the category feature.

• –

  ML approaches tend to be more efficient when applied to numeric inputs; hence, the need to encode categorical variables into numbers, for example, is recommended. We use one-hot encoding. In which each category value is represented by a binary vector, where only one element is 1 and the other is 0, as shown in Eq. (3). (3) $$One\text{-}Hot(C_i)=[b_1,b_2,b_3,\dots ,b_n]where{b}_j=\{\begin{array}{c}1ifj=i\\ 0otherwise\end{array}$$

• where C is the categorical variable, $n$ is the total number of unique C, and $b$ is the binary vector.

• –

  In medical datasets for definition and classification, targets with imbalanced classes (e.g., rare diseases) are usually available, and for such cases, the adaptive synthetic sampling (ADASYN) is said to be an effective method for creating synthetic samples of the minority class (Khan et al., 2021; Davagdorj et al., 2020). ADASYN focuses on instances that are difficult to learn and modifies itself according to the distribution of the data. This method address the class imbalance problem without generating uniform oversampling, resulting in optimal representations of the minor classes. This method prevents bias toward common classes and improves rare disease analysis. The following steps for ADASYN algorithm:

• Let the training dataset be $D=\{(x_i,y_i){\}}_{i=1}^n$, where $x_i\in {\mathbb{R}}^d$ and $y_i\in \{0,1\}$, with class one as the minority class.

  • 1.

    Identify minority class:

  • Let $m$ be the number of minority class examples. Define $\beta$ as the desired balance level, and compute the total number of synthetic samples to be generated, as shown in Eq. (4). (4) $$G=(n_{\mathrm{m}\mathrm{a}\mathrm{j}}-n_{\mathrm{m}\mathrm{i}\mathrm{n}})\times \beta$$

  • where $n_{\mathrm{m}\mathrm{a}\mathrm{j}}$ and $n_{\mathrm{m}\mathrm{i}\mathrm{n}}$ are the number of majority and minority samples, respectively.

  • 2.

    Compute local imbalance per minority sample:

  • For each minority sample $x_i$, find its $k$-nearest neighbors ( $k=5$) from the entire dataset. Let ${\eta }_i$ be the ratio of majority samples among the $k$ neighbors of $x_i$, as shown in Eq. (5). (5) $${\eta }_i=\frac{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}\mathrm{o}\mathrm{f}\mathrm{m}\mathrm{a}\mathrm{j}\mathrm{o}\mathrm{r}\mathrm{i}\mathrm{t}\mathrm{y}\mathrm{s}\mathrm{a}\mathrm{m}\mathrm{p}\mathrm{l}\mathrm{e}\mathrm{s}\mathrm{a}\mathrm{m}\mathrm{o}\mathrm{n}\mathrm{g}k\mathrm{o}\mathrm{f}{x}_i}{k}.$$

  • Normalize the ratios to obtain $r_i$, by using Eq. (6). (6) $$r_i=\frac{{\eta }_i}{{\sum }_{j=1}^m{\eta }_j}.$$

  • 3.

    Determine number of synthetic samples per instance:

  • Calculate the number of synthetic samples to be generated for each $x_i$, using Eq. (7). (7) $$g_i=r_i\cdot G.$$

  • This ensures that more synthetic data is generated for those minority samples surrounded by more majority samples.

  • 4.

    Generate synthetic samples:

  • For each $x_i$, randomly select one of its $k$-nearest minority neighbors $x_{zi}$, and generate synthetic samples using the following interpolation, as shown in Eq. (8). (8) $$x_{\mathrm{n}\mathrm{e}\mathrm{w}}=x_i+\delta \cdot (x_{zi}-x_i),\delta \sim \mathcal{U}(0,1)$$

  • Repeat this for $g_i$ times to create the required number of synthetic samples for $x_i$. In this work ReA class increased from 516 to 2,946 samples.

Model training

In this section, after preprocessing, the data is moved to the data training, and different algorithms are practiced. The selection of the 12 ML models was based on their established performance in medical classification tasks that involved tabular data. These models cover a broad spectrum of algorithmic types, including ensemble learning (RF, XGBoost, LGBM, GBoost, AdaBoost, extra trees (ET)), probabilistic learning naive Bayes (NB), distance-based methods k-nearest neighbors (KNN), linear classifiers (SVM, SGD), decision tree (DT), and artificial neural networks (ANN). This diversity allows for comprehensive benchmarking across different algorithmic families. Each of these algorithms is trained on the dataset, and then the different performances are evaluated.

Let us assume that we have dataset D and this dataset will split into training and testing sets. Both sets have the characteristics $X=x1,x2,x3,\dots ,xi$, and we have the labels $Y=y1,y2,y3,\dots ,yj$. The formula becomes as follows:

(9) $$D=x_i,y_j|i=1,2,3,\dots ,n,j=1,2,3,\dots ,m.$$

In this work, D is divided into two sets, $D_{train}$ and $D_{test}$, using ratios of 70% $D_{train}$ and 30% $D_{test}$ for each class, as shown in Table 4.

For training set $D_{train}$, the aim is to optimize the model parameters to minimize the error between the prediction $\hat{y}$ and the actual value $y$. Equation (10) shows how our models are trained.

(10) $${\theta }^{\ast }=\underset{\theta }{\arg min}\psi (g_{\theta }(x),y)$$where $\theta$ represents the model parameters, $g_{\theta }$ represents the model prediction with the parameters, $\psi$ is the loss function for the model, $y$ is the actual label, and ${\theta }^{\ast }$ represents the parameter optimization that was obtained during training, as shown in Table 5.

To mitigate overfitting in GBoost and XGBoost we increased regularization $\lambda$ = 0.2. For ANN algorithm Added dropout (rate = 0.2) and implemented early stopping (patience = 100 iter). Finaly, reduced max tree depth from 6 to 4 and feature subsampling (50% per tree).

After parameter setting and training our model $D_{train}$, the model is tested through $D_{test}$. Equation (11) is used for model prediction for the test data.

(11) $${\hat{y}}_l=g_{\theta }(x_i).$$

To evaluate our model performance base on the prediction $\hat{y}$ different metrics were used, such as accuracy, precision, recall, F-score, kappa, MCC, and HL. All these metrics are explained in the next section. Algorithm 1 shows the steps of the learning model.

Evaluation metrics

In ML, it is essential to examine models to evaluate how well they might perform on a task or whether they are fit for purpose. Different metrics are used to evaluate models on the basis of their prediction accuracy, error types, and ability to extend the model beyond the data used for training. These metrics help to report how well a model performs on some prediction tasks. Table 6 shows the most common evaluation metrics as they relate to model performance (Abd et al., 2023).

FDOSM

This section introduces the decision-making method in which fuzzy logic is used, incorporating subjective opinions and benchmarking ML models. This method is called the FDOSM (Jassim, Abd & Omri, 2023). This mathematical method is used to solve multicriteria decision-making problems with a single decision. Figure 2 shows the three main stages of FDOSM model, which are the data input unit (which implements a decision matrix (DM)), the data transformation unit (which converts the DM into an opinion DM), and the data processing unit (which employs a fuzzy opinion matrix to determine the ranking for each alternative (ML algorithms)).

Dataset analysis

Various clinical features revealed several important aspects, including count and missing data, with descriptive statistics presented in Table 10. The age and gender features are complete without any missing values. Continuous characteristics, such as ESR, CRP, RF or anti-CCP, suffer deficits in values of more than 9% to 27%. Binary features such as HLA-B27, ANA, Anti-Ro, and Anti-La are also features with missing value proportions between 16% and 43%. On the other hand, there are clinical features such as Anti-dsDNA and Anti-Sm that have 43% missing values, which could pose challenges in modeling. Finally, C4 and C3, which have moderate missing values between 14% and 17%. Addressing the missing values, appropriately normalizing the continuous variables, and adjusting binary variables will be necessary to obtain proper and exact modeling, especially since many features are missing for some of the variables.

Figure 3 shows the PCA method in 3D, where points that indicate instances of different classes, corresponding to various diseases. The demonstrative graph represents a visual projection of the dataset on the first two (PC1, PC2), where the PA class appears to be isolated, whereas other pathological classes appear to blend in the middle to a large extent. This blend suggests that there may be some overlap in the features associated with these diseases and that additional or more advanced discrimination methods might be necessary to fully separate them.

Machine learning results

This section provides performance metrics in terms of precision, recall, F-score, kappa, Hamming loss, MCC, accuracy, and other metrics for several ML models for unbalanced datasets and balanced datasets. Table 11 shows the unbalanced dataset, indicating that GBoost, RF, XGBoost, ET, and the LGBM perform better. While, GBoost and the RF outperforming many other metrics exhibiting high precision and well-balanced F-scores with low Hamming loss, exceptional MCC, and kappa scores reflecting high prediction reliability. Others, such as DT, ANN, and SGD, are in the middle tier of the models with an average performance level. However, SVM scores better on the reliability scale and poorer on the accuracy scale. In terms of KNN and NB, the two methods performed poorly, especially in terms of Hamming loss and lower accuracy. AdaBoost is the outlier, which performs hopelessly across all the metrics. In general, the ensemble models performed the best concerning this unbalanced dataset, whereas poorer performing models such as KNN, and NB, while AdaBoost did not perform as well.

Figure 4 shows a diagram of the ROC curve that enables us to perform a comparative assessment of the capabilities of different ML models that are based on distinguishing among classes. This evaluation takes the AUC results in the form of calculating the area under the curve. XGBoost, GBoost, AdaBoost, RF, and LGBM exhibit superior performance classifiably, with corresponding AUCs of approximately 0.99, whereas ET, ANN, and SVM also perform well, with AUCs of 98, 0.97, and 0.94, respectively. Conversely, KNN and SGD delivers the AUC of 0.92, this may be due to their imperviousness to the problem of high-dimensional data. Most apparently, NB, despite its simplifying assumption, records an efficient AUC of 0.91.

For the AUC measure, the DT performs lowest, with an AUC of 0.89, which is a likely overfitting metric. Additionally, all boosting algorithms appear to be more efficient than DT because of easier interactions of noise data. Some modifications in terms of feature scaling such as feature selection to reduce overfitting.

Table 12 outlines the evaluation measures of different combinations of ML models and highlights the level of overfitting. Each model in this case, XGBoost, DT, RF, LGBM, ET, and GBoost, showed overfitting tendencies to varying degrees, and this performance gap between training and validation scores was large, with DT recording of 20.24% drop in training performance. Such models cannot cope with overfitting without adopting regularization or pruning. In contrast, SVM, NB, AdaBoost, and SGD have no overfitting, and all training and validation scores are similar for such models, with some, such as AdaBoost and NB, being less effective in terms of performance. The overall performance of the ANN is balanced, whereby both training and validation do not perform outliers centered on the degree of overfitting that exists. Instead, tree-based models usually require tuning to overcome low performance, while tuning is likely to increase generalization. Further, whereby low-performing models may even are optimized or other measures taken without creating the risk of incorporating generalization into the validity. In this test, we can use Eq. (17) to calculate the percentage if overfitting occurs. It uses the difference between training and validation.

(17) $$Overf=\left(\frac{Train\text{-}Validation}{Train}\right)\ast 100.$$

Table 12 shows three overfitting statuses: yes, no, and mild. In the term yes, we have significant overfitting (overfitting % $\ge$ 10%), whereas mild refers to moderate overfitting (5% $\le$ Overfitting % < 10%). Finally, no refers to no overfitting (overfitting percentage < 5%).

Table 13 shows that the analysis of several ML models in a balanced dataset confirms that in models such as XGBoost, RF, GBoost, LGBM, and ET, the other models outperform in terms of all the metrics. XGBoost yields outstanding results, with the accuracy of 86.131%. Among the methods evaluated, GBoost was ranked as the best performing model, with 86.889% accuracy, indicating its strength. Moreover, linear models such as SGD, NB, and AdaBoost do not perform adequately, and with AdaBoost, the worst performance reflects accuracy of 34.397%. Overall, the performance of the ensemble models, especially GBoost, is superior and reliable; however, all the other models have unpredictable outcomes, and the performance is very poor because the balanced dataset is not favorable, particularly for the linear model AdaBoost. Note that after our dataset is balanced, the size of each class becomes ReA (2,946), RA (2,848), AS (2,127), SS (1,852), PA (1,783), N (1,604), and SLE (1,355).

The ROC curves in Fig. 5 clearly show that other ensemble models, such as XGBoost, RF, Gboost, AdaBoost, and LGBM, are the clear winners, as shown by their almost perfect AUC values, which indicates that such models can learn complex patterns and generalize well on various datasets. ANN performance is classifiable, with a corresponding AUC of 97, whereas SVM and KNN also perform well, with AUCs of 0.94 and 0.93, respectively. Most apparently, NB, despite its simplifying assumption, records an efficient AUC of 0.91. Conversely, KNN and SGD deliver an AUC of 0.90, this may be due to their imperviousness to the problem of high-dimensional data.

From the information in Table 14, the overfitting in various models is compared on the basis of their training and validation performance. XGBoost, DT, RF, LGBM, and ET exhibit heavy overfitting since they attain pleasant training scores, while their validation scores roll down and roll down a high margin, indicating that training is aimed at memorizing certain patterns and there are few opportunities to expose the models to new data. The KNN and gradient boosting models show moderate overfitting, meaning that for these methods, the overall performance is reasonable, although it could be increased with tuning or additional datasets to avoid overfitting. The SVM, ANN, SGD, AdaBoost, and NB methods were delayed or had no overfitting since the training and validation scores were reasonable.

FDOSM results

In this stage, three experts working with ML and evaluation were used. The experts have more than 6 years in this field. Tables 15, 16, and 17 present evaluations from three different experts on various machine learning models on the basis of multiple performance metrics.

Table 18 shows the benchmarking performance scores and the classifications by different ML models by three experts and a final aggregated consensus decision. The final decision column incorporates all the factors and includes the corresponding aggregate score and rank for each model. In the case of GBoost, the best among all the models with all three experts ranked 1st, with a score of 0.1333, which demonstrates perfect compliance with the evaluations. Additionally, the RF performs well and always holds second place, with a score of 0.1571 in each expert’s evaluation, indicating that there is no disagreement about how well it performs. LGBM has been placed in 3rd place in the aggregated score. XGBoost occupies with the final score of 0.2048 and ranks 4th, implying an average performance. In contrast, AdaBoost emerged as the least rated model 12th by all the experts, with a lower scaling of 0.8833; hence, it is the least efficient model assumed.

Figure 6 indicates that the performance of the GBoost algorithm exceeds that of the other algorithms with respect to the AUC values. In addition, the model performs perfect classification in a number of classes, which means that it is well optimized for the task and it is capable of determining the majority of these diseases. The grade of performance for AS is slightly lowers with highly accurate. The slight proximity of the curves to the lower zones of false positive rates suggests that there could be some problems in the early detection stage or in cases where some of the signs presented by one of the diseases are also exhibited by another illness. However, the results prove that there is a place for the GBoost algorithm with respect to the classification of these diseases, as it strongly varies in the diagnosis of ReA, SS, and SLE and differentiates healthy (normal) people from patients suffering from these autoimmune conditions.

Figure 7 shows the performance of the precision–recall (PR) curves of the GBoost algorithm in classes of rheumatic and autoimmune diseases. The accuracy or the PR score is helpful for skewed datasets, as it considers the relationship between precision and recall. The mean average precision for the entire GBoost algorithm is 0.95, which means that the model can dispose of positive predictions across most classes while making most positive predictions correct across all classes. The plot also includes lines of F1-scores, which illustrate the dependence of precision and recall on the threshold. The PR of SLE is 1.00, which essentially means that the disease diagnosis is perfect and has both recall and precision of nearly unity the entire time, that is, no false allegations or misses. The GBoost algorithm in the majority of classes provides good prediction of the majority of diseases with high precision and high recall, with the exception of AS, and this condition is likely related to symptoms of other diseases.

Evaluation analysis

This section focuses on two important key aspects, XAI and comparative analysis, to evaluate our framework. These keys show that our framework not only delivers accurate predictions but also provides transparency and interpretability, alongside a comprehensive comparison of various machine learning techniques.

XAI

XAI is an area of work that has rapidly expanded in recent years and focuses on improving the interpretability of ML models. The main objective of XAI is to design models in such a way that the reasons for its outputs can be clarified to the human being interfaced with the system. This can involve providing reasons for individual predictions or explaining why a model prediction is made in general and/or what features are relevant to the model’s outputs. In this work, use Local interpretable model-agnostic explanations (LIME) are a tool that makes ML models more understandable by providing default explanations of the decisions made by them.

In this work, three case studies were used to demonstrate the superiority and rationality of the proposed model. In Table 19, the following clinical characteristics are provided for the correct diagnosis of patients: ReA (case 1), PA (case 2), and N (case 3). In general, the characteristic weights reveal that there are different and sometimes conflicting roles of clinical markers in differentiating the pathological states of patients, particularly inflammatory markers and autoantibodies in certain cases, whereas they are less important in others.

Table 19:

Three cases from the testing set.

Feature	Case 1		Case 2		Case 3
	Weight	Value	Weight	Value	Weight	Value
Age	0.0047	66	0.0054	71	−0.0017	25
Gender	−0.0161	1	−0.0160	1	0.0169	0
ESR	−0.0214	24.003	−0.0691	33.538	0.1527	4.494
CRP	−0.0282	13.539	−0.0154	25.645	0.0025	13.299
RF	0.0506	14.674	0.0404	3.118	0.0448	10.020
Anti-CCP	0.0388	6.449	0.0388	13.132	0.0422	2.257
HLA-B27	−0.0243	1	−0.0250	1	−0.0228	1
ANA	−0.0669	1	−0.0701	1	−0.0723	1
Anti-Ro	−0.0331	0	0.0335	0	−0.0355	1
Anti-La	−0.0302	1	−0.0324	1	−0.0332	1
Anti-dsDNA	0.0166	0	−0.0148	1	0.0149	0
Anti-Sm	−0.0061	1	0.0049	0	−0.0031	1
C3	−0.0329	111.332	−0.0329	107.243	0.0173	147.178
C4	0.0041	39.314	0.0050	64.989	0.0052	38.807
True class	ReA		PA		N

DOI: 10.7717/peerj-cs.3096/table-19

Note:

For gender feature, Female will encoded as 0 and Male will be as 1. For the characteristics (HLA-B27, ANA, Anti-Ro, Anti-La, Anti-dsDNA and Anti-Sm) will Encode (Negative = 0 and Positive = 1).

Figure 8 shows the predicted probabilities of different diagnostic possibilities for a given medical case via the three exploratory diagnostic. In Case 1, the model predicts ReA with 100% confidence, ruling out AS, N, PA, etc., all of which are 0.00. The development of probabilities in Case 2 indicates that (83%) is the most likely diagnosis of PA, reducing the chances of AS to a further low (17%) and completely ruling out other states. In Case 3, the model has a mean prediction of (53%) SS and (46%) N. Overall, Case 1 produced the highest level of diagnostic accuracy, Case 2 the second highest, and Case 3 the broadest range of diagnostic potential.

Figure 9 shows the determination of the diagnosis in three cases via the GBoost algorithm based on clinical markers, where each of the three cases has a different threshold for one or more of the features as follows:

• –

  In case 1 (ReA), the model predetermines the ESR as an important parameter. When the ESR is less than 26.03, the model assumes lower values of CRP, indicating moderate inflammation. Furthermore, the absence of RF, Anti-CCP, Anti-dsDNA, and Anti-sm, for example, increases the strength of the model in diagnosing ReA. The importance of these features is illustrated by the size of the bars corresponding to them, that is, how much each is involve into the decision.

• –

  In case 2 (PA), the model differentiates patients above the CRP value of 19.06, where more inflammation is likely to be present in the PA. This value is also below the RF index of 11.27 and reflects the blue presence of Anti-CCP. These indicators, when added to the exclusion of HLA-B27 and other biomarkers of autoimmune diseases, are reasonable for the model to consider PA as the most likely diagnosis.

• –

  In case (SS), the model prioritizes the use of an ESR as one of the first steps, considering the erythrocyte sedimentation rate values to be of maximum $\le$14.19 as a very good predictor. Positive findings for anti-Ro and anti-La antibodies push the diagnosis even further toward SS. Lower CRP values and high C3 levels help reduce the number of diagnoses and exclude others.

Figure 10 shows the interpretability of GBoost algorithm concerning the classification of patients in (N or Not N), which is based on the clinical factors in use. Among Fig. 9 showed features, the most satisfactory feature for pushing the classification toward N is the ESR, as low ESR scores $\le$14.19 help improve the likelihood of a normal report, as evidenced by a relatively high importance score of 0.15. Additionally, reasonably important features are the factors RF and Anti-CCP, with cutoff values of $\le$11.27 and $\le$13.82, respectively. A normal diagnostic status is also supported by C3 levels within the range of 132.44 to 157.09, which score 0.02 in importance. However, the patient’s sex itself made a slight contribution to the N classification, with a score of 0.02. Anti-dsDNA levels $\le$0.06 are also a long way from promoting normalcy. C4 levels $\le$38.81 favor an N diagnosis with a very subtle importance score of 0.01. In contrast, high levels of ANA (>1.00) greatly help the classification of Not N so that it has an importance score of 0.07, which denotes an abnormal condition and explains why it is very useful.

Figure 11 shows how globally XAI important features are for the prediction of autoimmune and rheumatic diseases. The ESR is the most important factor, with an importance score of approximately 0.40, due to its close relationship with inflammation in autoimmune diseases. The second feature is RF, which is necessary for the diagnosis of rheumatism and it is where the largest gap occurs, followed by C3. C3 is important in diseases such as lupus because complement levels are altered under these conditions; hence, C3 becomes critical. Other second-level technical features, such as ANA, anti-Ro, anti-La, anti-dsDNA, and anti-SM autoantibodies, contribute less strongly but they are important in the diagnosis of certain autoimmune diseases, such as AS and SLE. Hypertension, which is driven by age and sex, is not critical, at least in measuring disease risk or poor health outcomes. Since clinical and laboratory biomarkers, especially those of inflammation and active tubular autoantibodies, are more common than are demographic biomarkers. This gives credence to two of our clinical tests in ESR, RF, C3, anti-CCP, and CRP, as first-line diagnostic tools for autoimmune diseases.

The global feature importance using GBoost model as shows in Fig. 11 reveals that ESR and RF are dominant predictors, aligning with current diagnostic guidelines where these markers are first-line tests for inflammatory arthritis. However, the model additionally highlights C3 as equally important as RF (scores 0.35 vs. 0.34), suggesting lab workflows might prioritize C3 testing more frequently than currently done a potentially practice changing insight.

Comparative analysis

The datasets in the ML approach often encounter the common challenge of considering class imbalance, where models overlearn a particular class and neglect minority classes, which makes their performance low on such classes. To overcome this, data balancing techniques, which involve oversampling, are applied. It is crucial to assess the effectiveness of such techniques by performing a performance comparison on the initial set, which was parched, and the resulting balanced set. The direct difference between the balanced and unbalanced metrics can be expressed as the absolute difference (AD), as shown in Eq. (18). A positive value means that there is an improvement within the bias correction method, and a negative value means that there is a reduction within the method. The percentage difference (PD) is used to represent the improvement or degradation of the performance after the dataset has been balanced over the values of the unbalanced performance and it is represented mathematically in Eq. (19).

(18) $$AD=Balance\mathrm{\_}dataset\text{-}Unbalance\mathrm{\_}dataset$$ (19) $$PD=\frac{AD}{Unbalanced\mathrm{\_}dataset}\ast 100.$$

Table 20 shows the performance of the different metrics of the ML models with a greater focus on the AD and PD between the balanced and unbalanced datasets. Models like KNN and SVM showed significant improvements, with KNN achieving a 15.14% precision PD, while RF and ET displayed minimal changes. The stochastic behavior of SGD was evident, with metrics fluctuating significantly, indicating inconsistent adaptation to dataset balancing. AdaBoost and NB demonstrated moderate enhancements, with PD values ranging from approximately 4% to 8% when handling highly enhanced dataset imbalances. The present results highlight that appropriately balancing the training datasets helps enhance model performance. However, the extent of the improvements is dependent on the models.

Table 20:

Performance comparison of ML models for Unbalanced and Balanced dataset, showing AD and PD for seven evaluation metrics.

Model	Precision		Recall		F-score		Kappa		HL		MCC		Accuracy
	AD	PD	AD	PD	AD	PD	AD	PD	AD	PD	AD	PD	AD	PD
XGBoost	3.304	3.97%	5.192	6.41%	4.458	5.45%	4.936	6.27%	−3.892	−21.91%	4.914	6.24%	3.892	4.73%
KNN	9.468	15.14%	9.107	14.69%	8.877	14.37%	7.099	12.07%	−5.554	−16.20%	7.490	12.70%	5.554	8.45%
DT	4.785	6.07%	4.746	6.02%	4.792	6.08%	4.321	5.72%	−3.375	−16.52%	4.314	5.71%	3.375	4.24%
RF	0.203	0.24%	5.050	6.27%	3.317	4.03%	2.891	3.61%	−2.136	−12.97%	2.886	3.59%	2.136	2.56%
SVM	12.333	19.64%	9.101	14.08%	9.780	15.40%	3.749	5.70%	−2.610	−9.21%	4.027	6.10%	2.610	3.64%
ANN	5.358	6.87%	5.361	6.98%	4.924	6.37%	2.870	3.80%	−2.149	−10.54%	3.155	4.18%	2.149	2.70%
SGD	2.202	3.06%	−1.956	−2.75%	−2.641	−3.85%	−3.728	−5.60%	3.296	11.81%	−3.348	−4.94%	−3.296	−4.57%
GBoost	0.803	0.93%	4.744	5.78%	3.259	3.90%	3.846	4.77%	−2.967	−18.45%	3.844	4.76%	2.967	3.54%
LGBM	2.755	3.30%	4.565	5.61%	3.819	4.65%	4.391	5.55%	−3.291	−19.04%	4.327	5.46%	3.291	3.98%
AdaBoost	3.539	4.45%	5.046	6.48%	4.546	5.80%	4.275	5.68%	−3.613	−17.44%	4.157	5.51%	3.613	4.56%
NB	4.117	6.17%	5.143	7.88%	4.723	7.17%	5.359	8.49%	−4.490	−12.93%	5.405	8.55%	4.490	6.88%
ET	−0.179	−0.21%	2.689	3.22%	1.355	1.59%	2.895	3.49%	−1.964	−11.83%	2.805	3.38%	1.964	2.36%

DOI: 10.7717/peerj-cs.3096/table-20

Figure 12 presents a comparative analysis of different ML procedures with respect to unbalanced and balanced datasets. The performance of most models generally increases when the datasets are balanced, especially for the minority class, and the precision-recall trade-off increases with increasing MCC and F-score for balanced datasets.

The least effective algorithm in this study is AdaBoost which recorded a very low accuracy ranging from 34% to 37%. This dismal outcome can be linked to a number of reasons: (i) highly correlated of features which was likely to confuse the weak learners and hence perform poorly, (ii) a restricted number of weak learners (set to 50 in this investigation) which may prove to be too few to learn the patterns that exist in the data, and (iii) the inclusion of noise in the form of irrelevant features in the dataset.

Table 21 achieve marginally higher accuracy (e.g., Wu et al., 2025; Wang et al., 2023; Jia et al., 2024), our study distinguishes itself in several key areas that make it more robust, generalizable, and clinically relevant. Most prior studies are limited to binary classification (e.g., RA vs. non-RA), whereas our proposed framework simultaneously classifies six distinct rheumatic and autoimmune diseases with normal case, providing a more realistic and useful tool for differential diagnosis. Our dataset comprises 12,085 real patient records, significantly surpassing synthetic datasets (Martins et al., 2023) and small image dataset (Nouira et al., 2024). This enhances the model’s applicability in actual clinical settings. Our proposed work use FDOSM which allows to choose best ML model from 12. We integrate LIME to provide local interpretability of predictions, a critical factor for adoption in healthcare environments where transparency is essential. With an accuracy of 86.89%, our proposed maintains high predictive performance while offering explainability and multi-class support, balancing performance with interpretability something rarely achieved together in the literature.

Table 21:

Summary of studies on disease diagnosis with ML methods.

Ref	Disease(s)	Method	Data size	Accuracy (%)
Martins et al. (2023)	ALS, MS, SLE, CD, AIH	LR, NB, SVM	Synthetic (1,000)	84
Wang et al. (2023)	LN	XGB	1,467	99
Park (2024)	SLE	GB, RF	24,990	88
Pham et al. (2024)	RA	XGB	5,600	81.6
Jia et al. (2024)	AS	SCJAYA-FKNN	UCI + AS	99.2
Nouira et al. (2024)	SS	CNN	225	99
Qiu et al. (2025)	RA	LightGBM	SHARE survey	90.2
Wu et al. (2025)	RA	RF, LR	2,106 + external	96.2
Our work	6 diseases with normal	12 ML models + FDOSM + LIME	12,085	86.89

DOI: 10.7717/peerj-cs.3096/table-21

This work addresses the pressing challenge of accurately diagnosing complex and overlapping autoimmune and rheumatic diseases, which are often misclassified due to the lack of specific biomarkers and the subjective interpretation of clinical signs. Traditional ML models may offer high accuracy, but they often lack transparency and reliability in critical domains like healthcare.

To tackle this, we proposed a hybrid framework that combines multiple ML models and FDOSM with XAI. The FDOSM enables a multi-criteria evaluation that goes beyond accuracy to include precision, recall, kappa, and other metrics, incorporating expert judgment in a rational and fair ranking process. Meanwhile, XAI enhances trust and clinical adoption by offering global and local explanations for each diagnosis, thereby supporting transparent and evidence-based medical decisions.

The theoretical contribution of this study lies in demonstrating how FDOSM and XAI can bridge the gap between predictive accuracy and interpretability. Our results show that GBoost algorithm is top performer according to FDOSM. Also, more interpretable when coupled with XAI. This dual strength significantly contributes to the emerging field of trustworthy AI in healthcare.

This work also has limitations such as; (i) The dataset was exclusively collected from Iraq (2019–2024), which may include (patient demographics, genetic, and healthcare protocols), (ii) investigating how the use of deep learning-based imputation may influence the results by remaining more robust, (iii) investigating whether any of the models, such as feature selection or feature reduction, could provide insight into which variables are more valuable in predicting the results of the model, and (iv) the absence of adversarial testing: the models were not evaluated against adversarial attacks. There are a number of possible fixes that might be used to overcome these limitations.

• –

  Future studies should incorporate datasets from multi-center, international cohorts to enhance generalizability.

• –

  Sophistication imputation methods, including denoising autoencoders are used to handle complex missing data patterns. These strategies may produce better imputed values and capture underlying data distributions.

• –

  To identify the most important variables for model prediction, t-SNE, PCA, or recursive feature elimination (RFE) is used.

• –

  Adversarial testing techniques, such as Projected Gradient Descent (PGD) or the Fast Gradient Sign Method (FGSM), are used to assess the model’s resistance to adversarial attacks.