Example 1: Palmdates data
Chemical variables
First we apply the method with the first two variables:
fibre and sorbitol to show the performance of
the package. The application of the method with two variables allows the
visualization of the hyperplane in two dimensions, in this case, a
straight line.
group = as.factor( c(rep("Spanish",11), rep("Other",10)) )
model1 = lpda(data = palmdates$conc[,1:2], group = group )First output of lpda is a matrix with the coefficients
of the hyperplane: \(a'x=b\) for
each pair-wise comparison. In this example:
model1$coef
#> Other-Spanish
#> [1,] -1.116131
#> [2,] -3.002923
#> [3,] -6.563646being -1.116 and -3.003 the coefficients of fibre and
sorbitol respectively and -6.564 the constant \(b\). We can plot the line on the points,
that represent the samples, with the following code:
plot(palmdates$conc[,1:2], col = as.numeric(group)+1, pch = 20,
main = "Palmdates example")
abline(model1$coef[3]/model1$coef[2], -model1$coef[1]/model1$coef[2], cex = 2)
legend("bottomright", c("Other","Spanish"),col = c(2,3), pch = 20, cex=0.8)
Predicted group with the model is obtained with predict
function. Confusion matrix is computed as follows. We observe that all
the samples are well classified.
pred = predict(model1)
table(pred$fitted, group)
#> group
#> Other Spanish
#> Other 10 0
#> Spanish 0 11By considering the five concentration variables, all samples are well classified as well. As a 5-dimensional plot can not be showed, we offer the possibility of seeing a plot that shows the situation of samples with respect to \(H\). Y-axis represents order in which they appear in the data matrix and X-axis distances of each sample to \(H\).
model2 = lpda(data = palmdates$conc, group = group )
plot(model2)
Another option is the application of lpda to the first
PCs scores. When data is highly correlated, two components are usually
preferred. In such case choosing as plot arguments
PCscores = TRUE we will visualize the first two PCs,
indicating in the axis the proportion of explained variance by each PC.
Moreover if there are two groups, as in this example, the optimal
hyperplane is also showed.
model3 = lpda(data = palmdates$conc, group = group, pca = TRUE, Variability = 0.7)
plot(model3, PCscores = TRUE, main = "PCA-Substances palmdates")
Spectra variables
When having data sets with more variables than individuals PCA is not
directly applicable. In lpda we have implemented the
possibility of dealing with such problem, working with the equivalences
between the PCA of the data matrix and the transposed matrix.
In palmdates$spectra there are 2050 very correlated
measurements as it is showed in the following figure:
Due to the high correlation the application of lpda to
all the spectra variables and to the first 2 PCs gives the same
solution: 0 prediction errors.
model4 <- lpda(data = palmdates$spectra, group = group)
pred = predict(model4)
table(pred$fitted, group)
#> group
#> Other Spanish
#> Other 10 0
#> Spanish 0 11model5 = lpda(data = palmdates$spectra, group = group, pca = TRUE, Variability = 0.9)
plot(model5, PCscores = TRUE, main = "Spectra palmdates")
Following example shows how predict the group of a test set that does not participate in the model. In this set we include two samples of each group.
test = c(10,11,12,13)
model6 = lpda(data = palmdates$spectra[-test,], group = group[-test], pca = TRUE,
Variability = 0.9)
pred = predict(model6, palmdates$spectra[test,])
pred$fitted
#> [1] "Spanish" "Spanish" "Other" "Other"We can also use lpdaCV function to compute the predicted
error in different test sets by crossvalidation with a specific model.
As explained before, two strategies are implemented: loo
(leave one out) and ktest that is specified in
CV argument and by default is loo. When
ktest is selected, ntest is the size of test
set (samples not used for computing the model, only to evaluate the
number of prediction errors) and R is the times the model
is computed and evaluated with different training and test sets.
lpdaCV(palmdates$spectra, group, pca = TRUE, CV = "loo")
lpdaCV(palmdates$spectra, group, pca = TRUE, CV = "ktest", ntest = 5, R = 10)CV results are so good due to the clear difference between the two groups in spectra data. We continue with CV functions in Example 2 that deals noisier data.
