A new multivariate data analysis model: constrained Naïve Bayes

Trabajo Fin de M´ aster A new multivariate data analysis model: constrained Na¨ıve Bayes Presented by: Mar´ıa Remedios Sillero Denamiel Supervisors: Dr. Rafael Blanquero Bravo, Universidad de Sevilla Dr. Emilio Carrizosa Priego, Universidad de Sevilla External Supervisor: Dr. Pepa Ram´ ırez Cobo, Universidad de C´adiz June 24, 2016 Agradecimientos En primer lugar quiero agradecer a los doctores Emilio Carrizosa Priego y Rafael Blanquero Bravo quienes, hace poco m´as de un a˜no, me brindaron la oportunidad de conocer qu´e es la investigaci´on matem´atica, un trabajo que requiere mucha dedicaci´on y constancia, las cuales me han sabido transmitir a la perfecci´on. El trabajo que aqu´ı se presenta me ha permitido conocer a la doctora Pepa Ram´ırez Cobo, quien adem´as de volcarse con ´este, ha sido un pilar fundamental para m´ı en este tiempo. Como alguien dir´ıa: “Las personas con las que trabajo son, sin duda, lo mejor de este trabajo”. A mis padres, Andr´es y Aurelia, y a mi hermana Roc´ıo, ellos son mi apoyo incondicional y los que me han ense˜nado a luchar por lo que quiero. Gracias Antonio J. Molero por animarme siempre y por estar siempre ah´ı. Gracias a todos por confiar en m´ı. 3 Contents Introduction 7 1 Motivation and Background 9 1.1 Definition of a classifier and performance measures . . . . . . . . . . . 9 1.2 Overall view of classification methods . . . . . . . . . . . . . . . . . . . 10 1.3 Benchmarkclassifiers............................ 11 1.3.1 Support Vector Machine . . . . . . . . . . . . . . . . . . . . . . 11 1.3.2 RandomForest ........................... 13 1.3.3 Na¨ıveBayes............................. 14 2 Na¨ıve Bayes classifier 15 2.1 Description ................................. 15 2.2 Discussion about the independence hypothesis . . . . . . . . . . . . . . 16 2.2.1 Maximal Information Coefficient (MIC) . . . . . . . . . . . . . . 17 2.3 Parameters’ estimation methodologies . . . . . . . . . . . . . . . . . . . 19 2.3.1 The frequentist versus the Bayesian paradigms . . . . . . . . . . 19 2.3.2 The Normal-Gamma model . . . . . . . . . . . . . . . . . . . . 21 2.3.3 A simulation example . . . . . . . . . . . . . . . . . . . . . . . . 25 2.4 Arealdatasetexample........................... 29 2.4.1 Univariate and bivariate analysis . . . . . . . . . . . . . . . . . 29 2.4.2 Comparison between the frequentist and Bayesian approaches . 30 3 The constrained Na¨ıve Bayes 43 3.1 Motivation and formulation of the optimization problem . . . . . . . . 43 3.2 An algorithm to solve the optimization problem . . . . . . . . . . . . . 45 3.3 Numericalexample ............................. 46 Conclusions 51 5 Introduction The classification problems are in vogue due to the fact that they are necessary in many filds of real life, where the datasets gather so much information. Some examples that reveal the importance of these problems are: the early detection of diseases, the granting of credit to a certain individual,. . . It seems clear that to work with good features, that is, attributes that are able to distinguish the different classes, is a key point. Over the years, the classification problems have been studied and many classifiers have been developed. Some examples of them are: Support Vector Machine ([5]), Random Forest ([3]), Na¨ıve Bayes ([6],[11],[8],[10]), between others. Once a classifier is defined, the question that everyone wants to know is how good is a given classification function. In order to solve this question, measures of effectiveness were defined ([16]). This work focuses on the Na¨ıve Bayes classifier. A thorough study of the behavior of the Na¨ıve Bayes classifier, the effect of assuming independence and the influence of the involved parameters estimation methods is presented. Finally, as a major contribution of this work, a different version of Na¨ıve Bayes classifier is presented, in which the estimation is made by imposing constraints on the effectiveness measures on the obtained classifier. The structure of the work is described next. In Chapter 1 an introduction to classification is made. In Chapter 2, Na¨ıve Bayes classification is explained in more detail, in addition to introduce a coefficient for the calculation of both linear and nonlinear relationships between variables ([15], [18]). Finally, the different parameters estimation methodologies will be compared. In Chapter 3 a study of the novel approach will be carried out. 7 Chapter 1 Motivation and Background 1.1 Definition of a classifier and performance measures A classifier is a function, called f, that associates input variable vectors x= (x1, . . . , xp)∈ Xto output classes y∈ {C1, . . . , CK}, where Xis the space of variables and xi∈R although they will be continuous and discrete. Moreover, it is considered that classes are categorical and respectively exclusive. The aim of statistical classification is to learn ffrom a labeled training dataset of Ninput-output pairs, (xn, yn), n= 1 . . . N, where xnare the individuals and yndenotes the class for such individual. Once the classifier has been obtained using the training set, its performance is measured ([16]) by using a test dataset so as to avoid the over-fitting. Note that both samples, training and test, can be assumed to be generated from the same population. Let the true class for a given individual be denoted xand the class for the same individual indicated by the classifier be denoted y. Although it is expected that they are the same for as many individuals as possible in the test dataset, nevertheless, the classifier will make errors. The probabilities that the classifier will return correct results on each of the classes are estimated in a frequentist manner by dividing the number of rightly classified individuals for a given class by the total number of individuals in that class in the test dataset. The probability P(y=Ck|x=Ck) is understood as the probability that the classifier designates class Ck, given that the individual belongs to class Ck. Now, P(y=Ck|x=Cl), k, l ∈ {1, . . . , K},k6=l, shows the probability that the classifier fails in classifying an individual actually from class Clas an individual in class Ck. 9 independent on the class), in practice, the interest is in computing the numerator, which can be written as the joint distribution p(x, Ck) = p(x1|x2. . . xn, Ck)p(x2|x3. . . xn, Ck). . . p(xn|Ck)π(Ck).(2.3) The key assumption of the NB, which makes it a so tractable classifier even for large values of p, is the independence of the features conditioned to the class, which implies that (2.3) can be simplified to p(x, Ck) = π(Ck) p Y i=1 p(xi|Ck). Finally, the probabilities of interest (2.2) will be computed according to p(Ck|x)∝π(Ck) p Y i=1 p(xi|Ck). Note that in order to implement the NB classifier, a prior distribution π(·) for the class as well as a probability distribution for the features conditioned to the class Xi|Ck need to be selected. Common choices for the last include the Normal and Multinomial distributions, for the continuous and discrete cases, respectively. Concerning the choice of the prior, it can be based on the researcher’s previous knowledge of the problem or, in the case of lack of prior knowledge, it can be set equal for all classes (that is, π(Ck)=1/K, for all k= 1, . . . , K). 2.2 Discussion about the independence hypothesis As previously commented, the independence assumption in the NB notably simplifies the computation of the conditional probabilities as in (2.2). However, this assumption is generally not true in practice and violated for many real datasets. The aim of this section is to discuss how this affects the performance of the classifier. As a toy example, consider three features X1, X2, X3such that X1and X2are independent, but X3is highly correlated with X1. The assumption of independence leads to decomposing the joint density as the product of the three marginals, and therefore the influence of X1(≃X3) is different than it should be. Hence, it is natural to think that such miscalculations may possibly deteriorate the classifier’s performance. However, the NB has proven to be comparable or superior to many well-known classification alternatives. In [6] this fact is explained as follows: “...although the individual class density estimates may be biased, this bias might not hurt the posterior probabilities as much, especially near the decision regions. In fact, the problem may be able to withstand considerable bias for the savings in variance such a naive assumption earns...” 16 Also, [8] provides another possible explanation about the over performance of the NB, which is based on the features’ selection process. Due to the measurement of variables is often quite expensive (for example, in biomedical contexts as metabolomics), a previous selection step is undertaken in order to eliminate redundant information. The result of such selection is a set of variables which tend to be weakly correlated. Next, we explore more in depth the effect of dependent features on the NB performance via a numerical example. Assume that the number of classes is only two (positive and negative), and four features (M1, M2, M3 and M4) are simulated according to a multivariate Normal distribution with means vectors given by µ+= (2,3,4,1), µ−= (3,6.5,2,3.5) and covariance matrix defined by M1M2M3M4 M1 M2 M3 M4     1.0 0.9 0 0 0.9 1.0 0 0 0 0 1.0 0 0 0 0 1.0     Note that features are linearly independent, except for M1 and M2 which are highly correlated. In order to train the NB classifier (where Normal conditional marginals are assumed), a sample of size equal to 1000 is used. For validation of the method, a sample of 300 observations is considered. In both cases, the proportion of the classes is equal. To implement the NB, the routine NaiveBayes from the library klaR is used (see the Appendix for the details). Table 2.1 shows the performance values (correct classification rates) for both classes, for all possible combinations of features. From the table, it can be observed how the results under the selection (M2, M3, M4) (that is, when one dependent variable is removed) are better than using the complete information provided by the four variables. This fact points towards the importance of a proper feature selection method that discards noisy information. Because of the independence assumption in the NB, such attributes’ selection process should produce a final set of independent features with high significance. Traditionally, the Pearson correlation coefficient has been used to measure the (linear) dependence between random variables. However, other type of correlations (nonlinear) may be present in the data and are ignored by the Pearson coefficient. The next section introduces an index, on which such feature selection process might be based, that takes into account nonlinear associations. 2.2.1 Maximal Information Coefficient (MIC) The Maximal Information Coefficient (MIC) (see [15],[18]) measures nonlinear relationships among random variables. Its definition arises from that of the concept of mutual 17 Combination CCR Positive Class CCR Negative Class M1 0.741 0.708 M2 0.954 0.96 M3 0.838 0.829 M4 0.888 0.899 M1, M2 0.935 0.952 M1, M3 0.865 0.896 M1, M4 0.905 0.895 M2, M3 0.968 0.993 M2, M4 0.973 0.967 M3, M4 0.935 0.959 M1, M2, M3 0.942 0.979 M1, M2, M4 0.966 0.948 M1, M3, M4 0.953 0.947 M2, M3, M4 0.987 0.993 M1, M2, M3, M4 0.974 0.98 Table 2.1: Results of each combination of variables. information (MI) (see, [9]), which quantifies the information about one variable Xthat is provided by a different variable Y. The MI coefficient performs by exploring if the products of marginal distribution p(X)p(Y) is similar to the joint distribution p(X, Y ). Specifically, the MI coefficient between two discrete variables Xand Yis defined as MI(X, Y ) = X y∈YX x∈X p(x, y) log p(x, y) p(x)p(y), where p(x) and p(y) denote the probability mass functions of Xand Y, and p(x, y) is the joint probability mass function of Xand Y. Similarly, for the continuous case, the coefficient is given by MI(X, Y ) = ZYZX p(x, y) log p(x, y) p(x)p(y)dxdy, where now p(·) refers to a density function. The MI coefficient presents two main problems. The first is from a computational viewpoint, since its calculation depends on a two dimensional smoothing. The second problem is that it is not bounded, which makes its interpretation difficult. After more than 50 years of development of the MI index, the MIC coefficient has been proposed to overcome the drawbacks of the MI. The MIC captures the relationships between two variables using a grid on the scatterplot of these two variables. To calculate the 18 MIC, all possible grids are explored up to a maximal grid resolution (depending on the sample size), and for each pair of integer values (x, y), the highest possible MI realizable by any x-by-ygrid is computed. Then, such MI values are normalized in the interval [0,1] in order to guarantee the comparison between all grids. To calculate the MIC coefficient, first the matrix M= (mx,y) is defined where mx,y is the largest normalized mutual information obtained by any x-by-ygrid, so that the maximum value in Mwill be the value of the MIC coefficient. Formally, let Ga grid and let IG be the mutual information of the probability distribution induced on the boxes of G, with the probability of a box proportional to the number of points dropping inside the box. The (x, y)-th entry of Mis denoted by mx,y, defined as mx,y = max IG/log min{x, y}, where this maximum is calculated over all x-by-ygrids G. The MIC index is defined as the maximum of mx,y between all ordered pairs (x, y) satisfying xy < B, where Bis dependent on the sample size. Note that the elements of Mtake values between 0 and 1, and consequently, the MIC does, too. Also, the MIC index satisfies the symmetry property, that is, MIC(X, Y ) = MIC(Y, X), because of the symmetry of the MI coefficient and because of IGdepends on the rank order of the data. Finally, it should be noted that to compute M, it is necessary to optimize over all grids, so it seems clear that computing tools will be need to calculate the MIC value. In particular, the software R computes the MIC coefficient through the command mine, included in the library miverva (see the Appendix for the details). 2.3 Parameters’ estimation methodologies As described in Section 2.1., in order to implement the NB classifier, a probability distribution for the features conditioned to the class needs to be selected, Xi|Ck∼Fθi,k (x). Once such model is set, then it needs to be estimated through sample data. Throughout this work, it will be assumed that Xi|Ckfollows a Normal distribution with unknown mean and variance denoted by θi,k =µi,k, σ2 i,k. In next section we explore two possible alternatives for estimating the model’s parameters. 2.3.1 The frequentist versus the Bayesian paradigms Two main approaches can be considered for statistical inference: the frequentist (or classic) and the Bayesian one. The first approach undertakes estimation just by taking into account the observed sample data and, basing on the concept of frequency, it 19 derives point estimates for the model parameters by maximizing the likelihood function. Therefore, under a Gaussian likelihood, the frequentist approach gives as parameters’ points estimates the sample mean and the sample variance. The second approach considers instead the possible prior knowledge that the researcher possesses about the problem, and this previous knowledge is updated (via the Bayes formula) by the likelihood function once data are observed. The use of the prior knowledge implies assigning a probability distribution (the so-called prior distribution) to the model parameters, which are treated as random variables (instead of unique, unknown values as in the frequentist approach). Such prior distribution will transform (using the likelihood function) into the posterior distribution, which is the object of interest for Bayesian statisticians. Formally, the Bayesian inference approach performs as follows. Suppose that before the experiment is undertaken, our prior distribution describing the model parameter θis π(θ).The data are coming from an assumed model (likelihood) which depends on the parameter and is denoted by f(x|θ).Bayes theorem updates the prior π(θ) to the posterior by accounting for the data x, π(θ|x) = h(x, θ) m(x)=f(x|θ)π(θ) m(x), where m(x) is a normalizing constant, m(x) = RΘf(x|θ)π(θ)dθ. Some of the advantages of the Bayesian statistics are pointed out next: •The uncertainty is expressed via the probability distribution. The statistical inference follows a conceptually simple recipe embodied in Bayes’ theorem. •Available prior information is coherently incorporated intro the statistical model describing the data. •The FDA (US Food and Drug administration) recommends the use of a Bayesian methodology in the design and analysis of clinical trials for medical devices. Some of the reasons are: –Valuable prior information is often available. –The use of prior information may alleviate the need for a larger sized trial. –Bayesian methods allow for great flexibility in dealing with missing data. –Bayesian models facilitate meta-analysis. For a detailed description of the Bayesian paradigm, we refer the reader to [19]. In this work, because of the versatility that priori distributions provide, we will assume a Bayesian framework. The next section introduces the Normal-Gamma model, as a useful tool for Bayesian inference of the Normal distribution. 20 2.3.2 The Normal-Gamma model Consider a pair of random variables (X, Y ) such that the conditional distribution X|Y is normal X|Y∼Nµ, 1 kY . Also, Yis distributed according to a Gamma probability model, Y|(α, β)∼Ga(α, β). Then, it is said that the pair (X, Y ) follows a Normal-Gamma distribution, noted by (X, Y )∼NG(µ, k, α, β). The joint probability density function can be found as f(x, y |µ, k, α, β) = βα√k Γ(α)√2πyα−1 2e−βy−ky(x−µ)2 2.(2.4) Figure 2.1 depicts the previous density function for an assortment of parameters’ values. Some properties of the model are as follows. Concerning the expected values of the random variables, they are given by E(X) = µ, E(Y) = α β. Regarding the variance values, they can be found as V(X) = β k(α−1), V (Y) = α β2. Finally, it can be proven that (2.4) is unimodal with the maximum attained at µ, α−0.5 β. The usefulness of the Normal-Gamma model in Bayesian statistics is that it is a conjugate model. This implies that, for a Normal likelihood (like the case considered in this work), a Normal-Gamma prior transforms into a Normal-Gamma posterior distribution, and therefore, the posterior does not need to be numerically calculated. This fact is stated by the next Theorem (see [12]). Theorem 2.3.1 Let Xfollow a Normal distribution, X∼N(µ, σ2)and let x= (x1, . . . , xn)denote a simple random sample of X. Assume that µand σ2follow a Normal-Gamma prior distribution, that is µ, 1 σ2∼NG(µ0, k0, α0, β0). 21 (a) NG(µ= 0.1, k = 2, α = 1, β = 1) (b) NG(µ= 0.1, k = 2, α = 3, β = 1) (c) NG(µ= 0.1, k = 2, α = 5, β = 1) (d) NG(µ= 0.1, k = 2, α = 5, β = 3) Figure 2.1: The probability density function of different Normal-Gamma models 22 Then, the posterior distribution of Xis also a Normal-Gamma distribution, given by µ, 1 σ2|x∼NG(µn, kn, αn, βn), where µn=k0µ0+n¯ x k0+n kn=k0+n αn=α0+n 2 βn=β0+1 2 n X j=1 (xj−¯ x)2+k0n(¯ x−µ0)2 2(k0+n). Proof: The likelihood function is l(x|µ, σ2) = n Y i=1 p(xj|µ, σ2) = 1 √2πσ2n exp (−1 2σ2 n X j=1 (xj−µ)2).(2.5) Let λdenote the so-called precision, λ=1 σ2. Then l(x|µ, σ2) = 1 (2π)n/2λn/2exp (−λ 2 n X j=1 (xj−µ)2).(2.6) Let ¯ xand s2denote the empirical mean and variance: ¯ x=1 n n X j=1 xj s2=1 n n X j=1 (xj−¯ x)2. Since n X j=1 (xj−¯ x)(µ−¯ x) = (µ−¯ x) n X j=1 xj!−n¯ x!= (µ−¯ x)(n¯ x−n¯ x)=0. 23 then, it is possible to rewrite the exponential term as follow n X j=1 (xj−µ)2= n X j=1 [(xj−¯ x)−(µ−¯ x)]2= = n X j=1 (xj−¯ x)2+ n X j=1 (¯ x−µ)2−2 n X j=1 (xj−¯ x)(µ−¯ x) = =ns2+n(¯ x−µ)2(2.7) Assume that (µ, λ) follows a Normal-Gamma conjugate prior distribution, then f(µ, λ|µ0, k0, α0, β0) = 1 Γ(α0) βα0 02π k01 2 λα0−1 2exp −λ 2[k0(µ−µ0)2+ 2β0], where µ0, k0, α0, β0are assumed to be known. Hence, the posterior density is f(µ, λ|x)∝f(µ, λ|µ0, k0, α0, β0)p(x|µ, λ) ∝λ1 2e−(k0λ)(µ−µ0)2 2λα0−1e−β0λλn 2e−(λ 2)Pn j=1(xj−µ)2 From (2.7) n X j=1 (xj−µ)2=n(µ−¯ x)2+ n X j=1 (xj−¯ x)2. Finally, it can be shown that k0(µ−µ0)2+n(µ−¯ x)2= (k0+n)(µ−µn)2+k0n(¯ x−µ0)2 k0+n, where µn=k0µ0+n¯ x k0+n. Therefore, k0(µ−µ0)2+ n X j=1 (xj−µ)2=k0(µ−µ0)2+ n X j=1 (xj−¯ x)2= = (k0+n)(µ−µn)2+k0n(¯ x−µ0)2 k0+n+ n X j=1 (xj−¯ x)2. Thus, 24 f(µ, λ|x)∝λ1 2e−λ 2(k0+n)(µ−µn)2λα0+n 2−1e−β0λe−(λ 2)Pn j=1(xj−¯ x)2e−(λ 2)k0n(¯ x−µ0)2 k0+n. In conclusion, the posterior density follows a Normal-Gamma distribution p(µ, λ|x) = NG(µ, λ|µn, kn, αn, βn) with µn=k0µ0+n¯ x k0+n kn=k0+n αn=α0+n 2 βn=β0+1 2 n X j=1 (xj−¯ x)2+k0n(¯ x−µ0)2 2(k0+n).  2.3.3 A simulation example Here, we estimate from a Bayesian viewpoint and using the Normal-Gamma model described in the previous section, the probability density in a Gaussian NB classification context. The model parameters are estimated by maximizing the posterior density (Maximum a posteriori estimates). In this case, the parameters are Normal-Gamma distributed a priori, and consequently, according to Theorem 2.3.1, they follow the same probability model (with different parameters) a posteriori. The maximum a posteriori estimate The maximum a posteriori (MAP) estimates are parameters’ point estimates obtained as the mode of the posterior distribution. Therefore, once the posterior distribution is known, then it needs to me maximized. Here, a simplification of the posterior density is found, under the Normal-Gamma setting. Assume first that feature Xiconditioned to the class Ckis Gaussian distributed Xi|Ck∼N(µi,k, σ2 i,k), 25 Sepal.Length 2.0 2.5 3.0 3.5 4.0 0.5 1.0 1.5 2.0 2.5 4.5 5.5 6.5 7.5 2.0 2.5 3.0 3.5 4.0 Sepal.Width Petal.Length 1234567 4.5 5.5 6.5 7.5 0.5 1.0 1.5 2.0 2.5 1234567 Petal.Width Iris Data Figure 2.3: Scatterplots for each pair of variables in the Iris dataset 32 SL SW PL PW SL SW PL PW     0.99 0.66 0.23 0.19 0.66 1 0.3 0.22 0.23 0.3 1 0.14 0.19 0.22 0.14 0.9     SL SW PL PW SL SW PL PW     0.99 0.40 0.62 0.29 0.40 0.99 0.22 0.3 0.62 0.22 1 0.28 0.29 0.3 0.28 1     SL SW PL PW SL SW PL PW     0.99 0.36 0.60 0.50 0.36 1 0.40 0.38 0.6 0.39 1 0.52 0.50 0.38 0.52 0.99     Note from the MIC values how the hypothesis of independence becomes unrealistic in this example (especially for some pairs of variables). Some details are given next, related to how the numerical experiment will be conducted. First, all possible combinations of variables will be considered. For each combination, we will proceed as follows. •The full dataset is separated in training dataset and validation samples with samples sizes given by the Table 2.7. Setosa Virginica Versicolor Training 35 35 35 Validation 15 15 15 Table 2.7: Distribution of samples in training and validation. •In order to get reliable results and to ensure that the results do not depend on a certain choice of the training set and validation set, the previous step will be repeated 10 times in a random way, so that the final performance measure will be the average/median of the measures for each fold on the validation set. Then, 33 –in the case of Classic (frequentist) NB classifier, as the parameters estimates are the sample mean and sample variance, we get them from the training sample. Once these estimations are set, then the correct classification rate for each class will be calculated on the validation dataset. The R command NaiveBayes will be used for this case. –in the case of the Bayesian NB classifier, where the parameters estimates are found as the solution to the maximization of the objective function (2.9) using the training sample, then the correct classification rate for each class will be calculated on the validation dataset. Here, similarly as in Section 2.3., the R command nmkb will be used to maximize the posterior density. Results We focus now on the classification rates under both estimation methods, shown by Tables 2.8-2.11. In view of the results obtained, it seems that, for the Iris dataset, in most of variable sets the average behavior is more or less similar under both approaches, although the frequentist estimation of parameters of the NB classifier gets more balanced values of correct classification rate for the three classes, while the Bayesian Inference of the parameters returns some cases that are more unbalanced in this aspect. As another observation, the best accuracy has been obtained, in both cases, for the combination Petal Width and Petal Length. It is important to note that these are independent in the three classes, because of the value 0.52 is the highest value that they take. This clearly points out to the already commented importance of the attributes selection process. Figures 2.4 and 2.5 depict the boxplots of the estimations for the precision parameter under the frequentist and Bayesian approaches, respectively, for Petal Length and Petal Width variables. Note from the first figure how the classic method leads to more variable results. Note also, how the Bayesian one seems more robust and its results are not very affected by the choice of the starting point in the maximization algorithm. Figures 2.6-2.7 are the analogous to the previous ones for the means. In this case, both estimation methods look more similar. We found it of interest to visually show how the prior density transformed into the posterior. Figures 2.8-2.9 depict such results for the case of the Petal Length variable conditioned to the Setosa class. Even though the densities belong to the same family of distributions it is clear the change in the modes. 34 Figure 2.4 Figure 2.5 35 Figure 2.6 Figure 2.7 36 Figure 2.8: Prior density for the Gaussian parameters of the Petal Length variable conditioned to the class Setosa. 37 Figure 2.9: Posterior density for the Gaussian parameters of the Petal Length variable conditioned to the class Setosa. 38 Sepal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 86.66 62 62.67 70.44 Median 86.67 60 66.67 71.11 Variation coefficient 0.14 0.18 0.18 0.08 Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 71.33 42.67 56 56.67 Median 73.33 43.33 56.67 57.78 Variation coefficient 0.13 0.17 0.22 0.1 Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 92 93.33 95.11 Median 100 93.33 93.33 94.44 Variation coefficient 0 0.75 0.07 0.02 Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 98 91.33 98 95.78 Median 100 90 100 95.56 Variation coefficient 0.03 0.06 0.03 0.03 Sepal Length & Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 98.67 63.33 74.67 78.89 Median 100 60 76.67 80 Variation Coefficient 0.03 0.14 0.14 0.04 Sepal Length & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 80 87.33 89.11 Median 100 86.67 86.67 87.78 Variation coefficient 0 0.16 0.12 0.04 Sepal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 99.33 91.33 95.33 95.33 Median 100 90 93.33 95.56 Variation coefficient 0.02 0.06 0.03 0.02 Sepal Width & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 82 89.33 90.44 Median 100 83.33 93.33 91.11 Variation coefficient 0 0.13 0.1 0.04 Sepal Width & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 99.33 89.33 94 94.22 Median 100 86.67 96.67 94.44 Variation coefficient 0.02 0.08 0.07 0.03 Table 2.8: Classic Na¨ıve Bayes. 39 Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 93.33 96.67 96.67 Median 100 93.33 100 96.67 Variation coefficient 0 0.07 0.05 0.03 Sepal Length & Sepal Width & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 72.67 87.33 86.67 Median 100 80 86.67 86.67 Variation coefficient 0 0.15 0.12 0.05 Sepal Length & Sepal Width & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 91.33 88.67 93.33 Median 100 90 86.67 94.44 Variation coefficient 0 0.06 0.09 0.03 Sepal Length & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 93.33 93.33 95.56 Median 100 93.33 93.33 95.56 Variation coefficient 0 0.07 0.06 0.02 Sepal Width & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 92.67 96.67 96.44 Median 100 90 100 96.67 Variation coefficient 0 0.07 0.05 0.03 Sepal Length & Sepal Width & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 92.67 94.67 95.78 Median 100 90 93.33 95.56 Variation coefficient 0 0.07 0.04 0.02 Table 2.9: Classic Na¨ıve Bayes 40 Sepal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 44 25.33 98 55.78 Median 46.67 26.67 100 56.67 Variation coefficient 0.3 0.35 0.05 0.09 Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 98 0 40 46 Median 100 0 36.67 44.44 Variation coefficient 0.03 0 0.02 0.06 Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 83.33 98.67 94 Median 100 86.67 100 94.44 Variation coefficient 0 0.14 0.03 0.04 Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98 69.33 89.11 Median 100 100 73.33 88.89 Variation coefficient 0 0.03 0.12 0.03 Sepal Length & Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 97.33 42 91.33 76.89 Median 100 40 93.33 75.56 Variation coefficient 0.04 0.2 0.08 0.05 Sepal Length & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 94 86.67 93.56 Median 100 93.33 86.67 93.33 Variation coefficient 0 0.07 0.11 0.03 Sepal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 60 100 86.67 Median 100 60 100 86.67 Variation coefficient 0 0.1 0 0.02 Sepal Width & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 99.33 62.67 87.33 Median 100 100 63.33 87.78 Variation coefficient 0 0.02 0.3 0.07 Sepal Width & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 90 98 96 Median 100 86.67 100 95.56 Variation coefficient 0 0.07 0.03 0.03 Table 2.10: Na¨ıve Bayes and Bayesian Inference. 41 Sepal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean - - - - Median - - - - Variation coefficient - - - - Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean - - - - Median - - - - Variation coefficient - - - - Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 96 88 94.67 Median 100 100 86.67 93.33 Variation coefficient 0 0.07 0.08 0.02 Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98 69.33 89.11 Median 100 100 73.33 88.89 Variation coefficient 0 0.03 0.12 0.03 Sepal Length & Sepal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean - - - - Median - - - - Variation coefficient - - - - Sepal Length & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98 76.67 91.56 Median 100 100 80 93.33 Variation coefficient 0 0.03 0.2 0.05 Sepal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 87.33 96.67 60.67 81.56 Median 90 100 66.67 82.22 Variation coefficient 0.16 0.05 0.18 0.07 Sepal Width & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98.67 72.67 90.44 Median 100 100 73.33 91.11 Variation coefficient 0 0.03 0.25 0.06 Sepal Width & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 92 98 63.33 84.44 Median 100 100 66.67 85.56 Variation coefficient 0.19 0.03 0.3 0.09 Table 3.1: Constrained Na¨ıve Bayes and Bayesian Inference. 48 Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98 88 95.33 Median 100 100 90 95.56 Variation coefficient 0 0.05 0.13 0.03 Sepal Length & Sepal Width & Petal Length CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 96 72 89.33 Median 100 96.67 73.33 88.89 Variation coefficient 0 0.05 0.22 0.06 Sepal Length & Sepal Width & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 90.67 92.67 54.67 79.33 Median 96.67 93.33 56.67 78.89 Variation coefficient 0.15 0.11 0.26 0.05 Sepal Length & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 98 58.67 85.56 Median 100 100 60 85.56 Variation coefficient 0 0.03 0.29 0.07 Sepal Width & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 99.33 98 73.33 90.22 Median 100 100 70 90 Variation coefficient 0.02 0.05 0.2 0.05 Sepal Length & Sepal Width & Petal Length & Petal Width CCR Setosa CCR Virginica CCR Versicolor Accuracy Mean 100 96.67 65.33 87.33 Median 100 100 63.33 86.67 Variation coefficient 0 0.05 0.28 0.06 Table 3.2: Constrained Na¨ıve Bayes and Bayesian Inference. 49 50 Conclusions Three main conclusions emerge from the study on the NB classifier that has been carried out in this work. First, NB classifier assumes independence of variables. When it is applied with dependent variables, the full set of variables may not yield the best results in terms of accuracy. This is illustrated by the example of the Iris dataset. In fact, Feature Selection, that has not been considered in this work beyond the complete enumeration of set of variables in the same example, is concluded to be necessary. Another main conclusion is that the different parameters estimation methods induce, of course, different classifiers. This work does not have as an objective to indicate which one is right (if the frequentist approach, if the Bayesian approach), but to show with an example that different methods induce different classifiers. Although, regarding the Bayesian approach, it remains to be done an exhaustive study of the sensitivity priors, because different a priori estimates will lead to different estimates values and therefore to different classifiers. The novelty of this work is the inclusion of constraints on the classical Bayesian approach. These constraints have provided a different classifier, with the difference that this new approach allows to control the performance measures, which is fundamental in many realworld classification problems. 51 52 Appendix: R routines Na¨ıve Bayes routines The language and environment R (see [14]) provides the following two packages relating to Na¨ıve Bayes classifier: •The package e1071, that provides the functions –naiveBayes: to perform Na¨ıve Bayes classification. –predict.naiveBayes: to obtain the predictions. •The package klaR that has the following function to perform Na¨ıve Bayes classification and to predict –NaiveBayes –predict.NaiveBayes In next subsections, let us explain in more detail the previous functions. Package “e1071” Function naiveBayes can be called in two different ways, that is, with different type of arguments. The first one would be: naiveBayes(x, y, laplace=0, subset, na.action=na.pass) •x: Numerical matrix or data-frame of numeric variables and/or categorical. •y: Vector of classes. •laplace: Laplace smoothing parameter to improve the estimates of the probabilities in the case of categorical data (by default, it does not apply). 53 •subset: Vector of indices specifying instances of the sample learning to use. •na.action: Action to perform in the presence of missing values (NA). And the second one is: naiveBayes(formula, data, laplace=0, subset, na.action=na.pass) •formula: Formula in the way class ∼x1+x2+. . . •data: Learning sample data-frame or contingency table with the frequencies resulting from the count. •laplace,subset,na.action: As before. On the other hand, the function naiveBayes returns an object with the next components: •apriori: Probabities of each class of the response variable. •tables: List of tables, one for each predictor variable. According the type of it, the information obtained is different: –Qualitative: the table contains the probabilities of each class of the output variable conditional on the different modalities of the predictor variables. –Quantitative: for each class the output variable, the table contains the mean and standard deviation of the predictor variable conditioned to that class. •levels: Output variable classes. •call: Sequence of the call to the naiveBayes function. Finally, once the classifier has been constructed according to the previous functions, predictions can be made: predict(object, newdata, type=class,threshold=0.001, eps=0) where •object: Object generated by the naiveBayes function. •newdata: Data-frame with the instances to predict. 54 •type: If raw is selected, the function returns the probability of each class of output variable (for each instance to predict). If class is selected, it only indicates the class with maximum probability for each instance to predict. •threshold: Value for which the probabilities are replaced below the value specified in eps. •eps: The lower probabilities than this value will be replaced by threshold. Package “klaR” There exists another function called NaiveBayes from klaR. The main differences between the naiveBayes function (e1071) and the NaiveBayes function (klaR) are the following: •NaiveBayes allows one to specify prior probabilities of classes of the response variable. In naiveBayes they are estimated from the learning sample. •NaiveBayes lets one use kernel function to estimate the density function of continuous variables. In naiveBayes always it is assumed that this kind of variables follows a Normal distribution. •NaiveBayes requires learning sample given as matrix or data-frame (it does not support contingency table). MIC routine R provides a function that performs the calculation of the MIC. The library miverva has a function called mine, which returns the calculation of the MIC for a data set given. The principal argument needed is x, a numeric vector, matrix or data frame (which is coerced to matrix). The rest of arguments are optional and they are explained in the manual provided by R. This function returns, among other outputs, the Maximal Information Coefficient (MIC). Optimization Routines In this section two different functions from two distinct packages will be explained. First, it will be described a function that is useful when unconstrained optimization problems are carried out. The second one allows one to include constraints to the 55 optimization problem and, moreover, both constraints as the objective function can be given as blackbox functions. Package “dfoptim” The library dfoptim has a function called nmkb. This function is an implementation of the Nelder-Mead algorithm for derivative-free optimization. It allows bounds over the parameters. The arguments of this function are: nmkb(par, fn, lower=-Inf, upper=Inf, control = list(), . . . ) where •par: A initial vector of parameters. Note that it must be between lower and upper bounds (if they exist). •fn: The nonlinear objective function to be optimized. •lower and upper: lower and upper bounds on the parameters. •control: a list with the control parameters such as tolerance, maximum number of objective function evaluations allowed, . . . •. . . : auxiliary arguments of the objective function. Package “crs” The library crs has a function called snomadr. This last function snomadr is an R interface to NOMAD (Nonsmooth Optimization by Mesh Adaptive Direct Search) ([1]), an open source software C++ implementation of the Mesh Adaptive Direct Search algorithm designed for constrained optimization of blackbox functions. This function needs the following input arguments: •eval.f: the function that returns the value of the objective function and the value of the constraints. •n: the number of variables. •bbin: a vector that indicates the type of the variables. •bbout: a vector that fixed how the constraints have to be considered during the optimization process. 56 •x0: the initial vector of parameters. •lb and ub: vectors with lower and upper bounds of the controls, respectively. •. . . : auxiliary arguments that will be needed to the objective and constraints functions. Note that this function has some additional arguments which are explained in the help of R. 57