A Bayesian Inference Based Computational Tool for Parametric and Nonparametric Medical Diagnosis

Medical diagnosis is the basis for treatment and management decisions in healthcare. Conventional methods for medical diagnosis commonly use established clinical criteria and fixed numerical thresholds. The limitations of such an approach may result in a failure to capture the intricate relations between diagnostic tests and the varying prevalence of diseases. To explore this further, we have developed a freely available specialized computational tool that employs Bayesian inference to calculate the posterior probability of disease diagnosis. This novel software comprises of three distinct modules, each designed to allow users to define and compare parametric and nonparametric distributions effectively. The tool is equipped to analyze datasets generated from two separate diagnostic tests, each performed on both diseased and nondiseased populations. We demonstrate the utility of this software by analyzing fasting plasma glucose, and glycated hemoglobin A1c data from the National Health and Nutrition Examination Survey. Our results are validated using the oral glucose tolerance test as a reference standard, and we explore both parametric and nonparametric distribution models for the Bayesian diagnosis of diabetes mellitus.


Introduction
Medical diagnosis is a critical process of accurately identifying pathological conditions in patients. The term "diagnosis" has its etymological origins in the ancient Greek word "διάγνωσις", signifying "discernment" [1]. Traditionally, diagnostic tests are used to divide individuals into two principal categories: those who are afflicted with a specific disease and those who are not. Notably, the probability distributions associated with quantitative diagnostic test outcomes often demonstrate some overlap between the diseased and nondiseased groups. To address this, numerical diagnostic thresholds or cut-off points have been formulated to provide a binary classification of these test outcomes [2]. Nevertheless, this introduces a certain measure of uncertainty into the diagnostic accuracy of those tests [3]. This dichotomous method represents a significant shift in medical decision-making by linking a continuum of evidence to binary clinical decisions such as to treat or not to treat [4].
Despite the evident efficiency of traditional diagnostic methods, they sometimes fail to capture the complexity and heterogeneity of disease presentations across diverse populations [5]. To address these limitations, our research focuses on implementing Bayesian inference to calculate the posterior probabilities associated with disease diagnosis [6][7][8][9]. Within this Bayesian paradigm, prior probabilities of disease are integrated with distributions of diagnostic measurands in both diseased and nondiseased populations. This approach enables the evaluation of the information conveyed via diagnostic measurements and the combination of data from multiple diagnostic tests, which may improve diagnostic accuracy and precision while introducing flexibility, adaptability, and versatility into the diagnostic process [10]. Furthermore, the Bayesian approach extends its utility beyond the medical field by offering a robust framework for quantifying uncertainty in various domains, thereby enriching its applicability in both diagnostic and prognostic contexts [11,12].
A considerable challenge in integrating Bayesian inference into medical diagnosis is the limited availability of literature detailing the statistical distributions of diagnostic variables in both pathological and non-pathological states [13].
The ubiquitous application of the normal distribution in clinical laboratory indicators is due, in part, to its mathematical simplicity, the foundational Central Limit Theorem, and a rich collection of statistical methods designed for Gaussian data [14]. However, the universal applicability of the normal distribution is subject to critique, especially when dealing with clinical measurands that exhibit skewness, bimodality, or multimodality [15]. Hence, while the normal distribution remains invaluable in statistical modeling, critical evaluation of its appropriateness for specific diagnostic measurands is necessary. This evaluation should be accompanied by an openness to adopt alternative statistical distributions when needed [16].
This foundational data is crucial for Bayesian inference, establishing the essential context against which new diagnostic measurements can be compared. The absence of such normative data could potentially compromise the reliability and validity of Bayesian diagnostic methods.
To address the complex issues related to Bayesian diagnosis and the selection of appropriate statistical distributions for diagnostic variables, we have developed the Bayesian Diagnosis program, an interactive software tool programmed in the Wolfram Language. This tool allows users to explore and compare both parametric and nonparametric distributions to calculate posterior probabilities for disease. It is designed to analyze datasets of measurements of two distinct diagnostic tests, performed on both diseased and nondiseased populations.

The Program
Bayesian Diagnosis was developed using Wolfram Mathematica ® Ver. 13.3 (Wolfram Research, Inc., Champaign, IL, USA (2023)). This interactive program consists of three primary modules with eighteen submodules. It allows the calculation, plotting, and comparison of Bayesian posterior probabilities of disease for two diagnostic tests, assuming two sets of alternative parametric and nonparametric distributions of the measurements of those tests in diseased and nondiseased populations (refer to Figure 1 and to Supplementary File S1). It is freely available as a Wolfram Notebook (.nb) (Supplementary File: BayesianDiagnosis.nb). It can be run on Wolfram Player ® or Wolfram Mathematica ® (refer to Appendix B).

Datasets
As datasets are considered tuples of data. Although the program includes four datasets of measurements, one for each diagnostic test, applied to a diseased and a nondiseased population, these can be replaced by other appropriate datasets selected by the user (refer to Appendix B). Therefore, it can be used for any combination of diagnostic tests and diseases.

Bayesian Diagnostic Approach
The Bayesian diagnostic approach is a cornerstone in statistical inference and particularly useful in medical diagnosis [6,17,18]. The approach relies on Bayes' theorem [7]. For effective implementation of the Bayesian diagnostic method, knowledge concerning the statistical distributions of the measurements of the diagnostic tests is essential [14].
Bayes theorem is presented in Appendix A.

Datasets
As datasets are considered tuples of data. Although the program includes four datasets of measurements, one for each diagnostic test, applied to a diseased and a nondiseased population, these can be replaced by other appropriate datasets selected by the user (refer to Appendix B). Therefore, it can be used for any combination of diagnostic tests and diseases.

Bayesian Diagnostic Approach
The Bayesian diagnostic approach is a cornerstone in statistical inference and particularly useful in medical diagnosis [6,17,18]. The approach relies on Bayes' theorem [7]. For effective implementation of the Bayesian diagnostic method, knowledge concerning the statistical distributions of the measurements of the diagnostic tests is essential [14].
Bayes theorem is presented in Appendix A.

Parametric Distributions
Parametric statistics assume that dataset data comes from a population that can be adequately modeled with a probability distribution that has a fixed set of parameters [19]. The parametric distributions provided by the program are the following:

Parametric Distributions
Parametric statistics assume that dataset data comes from a population that can be adequately modeled with a probability distribution that has a fixed set of parameters [19]. The parametric distributions provided by the program are the following:

Copula Distributions
The copula distributions of the program are bivariate, with a bivariate normal distribution with correlation ρ as kernel, and univariate normal, lognormal and gamma marginals.
The probability density functions (PDFs) of the parametric distributions are mathematically defined in Appendix A.

Nonparametric Distributions
Conversely, nonparametric models were also employed, which do not make a priori assumptions about the distribution's mathematical form [20]. These are particularly useful for exploratory data analysis and are implemented as shown in Appendix A.

Histograms
A histogram is the graphical representation of the distribution of a dataset as a series of bins.
The program plots histograms of the provided datasets.

Kernel Density Estimators (KDEs)
In contrast to histograms, a KDE generates a continuous and smooth estimate of the underlying PDF by summing the contributions of kernel functions centered at each data point.
KDEs offer a flexible nonparametric approach to density estimation, allowing for a better representation of the underlying distribution of the data.
The program provides univariate and bivariate Gaussian KDEs. The bivariate KDEs use radial-type kernels.

Interface of the Program
The Bayesian Diagnosis program is equipped with an intuitive tabbed user interface (refer to Figure 2). This design facilitates seamless navigation through its various modules and submodules. Users have the capability to input and modify a range of parameters, including prior probabilities and measurement parameters. Additionally, the interface allows for the selection of both parametric distributions and KDEs pertinent to medical diagnosis (refer to Appendix C and Supplementary File S1).

Input Parameters Prior Probability
The user initiates the diagnostic evaluation by specifying the prior probability of disease occurrence in the population under study. This serves as a foundational measure for subsequent analyses.

Parametric Distributions
To facilitate a diagnostic model, the program allows for the definition of various parametric distributions for both the diseased and nondiseased populations across two diagnostic tests.

1.
Distribution Selection: The user selects the type of distribution from a predefined list: Gamma Distribution.

2.
Statistical Parameters: For each chosen distribution, the user defines the mean µ and standard deviation σ of the measurand in the respective population.

3.
Correlation Coefficients: The user specifies the correlation coefficients ρ between the measurands of the first and second diagnostic tests for both diseased and nondiseased populations.

KDEs
Alternatively, the user can opt to define the KDEs for the measurands in both diseased and nondiseased populations across the two tests:

1.
Bandwidth Parameter: For each KDE, the user defines the bandwidth parameter h.

2.
Correlation Coefficients: As with parametric distributions, the user defines the correlation coefficients ρ between the measurands of the two diagnostic tests.

Input Parameters
Prior Probability The user initiates the diagnostic evaluation by specifying the prior probability of disease occurrence in the population under study. This serves as a foundational measure for subsequent analyses.

Parametric Distributions
To facilitate a diagnostic model, the program allows for the definition of various parametric distributions for both the diseased and nondiseased populations across two diagnostic tests.

KDEs
Alternatively, the user can opt to define the KDEs for the measurands in both diseased and nondiseased populations across the two tests: 1. Bandwidth Parameter: For each KDE, the user defines the bandwidth parameter ℎ. 2. Correlation Coefficients: As with parametric distributions, the user defines the correlation coefficients ρ between the measurands of the two diagnostic tests.

Output Specifications Visualizations
The program generates a series of plots designed to elucidate various diagnostic metrics and statistics:

1.
Posterior Probability of Disease: Plots are generated to show the posterior probability of disease for each measurand and their combination.

2.
PDFs: Univariate PDFs for each measurand and the bivariate PDF of their combination are plotted. An option to overlay histograms on these plots is also provided. 3.
Probability-Probability (P-P) Plots: Similar to Q-Q plots, P-P plots are generated for further assessment of the distribution of each measurand [21].
The descriptions of the Q-Q and P-P plots are presented in the Supplementary File S2.

1.
Population Statistics: The program tabulates key statistical metrics such as mean, median, standard deviation, skewness, kurtosis, and prior probability for each userdefined distribution and dataset. For each bivariate distribution of the two measurands in diseased and nondiseased populations, the correlation coefficients are calculated and displayed.

2.
Posterior Disease Probabilities: For a user-defined pair of test measurement values, the program computes and presents the posterior probabilities for disease for each measurand and their combination.
By providing this comprehensive set of input parameters and output specifications, the program offers a robust platform for exploring the Bayesian diagnosis of disease using either parametric distributions or KDEs of medical diagnostic measurands.

Illustrative Application
To demonstrate the application of the program, fasting plasma glucose (FPG) was used as the first measurand and glycated hemoglobin A1c (HbA1c) as the second measurand for Bayesian diagnosis of diabetes mellitus. The oral glucose tolerance test (OGTT) was used as the reference diagnostic method. A diagnosis of diabetes was confirmed if the plasma glucose (PG) value was equal to or exceeded 200 mg/dL, measured two hours after oral administration of 75 g of glucose [22], during an OGTT (2-h PG). It is noteworthy that the study population was confined to individuals aged between 40 and 60 years, a decision informed by the well-documented strong correlation between age and the prevalence of diabetes [23].
National Health and Nutrition Examination Survey (NHANES) data from participants was retrieved for the period from 2005 to 2016 [24] (n = 60,936). NHANES is a series of studies designed to evaluate the health and nutritional status of adults and children in the United States.
Participants with a 2-h PG measurement of ≥ 200 mg/dL were considered diabetic (n = 211).
Descriptive statistics, including the mean, median, and standard deviation, were computed for each dataset and correlation coefficients for their combination. Univariate distributions were employed to model the distributions of FPG and HbA1c and bivariate distributions to model the joint distribution of FPG and HbA1c. Likelihoods and posterior probabilities were estimated for FPG, HbA1c and their combination.
The prior probability of diabetes was estimated as follows: The statistics of the dataset are presented in Table 1.

Results
Using the settings of Table 2, the program generated the plots of     Table  2.  Table 2.
tions, and of the respective datasets (NHANES datasets), with the settings of the program in Table  2. Figure 15 shows a table of prior and posterior probabilities for disease (diabetes) for values of FPG equal to 126 mg/dL and of HbA1c equal to 6.5%, the established thresholds of the two measurands for the diagnosis of diabetes [22], assuming parametric and KDE distributions.  Table 2.  Table 2.

Reevaluation of Traditional Diagnostic Methods
The KDEs smoothing bandwidth was set to double that given with Silverman's rule of thumb [26,27]. Figures 3 and 4 show the plots of the posterior probability of diabetes versus FPG and HbA1c, respectively. The curves of the parametric distributions are smooth double sigmoidal, while the curves of the nonparametric distributions are multimodal. Figure 5 shows the plot of the posterior probability of diabetes versus FPG and HbA1c combined. The surface of the parametric distribution is smooth, while the surface of the nonparametric distribution is multimodal. Figures 6-9 show the PDF of FPG and HbA1c in diabetic and nondiabetic patients and the histograms of the respective NHANES datasets. It is visually evident that the nonparametric distributions fit the datasets better, especially in diabetic patients. Figures 10-13 show the Q-Q plots of the parametric and KDE distributions of FPG and HbA1c in diabetic and nondiabetic patients versus the respective NHANES datasets. The plots show clearly that the nonparametric distributions fit the datasets better, especially in diabetic patients. Figure 14 shows a table with the descriptive statistics of FPG and HbA1c in diabetic patients and nondiabetic patients, assuming parametric and KDE distributions, and of the respective NHANES datasets. The data, including the loglikehood values, supports the hypothesis that the nonparametric distributions fit the datasets better, especially in diabetic patients. Figure 15 shows a table of prior and posterior probabilities for disease (diabetes) for values of FPG equal to 126 mg/dL and of HbA1c equal to 6.5%, the established thresholds of the two measurands for the diagnosis of diabetes [22], assuming parametric and KDE distributions.

Reevaluation of Traditional Diagnostic Methods
The findings of the present study highlight the importance of considering incorporating Bayesian methods in medical diagnosis and management. Conventional approaches based on rigid diagnostic criteria are often unable to account for the intricate relationships between disease pathology and diagnostic procedures thus limiting personalized patient care options. [28]. In stark contrast, Bayesian methodologies offer a framework that enhances diagnostic precision through a more comprehensive probabilistic assessment [5]. This Bayesian foundation, therefore, could serve as an enabler for tailored medical interventions, echoing similar arguments in existing literature advocating for individualized medicine [29].
The study population was confined to individuals aged between 40 and 60 years. This restriction allowed for a more homogeneous prior probability, thereby reducing the impact of age-specific variations in the prevalence of the condition under study.
Even though the KDEs from our illustrative application, as parameterized in Table 2, provide only an approximate fit to the NHANES datasets for FPG and HbA1c measurements, the posterior probabilities for diabetes delineated in Figure 15 suggest a limited concordance between the classification criteria of diabetes derived from the OGTT, HbA1c, and FPG tests [22], as found previously in existing literature [30].

Challenges and Considerations in Bayesian Analysis for Disease Diagnosis
Despite the evident merits of Bayesian analytics in medical diagnostics, it is paramount to address the intrinsic challenges associated with this methodological shift. One such issue resides in the limited availability of scholarly publications that provide a comprehensive statistical exploration of the measurands in both the diseased and nondiseased populations [31].

1.
Over-dependence on Prior Probabilities: The scarcity of empirically derived distributions amplifies reliance on prior probabilities, thereby inducing distortions in the calculation of posterior probabilities. This could result in suboptimal clinical judgments and potentially inaccurate diagnoses [32].

2.
Elevated Uncertainty: Insufficient data contributes to broader confidence intervals in the computed posterior probabilities, which, in turn, could exacerbate clinical indecisiveness [33].

3.
Risk of Bias: The introduction of systemic bias due to unrepresentative datasets could compromise the fidelity of Bayesian calculations [7].

4.
Imperative for Collaborative Research: More coordinated research is needed, including multi-center studies, meta-analyses, and open-access databases-to accumulate and disseminate data essential for effective Bayesian diagnosis [34].

5.
Exploration of Alternative Methodologies: Given the lack of comprehensive data, the utility of combining Bayesian methods with other statistical and computational techniques or diagnostic modalities becomes increasingly pertinent [35,36].

Parametric Versus Nonparametric Bayesian Models
In the context of diagnosing diabetes mellitus through FPG and HbA1c levels, our computational tool revealed that nonparametric Bayesian models typically produce a better fit to data distributions, corroborating existing literature that emphasizes the robustness of nonparametric techniques in capturing complex data distributions [26,37].

Multimodal Versus Double Sigmoidal Bayesian Probability of Disease Curve
The nonparametric Bayesian probabilities for disease exhibited multimodal patterns, in contrast to the bimodal, double sigmoidal curves generated by parametric models.

Multimodal Curve
Potential Causes: (a) Complex Pathophysiology: Multiple etiological pathways may influence the same measurand in divergent ranges, adding layers of complexity to diagnostic processes [13]. (b) Diagnostic Confounders: External variables affecting the measurand could compromise its efficacy as a standalone diagnostic criterion [38]. (c) Population Subgroups: The existence of demographically or genetically distinct subgroups within the studied population could also account for the observed multimodality [39].
(d) Statistical Artifacts: Demographically or genetically distinct subgroups may be a factor contributing to observed multimodal distributions [39].
Theoretical Implications: Multimodal distributions present a clinical conundrum, compelling healthcare providers to potentially employ additional diagnostic tests or even alternative methodologies [13].

Double Sigmoidal Curve
A curve composed of two mirrored sigmoid functions, one delineating the probability behavior for lower measurand values and the other for higher values-presents a compelling intricacy in the realm of diagnostic statistics and medical decision-making. Interpretation (a) Two Zones of Risk: Such a curve suggests that the risk of the disease is heightened both at low and high extremes of the measurand but reduced in a middle "safe zone." (b) Multifactorial Etiology: This might reflect a situation where both deficiency and excess of a particular biological factor contribute to disease risk. For example, both low and elevated levels of hormones may pose challenges to physiological homeostasis.

Clinical and Diagnostic Implications
(a) Threshold Decision-making: Unlike a single sigmoid curve, where one threshold may be adequate for diagnosis, the double sigmoid may necessitate multiple thresholds, defining a "safe zone" for the measurand. (b) Treatment Strategies: Clinicians must be cautious when intervening based on such a measurand, as moving the measurand too far in either direction could heighten risk. (c) Population Stratification: This curve shape might imply that different sub-populations or disease subtypes could be better distinguished with additional tests or measurements.

Shortcomings of This Study
The main shortcomings of this study were the following: 1.
The OGTT was used as a reference diagnostic method for diabetes mellitus. The diagnostic threshold for 2-h PG was established in relation to the risk of diabetic retinopathy, a microvascular complication of diabetes mellitus [40]. However, glucose tolerance is influenced by complex interactions of factors, both physiological and environmental, which pose significant implications for clinical diagnosis and research. The considerations that could affect glucose tolerance and, therefore, the interpretation of the 2-h PG measurement, include the following: (a) Age and Gender: Age and gender are significant variables in glucose tolerance. Insulin sensitivity often decreases with age, resulting in higher 2-h PG levels [41]. Gender differences, particularly related to hormonal changes in females, could also affect glucose metabolism [42].
Diurnal Variability: Glucose tolerance is subject to diurnal variation, which could affect the 2-h PG test outcomes. Insulin sensitivity is generally higher in the morning than in the evening [43]. (c) Physical Activity: Exercise improves insulin sensitivity and therefore could affect glucose tolerance tests. The timing and intensity of physical activity could have a direct influence on the 2-h PG results [44]. (d) Dietary Patterns: Short-term and long-term dietary habits, including the macronutrient composition of the diet, may alter the body's glucose and insulin response [45]. (e) Stress and Emotional States: The acute stress response includes a transient rise in glucose levels as a result of catecholamine release, potentially affecting the 2-h PG test [46]. (f) Medications: Certain medications like corticosteroids, antipsychotics, and diuretics affect glucose metabolism, thereby influencing 2-h PG test outcomes [47].
(g) Genetic Factors: Genetic predispositions influence glucose tolerance, and not accounting for this introduce variability in the 2-h PG test [48].

2.
The lognormal distributions and the KDE, as parameterized in Table 2, fitted only approximately to the NHANES datasets of FPG and HbA1c measurements. It is well known that biological measurands, such as FPG and HbA1c, may not follow textbook statistical distributions like normal or lognormal distributions. Numerous papers have noted the skewness or kurtosis in the distribution of metabolic variables, urging the use of flexible statistical models [49,50]. The program presented in this work provides 29 different types of parametric and nonparametric plots. None of the above-mentioned programs provide this range of plots without advanced statistical programming.

Conclusions and Future Directions
The intricacies of the double-sigmoid curve and multimodal distributions introduce a new frontier in personalizing healthcare provision. While smoother relationships between measurements and Bayesian probability facilitate clinical interpretability, multimodal distributions might serve as sentinel indicators of underlying complexities or methodological shortcomings, thus providing a useful tool in the field of medical diagnosis.
As a pivotal next step, future research should aim to validate the utility and reliability of the Bayesian inference based method applied in this study through real-world clinical trials, in addition to extending its application to include more diagnostic modalities. The aim is to combine this approach with existing clinical protocols, thereby optimizing the diagnostic precision and consequently improving patient outcomes.
In addition to its potential for clinical applications, the computational tool developed for this study could hold considerable promise as an educational and research adjunct. By facilitating the analysis of Bayesian probabilities in disease diagnosis, it serves as an invaluable resource for both medical practitioners in training and experienced researchers in the field. Its modular design and user-friendly interface make it easily adaptable to various research settings and educational curricula, thereby accelerating the adoption and dissemination of Bayesian approaches in medical statistics and diagnostics.

Informed Consent Statement:
Written consent was obtained from each subject participating in the survey.

Data Availability Statement:
The data presented in this study is available at https://wwwn.cdc. gov/nchs/nhanes/default.aspx (accessed on 4 September 2023).

Conflicts of Interest:
The authors declare no conflict of interest.  Therefore, for the possibly multivariate parameter θ:

Abbreviations
where L D (θ|z) and f D (z|θ) denote the likelihood function and the PDF in the presence of the disease, while L D (z|θ) and f D (z|θ) denote the respective functions in the absence of the disease. The univariate normal distribution or Gaussian distribution is a continuous probability distribution of a random variable X. The general form of its PDF is: where the parameter µ is the mean of X, while σ is its standard deviation [52].
The bivariate normal distribution or Gaussian distribution is a continuous probability distribution of two normally distributed random variables X and Y. The general form of its PDF is: where µ X and µ Y are the means of X and Y, σ X and σ Y are their standard deviations, and ρ their correlation coefficient [52].
The univariate lognormal distribution is a continuous probability distribution of a random variable X whose logarithm is normally distributed. The general form of its PDF is: where m is the mean and s the standard deviation of ln(X) [52].
If µ and σ are the mean and the standard deviation of X, we have: Therefore, The bivariate lognormal distribution is a continuous probability distribution of two lognormally distributed variables X and Y. If m X and m Y are the means of ln(X) and ln(Y), s X and s Y their standard deviations, and r their correlation coefficient, the general form of its PDF is [52]: We have where µ X and µ Y are the means of X and Y, σ X and σ Y are their standard deviations and ρ their correlation coefficient. Therefore, The univariate Gamma distribution is a continuous probability distribution of a random variable X. The general form of its PDF is: where k is a shape parameter, θ a scale parameter and Γ(u) the gamma function [52]. The mean µ and the standard deviation σ of X, are calculated as following: The bivariate Gamma distribution is a continuous probability distribution of two variables X and Y. The copula version of its PDF is: and k X , k Y are shape parameters, ϑ X , ϑ Y are scale parameters, and ρ the correlation coefficient of X and Y.
If µ X and µ Y are the means of X and Y, σ X and σ Y their standard deviations, and ρ their correlation coefficient, it can be shown that:

. Copulas
If µ X and µ Y are the means of the variables X and Y, σ X and σ Y their standard deviations, and ρ their correlation coefficient, it can be shown that the bivariate PDFs of the other copulas of the program are defined as follows: Appendix A.4.1. X: Normally Distributed-Y: Lognormally Distributed