I perform PLS with pls package in R using plsr function.

Why am I unable to get the same predicted Y value as when I use the predict function as when I dot product the regression coefficients with the test data. I have read from other answers that I should be able to do this, especially if the method is set to simpls.

EDIT: I did some further testing and it seems like there is a bias/intercept constant that I have not added in. However, from elsewhere I read that with mean centering, the intercept should be 0 and I believe plsr has already mean centered all my data by default. I am wondering why is there still an intercept and also how can I find this intercept constant from the plsr model?

EDIT again: I have found the coef(,intercept=TRUE) function! However, my question still remains, why is there an intercept when my data has been mean centered (Both X and Y)?

Many thanks for your help!! Let me know if I can clarify anything at all!


pls::plsr centers both $\mathbf X$ and $\mathbf Y$, and the corresponding intercepts are in $Xmeans and $Ymeans.

So in order to predict using the coefficients that map $\mathbf{Y_c} = \mathbf { X_c} \mathbf B$, you need to

  1. center $\mathbf X$: $\mathbf X_c = \mathbf X - \bar x$ ($\bar x$ is $Xmeans)
  2. matrix-multiply by $\mathbf B$: $\mathbf Y_c = \mathbf X_c \mathbf B$
  3. "uncenter" $\mathbf Y$: $\mathbf Y = \mathbf Y_c + \bar y$ (= add $Ymeans)
library (pls)
yarn.pls <- plsr(density ~ NIR, 6, data = yarn)

X_c <- yarn$NIR [1,, drop = FALSE] - yarn.pls$Xmeans
Y_c <-  (X_c %*% yarn.pls$coefficients [,1,])
Y_c + yarn.pls$Ymeans


      1 comps  2 comps  3 comps  4 comps  5 comps  6 comps
[1,] 90.53581 90.49171 99.13326 98.94135 99.39108 99.54403

which is the same as (up to some small numerical error)

predict (yarn.pls, yarn [1,]) [,1,]

After putting the 3 steps above together, we can separate the calculations involving the centering:

$$\mathbf Y = (\mathbf X - \bar x) \mathbf B + \bar y\\ = \mathbf X \mathbf B - \bar x \mathbf B + \bar y \\ = \mathbf X \mathbf B + \beta_0$$

with $\beta_0 = - \bar x \mathbf B + \bar y$

leaving us with a single intercept, which is the intercept from coef (intercept = TRUE)

This way:

coef1 <- coef (yarn.pls, intercept = TRUE)
X1 <- cbind (1, yarn$NIR [1,, drop = FALSE])
X1 %*% c 

we get the same prediction as above using all 6 components:

[1,] 99.54403
  • $\begingroup$ thank you for such clarity!!! but just to solidify what I have in mind, adding the intercept is what you mean by 'uncenter' Y right? $\endgroup$ – sunnydk Dec 3 '19 at 21:52
  • $\begingroup$ @sunnydk: should be clear now. One important point is that we have 3 different intercepts: two for X and Y in the "low-level" formulation, and one "high-level" intercept that combines the two (with the coefficients): for PLS it is always important to spell out which intercept is meant. $\endgroup$ – cbeleites unhappy with SX Dec 3 '19 at 22:02
  • $\begingroup$ I was just going to edit my previous comment; but yes! combining your answer with theGD's, I was able to work out that the intercept is made up of both ymeans and coef times xmeans. Thank you very much!! $\endgroup$ – sunnydk Dec 3 '19 at 22:04

For PLS, mean-centering is equivalent to adding a constant "1" as a variable to every observation. In other words, it is matter of algebra to account for mean-centering with an intercept term. One can also account for scaling too. See below:

accounting for centering and scaling

Sorry for bad notation. So, $n$ is number of variables, $B$ is regression coefficients obtained by mean-centered and scaled data, $s$ stands for stddev and $\bar{x}$ and $\bar{y}$ are the means. $B_0$ here is the intercept term. Ignoring $s$ terms gives the version for mean-centering only case.

  • $\begingroup$ thank you for reminding me about the constant '1'. However, do you mind elaborating a little more about what you mean by 'equivalent' (perhaps intuitively if possible, I can follow the above algebra but I am just trying to make the link in my mind between the bias constant and mean-centering)? Also, should your Bnewi be equals to (sy*Bi)/sxi instead, since xi is from the X matrix? $\endgroup$ – sunnydk Dec 3 '19 at 21:26
  • $\begingroup$ In PLS, the aim is to find $B$ that maximizes covariance. If you add 1's as a variable, the best (according to aim of PLS) $B$ for that variable ($B_0$ or intercept) is the one that the rest of the model will vary around, basically the average. In other words, since the unbiased estimate (very vaguely means an estimate where no error is "expected") for intercept is the average, you can exploit this directly (mean-centering). Alternatively, you can let the algorithm find the best estimate for the intercept and this corresponds to adding 1s. $\endgroup$ – theGD Dec 3 '19 at 22:19
  • $\begingroup$ thank you for your answer!! there is quite a bit to unpack here mentally, but combining your answer and cbeleities' have helped me make much more sense already :) $\endgroup$ – sunnydk Dec 3 '19 at 22:29

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.