1### R code from vignette source 'pls-manual.Rnw'
2
3###################################################
4### code chunk number 1: pls-manual.Rnw:67-69
5###################################################
6pdf.options(pointsize=10)
7options(digits = 4)
8
9
10###################################################
11### code chunk number 2: pls-manual.Rnw:304-305
12###################################################
13library(pls)
14
15
16###################################################
17### code chunk number 3: pls-manual.Rnw:332-335
18###################################################
19data(yarn)
20data(oliveoil)
21data(gasoline)
22
23
24###################################################
25### code chunk number 4: pls-manual.Rnw:343-349
26###################################################
27par(mar = c(2, 4, 0, 1) + 0.1)
28matplot(t(gasoline$NIR), type = "l", lty = 1, ylab = "log(1/R)", xaxt = "n")
29ind <- pretty(seq(from = 900, to = 1700, by = 2))
30ind <- ind[ind >= 900 & ind <= 1700]
31ind <- (ind - 898) / 2
32axis(1, ind, colnames(gasoline$NIR)[ind])
33
34
35###################################################
36### code chunk number 5: pls-manual.Rnw:358-360
37###################################################
38gasTrain <- gasoline[1:50,]
39gasTest <- gasoline[51:60,]
40
41
42###################################################
43### code chunk number 6: pls-manual.Rnw:363-364
44###################################################
45gas1 <- plsr(octane ~ NIR, ncomp = 10, data = gasTrain, validation = "LOO")
46
47
48###################################################
49### code chunk number 7: pls-manual.Rnw:369-370
50###################################################
51summary(gas1)
52
53
54###################################################
55### code chunk number 8: pls-manual.Rnw:379-380 (eval = FALSE)
56###################################################
57## plot(RMSEP(gas1), legendpos = "topright")
58
59
60###################################################
61### code chunk number 9: pls-manual.Rnw:385-387
62###################################################
63par(mar = c(4, 4, 2.5, 1) + 0.1)
64plot(RMSEP(gas1), legendpos = "topright")
65
66
67###################################################
68### code chunk number 10: pls-manual.Rnw:406-407 (eval = FALSE)
69###################################################
70## plot(gas1, ncomp = 2, asp = 1, line = TRUE)
71
72
73###################################################
74### code chunk number 11: pls-manual.Rnw:413-415
75###################################################
76par(mar = c(4, 4, 2.5, 1) + 0.1)
77plot(gas1, ncomp = 2, asp = 1, line = TRUE)
78
79
80###################################################
81### code chunk number 12: pls-manual.Rnw:429-430 (eval = FALSE)
82###################################################
83## plot(gas1, plottype = "scores", comps = 1:3)
84
85
86###################################################
87### code chunk number 13: pls-manual.Rnw:435-436
88###################################################
89plot(gas1, plottype = "scores", comps = 1:3)
90
91
92###################################################
93### code chunk number 14: pls-manual.Rnw:444-448
94###################################################
95par(mar = c(4, 4, 0.3, 1) + 0.1)
96plot(gas1, "loadings", comps = 1:2, legendpos = "topleft",
97     labels = "numbers", xlab = "nm")
98abline(h = 0)
99
100
101###################################################
102### code chunk number 15: pls-manual.Rnw:462-463
103###################################################
104explvar(gas1)
105
106
107###################################################
108### code chunk number 16: pls-manual.Rnw:469-472 (eval = FALSE)
109###################################################
110## plot(gas1, "loadings", comps = 1:2, legendpos = "topleft",
111##      labels = "numbers", xlab = "nm")
112## abline(h = 0)
113
114
115###################################################
116### code chunk number 17: pls-manual.Rnw:481-482
117###################################################
118predict(gas1, ncomp = 2, newdata = gasTest)
119
120
121###################################################
122### code chunk number 18: pls-manual.Rnw:486-487
123###################################################
124RMSEP(gas1, newdata = gasTest)
125
126
127###################################################
128### code chunk number 19: pls-manual.Rnw:587-588
129###################################################
130dens1 <- plsr(density ~ NIR, ncomp = 5, data = yarn)
131
132
133###################################################
134### code chunk number 20: pls-manual.Rnw:592-594
135###################################################
136dim(oliveoil$sensory)
137plsr(sensory ~ chemical, data = oliveoil)
138
139
140###################################################
141### code chunk number 21: pls-manual.Rnw:615-617
142###################################################
143trainind <- which(yarn$train == TRUE)
144dens2 <- update(dens1, subset = trainind)
145
146
147###################################################
148### code chunk number 22: pls-manual.Rnw:621-622
149###################################################
150dens3 <- update(dens1, ncomp = 10)
151
152
153###################################################
154### code chunk number 23: pls-manual.Rnw:652-653
155###################################################
156olive1 <- plsr(sensory ~ chemical, scale = TRUE, data = oliveoil)
157
158
159###################################################
160### code chunk number 24: pls-manual.Rnw:658-659
161###################################################
162gas2 <- plsr(octane ~ msc(NIR), ncomp = 10, data = gasTrain)
163
164
165###################################################
166### code chunk number 25: pls-manual.Rnw:664-665 (eval = FALSE)
167###################################################
168## predict(gas2, ncomp = 3, newdata = gasTest)
169
170
171###################################################
172### code chunk number 26: pls-manual.Rnw:717-721 (eval = FALSE)
173###################################################
174## ncomp.onesigma <- selectNcomp(gas2, method = "onesigma", plot = TRUE,
175##                               ylim = c(.18, .6))
176## ncomp.permut <- selectNcomp(gas2, method = "randomization", plot = TRUE,
177##                             ylim = c(.18, .6))
178
179
180###################################################
181### code chunk number 27: pls-manual.Rnw:735-740
182###################################################
183par(mfrow = c(1,2))
184ncomp.onesigma <- selectNcomp(gas1, "onesigma", plot = TRUE,
185                              ylim = c(.18, .6))
186ncomp.permut <- selectNcomp(gas1, "randomization", plot = TRUE,
187                            ylim = c(.18, .6))
188
189
190###################################################
191### code chunk number 28: pls-manual.Rnw:760-763
192###################################################
193gas2.cv <- crossval(gas2, segments = 10)
194plot(MSEP(gas2.cv), legendpos="topright")
195summary(gas2.cv, what = "validation")
196
197
198###################################################
199### code chunk number 29: pls-manual.Rnw:818-820 (eval = FALSE)
200###################################################
201## plot(gas1, plottype = "coef", ncomp=1:3, legendpos = "bottomleft",
202##      labels = "numbers", xlab = "nm")
203
204
205###################################################
206### code chunk number 30: pls-manual.Rnw:824-827
207###################################################
208par(mar = c(4, 4, 2.5, 1) + 0.1)
209plot(gas1, plottype = "coef", ncomp=1:3, legendpos = "bottomleft",
210     labels = "numbers", xlab = "nm")
211
212
213###################################################
214### code chunk number 31: pls-manual.Rnw:859-861
215###################################################
216par(mar = c(4, 4, 0, 1) + 0.1)
217plot(gas1, plottype = "correlation")
218
219
220###################################################
221### code chunk number 32: pls-manual.Rnw:930-931
222###################################################
223predict(gas1, ncomp = 2:3, newdata = gasTest[1:5,])
224
225
226###################################################
227### code chunk number 33: pls-manual.Rnw:943-944
228###################################################
229predict(gas1, comps = 2, newdata = gasTest[1:5,])
230
231
232###################################################
233### code chunk number 34: pls-manual.Rnw:961-962
234###################################################
235drop(predict(gas1, ncomp = 2:3, newdata = gasTest[1:5,]))
236
237
238###################################################
239### code chunk number 35: pls-manual.Rnw:1000-1001 (eval = FALSE)
240###################################################
241## predplot(gas1, ncomp = 2, newdata = gasTest, asp = 1, line = TRUE)
242
243
244###################################################
245### code chunk number 36: pls-manual.Rnw:1007-1009
246###################################################
247par(mar = c(4, 4, 2.5, 1))
248predplot(gas1, ncomp = 2, newdata = gasTest, asp = 1, line = TRUE)
249
250
251###################################################
252### code chunk number 37: pls-manual.Rnw:1049-1050
253###################################################
254pls.options()
255
256
257###################################################
258### code chunk number 38: pls-manual.Rnw:1056-1057
259###################################################
260pls.options(plsralg = "oscorespls")
261
262
263###################################################
264### code chunk number 39: pls-manual.Rnw:1210-1219
265###################################################
266X <- gasTrain$NIR
267Y <- gasTrain$octane
268ncomp <- 5
269cvPreds <- matrix(nrow = nrow(X), ncol = ncomp)
270for (i in 1:nrow(X)) {
271    fit <- simpls.fit(X[-i,], Y[-i], ncomp = ncomp, stripped = TRUE)
272    cvPreds[i,] <- (X[i,] - fit$Xmeans) %*% drop(fit$coefficients) +
273        fit$Ymeans
274}
275
276
277###################################################
278### code chunk number 40: pls-manual.Rnw:1222-1223
279###################################################
280sqrt(colMeans((cvPreds - Y)^2))
281
282
283