PLS Generalized Linear Regression and Kernel Multilogit Algorithm (KMA) for Microarray Data Classification

TL;DR

This paper introduces extensions of PLSGLR combined with logistic regression and LDA for microarray data classification, compares them with classical methods, and demonstrates that KMA achieves the lowest error rates.

Contribution

It develops new PLSGLR-based classifiers and compares their performance with existing methods, highlighting KMA's superior accuracy on microarray data.

Findings

KMA has the lowest classification error rates among tested methods.

Extensions of PLSGLR improve classification performance.

KMA outperforms classical classifiers on both preprocessed and raw data.

Abstract

We implement extensions of the partial least squares generalized linear regression (PLSGLR) due to Bastien et al. (2005) through its combination with logistic regression and linear discriminant analysis, to get a partial least squares generalized linear regression-logistic regression model (PLSGLR-log), and a partial least squares generalized linear regression-linear discriminant analysis model (PLSGLRDA). These two classification methods are then compared with classical methodologies like the k-nearest neighbours (KNN), linear discriminant analysis (LDA), partial least squares discriminant analysis (PLSDA), ridge partial least squares (RPLS), and support vector machines(SVM). Furthermore, we implement the kernel multilogit algorithm (KMA) by Dalmau et al. (2015)and compare its performance with that of the other classifiers. The results indicate that for both un-preprocessed and preprocessed data, the KMA has the lowest classification error rates.

Tables

Table 1. Table 1: Proposed strategy

Steps

Step 1: Dimension reduction In this stage, we propose to use PLSGLR to project the high dimensional data to a low dimension space thus resulting in new components (latent variables), which preserve the information in the intrinsic structure of the data. Step 2: Use of latent variables for classification Analyze the obtained latent variables with the classical statistical classifiers: i PLSGLR components with logistics regression to get the PLSGLR-logistic model denoted as (PLSGLR-log) ii PLSGLR components with linear discriminant analysis model to get PLSGLR-Linear Discriminant Analysis model denoted as (PLSGLRDA)

Table 2. Table 2: Rate of classification error for the different methods when applied to the un-preprocessed data set

DATA	PLSGLR-log	PLSGLRDA	KNN	LDA	PLSDA	RPLS	SVM	KMA
Colon	38.3	31.7	60.0	25.0	11.7	15.0	18.3	1.7

Table 3. Table 3: classification error rates for the different methods when applied to the preprocessed data set

DATA	PLSGLR-log	PLSGLRDA	KNN	LDA	PLSDA	RPLS	SVM	KMA
Colon	16.4	13.3	26.7	15.0	11.7	11.7	14.8	11.2

Equations

$$\begin{array}{rl}\text{Pr}(\mathbf{y}& =j|\mathbf{x})=\frac{\text{exp}\{f(\mathbf{x};{\mathbf{\theta }}_j)\}}{1+\sum _{i=1}^{C-1}\text{exp}\{f(\mathbf{x};{\mathbf{\theta }}_i)\}},j=\{1,2,\dots ,C-1\}\\ \text{Pr}(\mathbf{y}& =C|\mathbf{x})=\frac{1}{1+\sum _{i=1}^{C-1}\text{exp}\{f(\mathbf{x};{\mathbf{\theta }}_i)\}},\end{array}$$

$${\vartheta }_j=\text{logit}(t_j)=\text{log}\frac{t_j}{t_C},j=\{1,2,\dots ,C-1\}$$

$$\begin{array}{rl}\hat{\mathbf{\Theta }}& =\underset{\Theta }{\mathrm{argmin}}\mathbf{U}(\mathbf{\Theta })\\ U(\mathbf{\Theta })& =\Vert \mathbf{\vartheta }-\mathbf{X\Theta }{\Vert }_F^2+\lambda \Vert \mathbf{\Theta }{\Vert }_F^2,\end{array}$$

$$U(\mathbf{\Gamma })=\Vert \mathbf{\vartheta }-\mathbf{X}{\mathbf{X}}^T\Gamma {\Vert }_F^2+\lambda \Vert {\mathbf{X}}^T\mathbf{\Gamma }{\Vert }_F^2$$

$$\begin{array}{rl}{t}_j^{new}& =\frac{\text{exp}\{{\vartheta }_j^{new}\}}{1+\sum _{i=1}^{C-1}\text{exp}\{{\vartheta }_j^{new}\}},j=\{1,2,\dots ,C-1\}\\ t_J^{new}& =\frac{1}{1+\sum _{i=1}^{C-1}\text{exp}\{{\vartheta }_j^{new}\}}.\end{array}$$

$$\displaystyle g(\theta)=\sum_{h=1}^{m}c_{h}\Big{(}\sum_{j=1}^{p}w_{hj}^{*}\mathbf{x}_{j}\Big{)},$$

$$BS{S}_j/WS{S}_j=\frac{{\sum }_i{\sum }_{k}I(y_i=k)({\bar{x}}_{kj}-{\bar{x}}_{.j}{)}^2}{{\sum }_i{\sum }_{k}I(y_i=k)(x_{ij}-{\bar{x}}_{kj}{)}^2}$$

Full text

1 Introduction

The field of genomics has witnessed a tremendous increase in the amount of data generation due to biotechnological advances like microarrays and next-generation sequencing platforms. These biotechnological advances have made it possible to simultaneously monitor expression levels for thousands of genes, and thus help in solving particular problems related to the identification of molecular variants in them, and their relation to the classification, diagnosis, prognosis and treatment of different conditions. The high dimensional data generated from microarray technology involve many thousands of genes measured simultaneously, a different microarray for each individual. This definitely introduces some noise and unwanted variations that might stem from technical or unknown sources.

In a microarray experiment let $n$ and $p$ be the numbers of the samples and genes respectively, so that the generated data is a $n\times p$ matrix. The main challenge with these technologies is that the resultant data are noisy due to biological and technological variations, and at the same time they usually are high dimensional, i.e., they have more variables than cases due to a low sample size, so $n<<p$. This condition makes the direct application of most classical statistical methodology implausible, leading researchers to propose new solutions for this type of problem.

Normally before the down stream analysis of the data generated from DNA microarrays, a preprocessing and normalization stage is performed to remove the noise, filtering out the genes with low expression values, addressing missing values, and standardizing the data via a log-transformation. One of the most used preprocessing procedures for microarray data was proposed by Dudoit et al. (2002), which entails three basic steps, namely: thresholding, filtering out of genes outside of a range of minimum/maximum intensities, and finally, standardization of the expression values by a log transformation (see(Alshamlan et al. (2013); Dudoit et al. (2002))).

This work considers classification problems for microarray data sets under two conditions: un-preprocessed and preprocessed. In the un-preprocessed data all genes available in the study are included, while in the preprocessed only the subset of genes believed to play important roles in the biological problem of interest are used. We extend the Partial Least Squares Generalized Linear Regression (PLSGLR) algorithm of Bastien et al. (2005) by combining it with Logistic Regression, to give PLSGLR-log, and with Linear Discriminant Analysis to come up with PLSGLRDA. Furthermore, we compare their performance with that of the kernel multilogit algorithm (KMA) proposed by (Dalmau et al., 2015), and of the classical methods: the k-Nearest Neighbour (KNN), Ridge Partial Least Squares (RPLS), Partial Least Squares-Linear Discriminant Analysis (PLSDA), the usual Linear Discriminant Analysis (LDA) and the Support Vector Machines (SVM), when applied to a set of microarray data, referred to in this work as the Colon data set (by (Alon et al., 1999)). We evaluate the classifiers with regard to their classification error rates in this data set and compare them.

Our work addresses problems similar to many studies involving classification in microarrays, with typically high dimensional data and low numbers of samples (or subjects). Following a two stage strategy, many involve the use of the original PLS to build the components, even though the response variables are discrete, for example the analysis of (Nguyen and Rocke, 2002a, b); this is intuitively not correct since the original PLS is an algorithm best suited for continuous response variables. And in almost all of the procedures a variable (gene) selection step is implemented, with an accompanying computing cost. This paper describes a procedure suitable for categorical data, and its performance is studied with and without the gene selection step, and compared to that of each of the other classifiers used. An additional advantage of our approach is that the PLSGLR can deal with missing values, unlike the original PLS, commonly used in the literature.

The proposed two stage strategy for the classification problem is described as follows.

To the best of our knowledge, the proposed combination of PLS generalized linear regression algorithm with logistic and discriminant analysis has not been used before in cases where $n<<p$. The PLS generalized linear regression algorithm is simple, and a good performance when compred to the classical methods would make it an attractive alternative.

2 Kernel multilogit algorithm (KMA)

The KMA was recently proposed by Dalmau et al. (2015). This algorithm works by first transforming a categorical response variable to a continuous one via a multilogit transformation. The transformed variable is then used with the explanatory variables in a regression model for classification and prediction. Finally, the new predicted variables are transformed back using the inverse multilogit function to the original space to enable classification.

Let the response variable vector $\bm{y}$ be categorical with class labels $\{1,2,\ldots,C\}$. To classify a discrete variable from predictor variables $\bm{x}$, the first step is to transform the response variable $\bm{y}$ into a new space using the multilogit function. The multinomial logit model with $C$ as the reference category can be given as

[TABLE]

where $f(\bm{x};\bm{\theta}_{i})=\bm{x}^{T}\bm{\theta}_{i}$. The expected value of $\bm{y}$ being a multinomial random variable is given by $E(\bm{y}|\bm{x})=[\text{Pr}(\bm{y}=1|\bm{x}),\text{Pr}(\bm{y}=2|\bm{x}),\ldots,\text{Pr}(\bm{y}=C|\bm{x})]^{T}$. Now, denoting $\bm{t}=E(\bm{y}|\bm{x})$, the original response variable $\bm{y}$ is not used but instead a transformed version $\bm{\vartheta}=\text{logit}(\bm{t})$ is used. The logit transformation is done with $C$ as the reference category as follows

[TABLE]

where $\vartheta_{j}\in\bm{\vartheta},t_{j}\in\bm{t}$.

In the second step a parametric linear model is proposed and its parameter estimates can be obtained via the standard Bayesian formula $\text{Pr}(\bm{\vartheta}|\bm{x})=\text{Pr}(\bm{x}|\bm{\vartheta})\text{Pr}(\bm{\vartheta})/\text{Pr}(\bm{x})$ where $\text{Pr}(\bm{\vartheta}|\bm{x})$ is the posterior probability distribution, $\text{Pr}(\bm{x}|\bm{\vartheta})$ is the likelihood function and $\text{Pr}(\bm{x})$ is the normalization constant, assuming that $\bm{\vartheta}\in\mathbb{R}^{C-1}$ for a given $\bm{x}\in\mathbb{R}^{m}$ follows a multivariate normal distribution $\bm{\vartheta}|\bm{x}\thicksim\mathcal{N}(\Theta^{T}\mathbf{x},\alpha^{2}\mathbf{I})$, $\Theta\in\mathbb{R}^{m\times C-1}$,$\text{Pr}(\bm{\vartheta}|\bm{x})$ is also multivariate normally distributed. Furthermore, the prior parameters are assumed to follow a normal distribution, i.e. $\bm{\theta}\sim\mathcal{N}(\mathbf{0},\beta^{2}\mathbf{I})$ where $\beta$ is known. The parameter matrix $\mathbf{\Theta}$ is thus estimated by optimizing an equivalent of a regularized least squares function

[TABLE]

where $\bm{\vartheta}=[\bm{\vartheta}^{(i)}]_{i=1,2,\ldots,n}^{T},\ \mathbf{X}=[\bm{x}^{(i)}]_{i=1,2,\ldots,n}^{T},\ \|.\|_{F}$ is the Frobenius norm of a matrix and $\lambda$ is the regularization parameter. The result is a closed form estimate given by $\mathbf{\hat{\Theta}}=(\mathbf{X}^{T}\mathbf{X}+\lambda\mathbf{I})^{-1}\mathbf{X}^{T}\bm{\vartheta}.$ To capture non-linearities which may be present, a dual representation $\mathbf{\Theta}=\mathbf{X}^{T}\bm{\Gamma}$ is taken so that

[TABLE]

then $U(\mathbf{\Gamma})$ is optimized to get $\hat{\bm{\Gamma}}=(\mathbf{K}+\lambda\mathbf{I})^{-1}\bm{\vartheta}$, where $\mathbf{K}=\mathbf{XX}^{T}$ is the Gram matrix, $K_{ij}=\langle\bm{x}^{(i)},\bm{x}^{(j)}\rangle+1$. However a more general kernel $K_{ij}=((\phi(\bm{x}^{(i)}),\phi(\bm{x}^{(j)}))$ where $\phi(\cdot)$is a nonlinear mapping, is preferred in practice.

The final step of the algorithm involves prediction/classification given a new set of response variables $\bm{x}^{new}$. This entails estimation of $\bm{\vartheta}^{new}$ by $\bm{\vartheta}^{new}=\hat{\mathbf{\Gamma}}^{T}\hat{\bm{x}}^{new}$, but $\hat{\bm{x}}^{new}=K((\phi(\bm{x}^{(i)}),\phi(\bm{x}^{(new)}))$. The computed $\bm{\vartheta}^{new}$ is used to estimate $\mathbf{t}^{new}$ by using $\bm{t}^{new}=\text{logit}^{-1}(\bm{\vartheta}^{new})$. The inverse of a logit function is given by

[TABLE]

The class labels associated with $\bm{x}^{new}$ are then computed using the estimated conditional distribution by finding the components that maximize those of $\bm{t}^{new}$ i.e. using the Bayes rule. The computed $\mathbf{t}^{new}$ is then used to get the class label ($\bm{\hat{y}}^{new}$) of the new data; for details see Dalmau et al. (2015).

3 Partial least squares (PLS) and some of its applications in genomics

PLS is a very useful approach because it is able to analyze data with many, noisy, collinear as well as incomplete variables. PLS is usually utilized in data reduction when there is multicollinearity or when the data have more variables than the number of samples. Essentially, the PLS aims at maximizing the covariance between the response variables $\bm{Y}$ and the predictors $\bm{X}$, i.e., $cov(\bm{X}^{T}\bm{Y})$ of highly multidimensional data by finding a linear subspace of the explanatory variables (Wold et al., 2001; Höskuldsson, 1988). Some literature on PLS can be found in Wold et al. (2001, 1984); andHosk1988, among others.

The research on PLS is still very active due to its ability to address problems associated with the high dimensional data such as multicollinearity and high dimensionality, among others. In the recent past, PLS has been utilized predominantly in high dimensional data in different fields like chemometrics and the “omics” like genomics, proteomics, metabolomics, and many other fields that generate large amounts of data, like spectroscopy (Gromski et al., 2015). Recent applications of PLS in microarray studies include Huang et al. (2013), who applied PLS regression (PLSR) in breast cancer intrinsic taxonomy, for classification of distinct molecular sub-types by using PAM50 signature genes as predictive variables in PLS analysis and the latent binary gene component analyzed by a logistic regression for each molecular sub-type. Also, Telaar et al. (2013) extended the notion of PLS-discriminant analysis (PLS-DA) to Powered PLS-DA (PPLS-DA), introducing a ‘power parameter’ maximised towards the correlation between the components and the group-membership, thereby achieving a minimal classification error. Furthermore, Xi et al. (2014) discussed the PLS-DA with applications to metabolites data. Other articles involving the usage of PLS include: Dong et al. (2014) who used PLS to investigate the underlying mechanism of the post-traumatic stress disorder (PTSD) using microarray data; Gusnanto et al. (2013), who made gene selection based on partial least squares and logistic regression random-effects (RE) in classification models; gene selection involving PLS was also done by Wang et al. (2015). The sparse PLS has also been utilized by many researchers; for instance, Chun and Keles (2009); Lee et al. (2011) and Chung2010 provided an efficient algorithm for the implementation of sparse PLS for variable selection in high dimensional data. Furthermore, Lê Cao et al. (2008) used sparse PLS for variable selection when integrating omics data. They implemented sparsity via lasso penalization of the PLS loading vectors when computing the singular value decomposition.

4 PLS generalized linear regression algorithm

In this section, we present an algorithm that can be applied to any Generalized Linear Regression which was developed by Bastien et al. (2005). Consider the response data $\mathbf{y}$ with the explanatory variables $\mathbf{x_{1},\ldots,x_{p}}$; then a PLS-General Linear Regression (PLSGLR) can be written as

[TABLE]

where $\theta$ a conditional expectation of the variable $\mathbf{y}$ if its distribution is continuous, or a vector of probabilities if the variable $\mathbf{y}$ follows a discrete distribution with a finite support, while $g(.)$ is the link function chosen according to the probability distribution of $\mathbf{y}$. The PLS components are given by $t_{h}=\sum_{j=1}^{p}w_{hj}^{*}\mathbf{x}_{j},j=1,\ldots,p,h=1,\ldots,m$. To compute the PLS components, let $\mathbf{X}=\mathbf{x_{1}\ldots,x_{p}}$ be a matrix of $p$ centred explanatory variables $\mathbf{x}_{j}$’s. The key objective is to determine $m$ orthogonal PLS components defined as a linear combination of the $\mathbf{x}_{j}$’s. The algorithm is presented as follows:

Computation of the first PLS component $t_{1}$ : First, the regression coefficients $a_{1j}$ of $\mathbf{x}_{j}$ are computed using the usual GLM procedure of $\mathbf{y}$ on $\mathbf{x}_{j},j=1\ldots p$. The column vector $\mathbf{a}_{1}$ which contains $a_{1j}$ is then normalized : $\mathbf{w}_{1}=\mathbf{a}_{1}/\lVert\mathbf{a}_{1}\rVert$. Finally, the component $t_{1}$ is computed as $\mathbf{t}_{1}=\mathbf{Xw}_{1}/\mathbf{w}_{1}^{\prime}\mathbf{w}_{1}$. 2. 2.

Computation of the second PLS component $t_{2}$ : Involves the computation of the linear model coefficients $a_{2j}$ of $\mathbf{x}_{j}$ in the GLM setting of $\mathbf{y}$ on $t_{1}$ and $\mathbf{x}_{j},j=1,...,p$. Since the main idea of PLS is to create the orthogonal components $t_{2}$, the component $t_{1}$ is added as a variable in estimating $\mathbf{y}$ on $t_{1}$ and $\mathbf{x}_{j},j=1,...,p$. This is because the structure of PLSGLR does not allow the residuals of $y$ to be obtained in each iteration that would aid in construction of orthogonal components. The column vector $\mathbf{a}_{2}$ which contains $a_{2j}$ is normalized: $\mathbf{w}_{2}=\mathbf{a}_{2}/\lVert\mathbf{a}_{2}\rVert$ and thereafter, the residual matrix $\mathbf{X}_{1}$ is obtained via the regression of $\mathbf{X}$ on $t_{1}$. The use of residual matrix in the attainment of the next component ensures orthogonality between the different components. The component $t_{2}$ is calculated by $\mathbf{t}_{2}=\mathbf{X}_{1}\mathbf{w}_{2}/\mathbf{w}_{2}^{\prime}\mathbf{w}_{2}$. Finally, $\mathbf{t}_{2}$ is expressed in terms of $\mathbf{X}$ : $\mathbf{t}_{2}=\mathbf{Xw}_{2}^{*}$. 3. 3.

Computation of the $h^{th}$ PLS Component $t_{h}$: Consider the already computed components $t_{1},...,t_{h-1}$; the final component $t_{h}$ is computed by calculating the GLM coefficients $a_{hj}$ of $\mathbf{x}_{j}$ by fitting $\mathbf{y}$ on $t_{1},...t_{h-1}$ and $\mathbf{x}_{j},j=1,...,p$. Next, the column vector $\mathbf{a}_{h}$, which contains $a_{hj}$ is normalized as: $\mathbf{w}_{h}=\mathbf{a}_{h}/\lVert\mathbf{a}_{h}\rVert$. The residual matrix $\mathbf{X}_{h-1}$ of the regression of $\mathbf{X}$ on $t_{1},...,t_{h-1}$ is then computed. The use of the residual matrix and the previously obtained $t_{1},...t_{h-1}$ as covariables in calculating the GLM coefficients helps the creation of orthogonal components, as previously explained. The final component $t_{h}$ is thus computed as $\mathbf{t}_{h}=\mathbf{X}_{h-1}\mathbf{w}_{h}/\mathbf{w}_{h}^{\prime}\mathbf{w}_{h}$. Finally, $\mathbf{t}_{h}$ is expressed in terms of $\mathbf{X}$ : $\mathbf{t}_{h}=\mathbf{Xw}_{h}^{*}$.

Bastien et al. (2005) note that while computing the components $t_{h}$, the nonsignificant elements in $a_{h}$ can be set to zero in order to simplify calculations, since only the significant response variables are needed to build the PLS components. The number of $m$ components to be used can be determined through cross-validation or by hard thresholding. The iteration can be stopped once there are no more significant coefficients in $a_{h}$.

Consider $x_{h-1,i}$, a column vector of the transpose of the $i$th row of $X_{h-1}$; then $t_{hi}=x^{\prime}_{h-1,i}w_{h}/w^{\prime}_{h}w_{h}$ of the $i$th case on the component $t_{h}$. This is basically the slope of the fitted line of the univariate OLS linear regression without intercept for $x_{h-1,i}$ on $w_{h}$, which can be estimated even with some data missing. Consequently, the component is computed based on the available data. Therefore the PLSGLR algorithm by Bastien et al. (2005) effectively copes up with missing data.

5 Applications to real data sets

5.1 Some exploratory analysis

In this study, we will describe in detail the analysis of the Colon data by Alon et al. (1999), obtained from the R package plsgenomics, which consists of a $(62\times 2000)$ matrix giving the expression levels of 2000 genes for 62 colon tissue samples.

An exploratory analysis of the data was done in order to visualize the differences in the un-preprocessed and preprocessed microarray data sets. The preprocessing is done using the R package plsgenomics see https://rdrr.io/cran/plsgenomics/, that implements the recommendations of (Dudoit et al., 2002). To visualize the differences between the preprocessed and un-preprocessed data sets, we consider the pairs of box plots, relative log expression (RLE), and principal components analysis (PCA) plots presented in Figures 1,2, 3, 4 and 5 respectively.

The same pair of data sets is examined using RLE plots, to show how the preprocessed data compares with the un-preprocessed data set with regard to the batch effect or any other abnormality. The RLE plots have been extensively used in studies of microarray data to reveal the effectiveness of data normalization; for an example see Gagnon-Bartsch and Speed (2011). The RLE plots are simple yet very powerful in the visualization of data to detect unwanted variations. To understand how an RLE plot is constructed, consider a data matrix $\bm{X}_{p\times n}$ where $p$ is the number of genes while $n$ the number of microarray samples, and so the element of the data matrix $x_{ij}$ represents the $i^{th}$ gene in the $j^{th}$sample. The RLE plot is then constructed by first calculating the median across each of the $p$ rows, and then substracting the respective median across each row of the data matrix $\bm{X}$, i.e $(x_{ij}-\text{median}\ x_{i*})$. The median is used because it is robust and not affected by outliers. A box plot is then generated for each of the $n$ samples, and a good one will be centered around zero and its width (interquartile range) should be equal to or less than $0.2$ (see (Gagnon-Bartsch and Speed, 2011)).

Finally we compare the ease of classification between the un-preprocessed and preprocessed data. The simplest way to visualize the separability of categories in a given data set is the use of principal components analysis (PCA) plots.

Gagnon-Bartsch and Speed (2011) note that one of the key challenges of the removal of unwanted variation is the difficulty in distinguishing the unwanted variations from the biological variation of interest. Furthermore, they note that the most appropriate way to deal with unwanted variation depends on the final objective of the analysis, for instance: differential expression (DE), classification, or clustering.

5.2 Analysis of the un-preprocessed data

In this analysis, we compare the performance of our proposed model extensions PLSGLR-log, PLSGLRDA and the KMA Dalmau et al. (2015) to that of the classical methods when the data has neither been preprocessed nor variables been selected, thus testing the performance of the classification algorithms in the presence of noise. The performance of the methodologies is then compared using a 10 fold cross validation (10-CV) and the corresponding classification error rates are computed. The results are presented in Table 2.

A particular method is judged to be the “best” if it has a lower classification error rate relative to the other methods, otherwise it is a poor classifier. The results based on minimal cross validation classification error rates indicate that for the Colon data, the KMA emerges as the best, followed by PLSDA, and RPLS, while the worst were KNN and PLSGLR-log.

5.3 Analysis of preprocessed data

During the preprocessing of microarray data the feature selection step is usually performed. This is because out of the thousands of variables (genes) generated, only a handful may play an important role towards the biological problem of interest. The thousands of data points are likely to be noisy due to biological or technical reasons. Thus the feature selection extracts a subset of the genes that are most informative (optimum subset of features). This reduces the noise by removing irrelevant or redundant features (Awada et al., 2012; Dudoit et al., 2002). Most commonly used feature selection methods involve ranking the genes based on some value of a univariate statistic, like the $t$-statistic, the F-statistic, or the Wilcoxon and Kruskal-Wallis statistics. A cut-off point based on either the number of genes or the p-value is imposed, to determine the number of variables to be used. Dudoit et al. (2002) suggest a gene selection method based on ranking. This is achieved by finding the ratio of between-group to within-group sum of squares $(BSS/WSS)$ so that for a gene $j$,

[TABLE]

where $\bar{x}_{.j}$ and $\bar{x}_{kj}$ are the average expression levels of gene $j$ and across all samples in class $k$, respectively. The $p$ genes with the biggest ratio are selected. In this study, we adopted the Dudoit et al. (2002) method of feature selection.

The preprocessing and the gene selection were performed using the recommendations of Dudoit et al. (2002). The top $p$ genes were thus selected using Equation 6 for the implementation of the classification methods.

The classification error rates for the various methodologies when applied to the data under consideration are presented in Table 3.

The results indicate that KMA was the best, followed by RPLS, PLSDA. PLSGRDA performed equally well, while KNN emerged as the worst classifier, also in every comparison.

6 Summary and conclusions

In this study, two extensions of the PLSGLR were considered in addition to the KMA for a comparative study with some classical classification methodologies, namely KNN, LDA, PLSDA, RPLS and SVM, when applied to one commonly used microarray data set. The data were considered when un-preprocessed and when preprocessed. For both the un-preprocessed and preprocessed cases, the KMA emerged as a clear “winner” based on lower classification error rates. The KMA algorithm can therefore be recommended for classification problems involving noisy and non-noisy data. This could be due to the fact that the chosen kernels map the samples to a higher dimensional space, where they become linearly separable. This leads to a better classification ability by the KMA. Furthermore, the three new algorithms can therefore be considered as an addition to the existing literature for the microarray data classification problems.

Acknowledgements

We acknowledge the partial support from the Mexico’s Consejo Nacional de Ciencias y Tecnología (CONACyT) project number 252996.