1
$\begingroup$

General question about predictions in PLS, using the chemometrics package in R here for an example. The standard example for pls1_nipals is

data(PAC)
res <- pls1_nipals(PAC$X,PAC$y,a=5)

This returns the coefficient vector res$b. Now lets assume the first 200 of 209 objects (rows) of PAC$X,PAC$y are used for calibration and the rest is used for validation.

res <- pls1_nipals(PAC$X[1:200,],PAC$y[1:200],a=5)
cbind(PAC$X[201:209,]%*%res$b,PAC$y[201:209])

leads to

    [,1]   [,2]
470103.8 486.81
457944.5 488.18
491483.9 495.01
487573.0 495.45
499962.8 497.66
487223.3 500.00
656960.3 501.32
661962.2 503.89
635861.9 503.91

The reason for the extreme difference is, that X and Y were centered. The centering of Y can be eliminated bei adding mean(PAC$[1:200]) to PAC$X[201:209,]%*%res$b. But how to correct the centering of X? First I thought by using colMeans(PAC$X[1:200,]), but it didn't work.

cbind((PAC$X[201:209,]
              -matrix(rep(colMeans(PAC$X[1:200,]),9),ncol=467,nrow=9,byrow=TRUE))
              %*%res$b+mean(PAC$y[1:200]),
      PAC$y[201:209])

Any ideas? Thanks.

$\endgroup$
  • 1
    $\begingroup$ Summary: I'm puzzled as well - I'd have thought your approach should work. In fact, it does work like that with pls::plsr models. I suspect you may have found a bug (in the normalization of the inner relation coefficients), and I emailed the package maintainer about it. I'll write an answer once I know more. $\endgroup$ – cbeleites supports Monica Aug 2 '16 at 14:38
2
$\begingroup$

You really found a bug, which was immediately fixed by the maintainer (version 1.4.1 on CRAN - available already as .tar.gz, win binaries will take a bit longer):

library (chemometrics)
library (pls)

train <- yarn [yarn$train,]

res <- pls1_nipals(train$NIR, train$density, a = 2)
x.centered <- scale (test$NIR, center = colMeans (train$NIR), scale = FALSE)
y.center <- mean (train$density)
y.hat <- x.centered %*% res$b + y.center

For comparison, set up the (kernel)PLS model using package pls:

res.pls <- plsr (density ~ NIR, data = train, ncomp = 2)
y.pls <- predict (res.pls, data = yarn [! yarn$train])

now we have:

> cbind (y.hat, y.pls, y.hat - y.pls)
              density       density
110 53.09713 53.09713  0.000000e+00
22  52.74115 52.74115 -8.526513e-14
31  34.67483 34.67483 -8.526513e-14
41  37.48693 37.48693 -7.105427e-15
51  32.26085 32.26085 -3.552714e-14
61  21.23447 21.23447 -6.394885e-14
71  22.18097 22.18097 -4.973799e-14 

So basically the same (up to "numeric noise")

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.