間違って解釈している場合はあると思うが、公開してしまう。
set.seed(0)
library(pls)
##
## 次のパッケージを付け加えます: 'pls'
## 以下のオブジェクトは 'package:stats' からマスクされています:
##
## loadings
model <- plsr(octane~., data = gasoline, validation = "CV")
summary(model)
## Data: X dimension: 60 401
## Y dimension: 60 1
## Fit method: kernelpls
## Number of components considered: 53
##
## VALIDATION: RMSEP
## Cross-validated using 10 random segments.
## (Intercept) 1 comps 2 comps 3 comps 4 comps 5 comps 6 comps
## CV 1.543 1.350 0.3852 0.2467 0.2373 0.2324 0.2241
## adjCV 1.543 1.347 0.3796 0.2467 0.2360 0.2293 0.2200
## 7 comps 8 comps 9 comps 10 comps 11 comps 12 comps 13 comps
## CV 0.2113 0.2165 0.2242 0.2388 0.2543 0.2700 0.2802
## adjCV 0.2070 0.2129 0.2189 0.2332 0.2459 0.2597 0.2700
## 14 comps 15 comps 16 comps 17 comps 18 comps 19 comps 20 comps
## CV 0.2799 0.2764 0.2927 0.3091 0.3269 0.3369 0.3405
## adjCV 0.2701 0.2665 0.2812 0.2965 0.3122 0.3217 0.3241
## 21 comps 22 comps 23 comps 24 comps 25 comps 26 comps 27 comps
## CV 0.3362 0.3145 0.3062 0.3106 0.3095 0.3148 0.3057
## adjCV 0.3200 0.2991 0.2914 0.2955 0.2945 0.2993 0.2904
## 28 comps 29 comps 30 comps 31 comps 32 comps 33 comps 34 comps
## CV 0.3044 0.3067 0.3122 0.3103 0.3106 0.3092 0.3089
## adjCV 0.2891 0.2912 0.2963 0.2944 0.2947 0.2934 0.2930
## 35 comps 36 comps 37 comps 38 comps 39 comps 40 comps 41 comps
## CV 0.3096 0.3115 0.3107 0.3104 0.3108 0.3107 0.3108
## adjCV 0.2938 0.2955 0.2948 0.2945 0.2949 0.2948 0.2949
## 42 comps 43 comps 44 comps 45 comps 46 comps 47 comps 48 comps
## CV 0.3108 0.3108 0.3108 0.3108 0.3108 0.3108 0.3108
## adjCV 0.2948 0.2949 0.2949 0.2949 0.2949 0.2949 0.2949
## 49 comps 50 comps 51 comps 52 comps 53 comps
## CV 0.3108 0.3108 0.3108 0.3108 0.3108
## adjCV 0.2949 0.2949 0.2949 0.2949 0.2949
##
## TRAINING: % variance explained
## 1 comps 2 comps 3 comps 4 comps 5 comps 6 comps 7 comps 8 comps
## X 70.97 78.56 86.15 95.40 96.12 96.97 97.32 98.10
## octane 31.90 94.66 97.71 98.01 98.68 98.93 99.06 99.11
## 9 comps 10 comps 11 comps 12 comps 13 comps 14 comps 15 comps
## X 98.32 98.71 98.84 99.00 99.21 99.46 99.52
## octane 99.20 99.24 99.36 99.44 99.49 99.51 99.58
## 16 comps 17 comps 18 comps 19 comps 20 comps 21 comps 22 comps
## X 99.57 99.64 99.68 99.76 99.78 99.82 99.84
## octane 99.65 99.69 99.78 99.81 99.86 99.89 99.92
## 23 comps 24 comps 25 comps 26 comps 27 comps 28 comps 29 comps
## X 99.88 99.91 99.92 99.93 99.94 99.95 99.96
## octane 99.93 99.94 99.95 99.97 99.98 99.99 99.99
## 30 comps 31 comps 32 comps 33 comps 34 comps 35 comps 36 comps
## X 99.96 99.97 99.97 99.98 99.98 99.98 99.98
## octane 99.99 100.00 100.00 100.00 100.00 100.00 100.00
## 37 comps 38 comps 39 comps 40 comps 41 comps 42 comps 43 comps
## X 99.99 99.99 99.99 99.99 99.99 99.99 99.99
## octane 100.00 100.00 100.00 100.00 100.00 100.00 100.00
## 44 comps 45 comps 46 comps 47 comps 48 comps 49 comps 50 comps
## X 99.99 99.99 100 100 100 100 100
## octane 100.00 100.00 100 100 100 100 100
## 51 comps 52 comps 53 comps
## X 100 100 100
## octane 100 100 100
plot(model, "validation")
#comp数の選択
#selectNcomp()関数を用いて適切なcomp数を算出する
ncomp.onesigma <- selectNcomp(model, method = "onesigma", plot = TRUE, ylim = c(.1, 1))
ncomp.onesigma
## [1] 4
#選択されたcomp数でのcoefを出す
# model$coefficients[, , ncomp.onesigma] #長いのでコメントアウト
library(plsVarSel)
#さっきはcomp数が少ないほうが良いことも考慮してcomp数を算出したが、今回はmin値の所を使うことにする
#(なぜ?数学的なところはわかってない)
ncomp <- which.min(model$validation$PRESS)
#vipの算出。VIP関数を用いる
#プロットすると右の方にたくさんあるのがわかる
vip <- VIP(model, ncomp)
vip
## 900 nm 902 nm 904 nm 906 nm 908 nm 910 nm 912 nm 914 nm
## 0.2775927 0.3031609 0.3261057 0.3901998 0.4105373 0.4426345 0.4860273 0.4307042
## 916 nm 918 nm 920 nm 922 nm 924 nm 926 nm 928 nm 930 nm
## 0.3358803 0.2194551 0.1899314 0.2102109 0.2134226 0.2068380 0.1904668 0.1790714
## 932 nm 934 nm 936 nm 938 nm 940 nm 942 nm 944 nm 946 nm
## 0.1730036 0.1723810 0.1687137 0.1717738 0.1749492 0.1812564 0.1913403 0.1882833
## 948 nm 950 nm 952 nm 954 nm 956 nm 958 nm 960 nm 962 nm
## 0.2004638 0.1979627 0.1995020 0.1830608 0.1882648 0.2041717 0.2051845 0.1962603
## 964 nm 966 nm 968 nm 970 nm 972 nm 974 nm 976 nm 978 nm
## 0.2174082 0.2075570 0.2104023 0.2200553 0.2225669 0.2213238 0.2349984 0.2308110
## 980 nm 982 nm 984 nm 986 nm 988 nm 990 nm 992 nm 994 nm
## 0.2565421 0.2667910 0.2587467 0.2608035 0.2655957 0.2605676 0.2798569 0.2815419
## 996 nm 998 nm 1000 nm 1002 nm 1004 nm 1006 nm 1008 nm 1010 nm
## 0.2823179 0.2872692 0.2713798 0.2778332 0.2850380 0.2869408 0.2872291 0.2957501
## 1012 nm 1014 nm 1016 nm 1018 nm 1020 nm 1022 nm 1024 nm 1026 nm
## 0.2983284 0.3068966 0.3112595 0.3068846 0.2984879 0.2822172 0.2643655 0.2580385
## 1028 nm 1030 nm 1032 nm 1034 nm 1036 nm 1038 nm 1040 nm 1042 nm
## 0.2393238 0.2178560 0.2056439 0.1960317 0.1909958 0.1905676 0.1905706 0.1913541
## 1044 nm 1046 nm 1048 nm 1050 nm 1052 nm 1054 nm 1056 nm 1058 nm
## 0.1966565 0.2012427 0.1955984 0.2040601 0.2032062 0.2034688 0.2059207 0.2107631
## 1060 nm 1062 nm 1064 nm 1066 nm 1068 nm 1070 nm 1072 nm 1074 nm
## 0.2157429 0.2254961 0.2277107 0.2341726 0.2380540 0.2425892 0.2384115 0.2379669
## 1076 nm 1078 nm 1080 nm 1082 nm 1084 nm 1086 nm 1088 nm 1090 nm
## 0.2284931 0.2158308 0.2089843 0.2022033 0.2063344 0.2045436 0.2033534 0.1988961
## 1092 nm 1094 nm 1096 nm 1098 nm 1100 nm 1102 nm 1104 nm 1106 nm
## 0.1996389 0.1931753 0.1919456 0.1784152 0.1747307 0.1836130 0.1737369 0.1793558
## 1108 nm 1110 nm 1112 nm 1114 nm 1116 nm 1118 nm 1120 nm 1122 nm
## 0.1707064 0.1615564 0.1452457 0.1354978 0.1233921 0.1183963 0.1176019 0.1296653
## 1124 nm 1126 nm 1128 nm 1130 nm 1132 nm 1134 nm 1136 nm 1138 nm
## 0.1528608 0.1824765 0.2433137 0.3068752 0.3962746 0.5290200 0.6720560 0.8161248
## 1140 nm 1142 nm 1144 nm 1146 nm 1148 nm 1150 nm 1152 nm 1154 nm
## 0.9211033 0.9737832 0.9677522 0.9300330 0.8592872 0.7829940 0.6435099 0.4430165
## 1156 nm 1158 nm 1160 nm 1162 nm 1164 nm 1166 nm 1168 nm 1170 nm
## 0.2320939 0.2124867 0.3819526 0.4789085 0.5398725 0.5392261 0.5078064 0.4867062
## 1172 nm 1174 nm 1176 nm 1178 nm 1180 nm 1182 nm 1184 nm 1186 nm
## 0.5214534 0.5573527 0.6308891 0.6748382 0.7348030 0.7915742 0.9056629 1.0966365
## 1188 nm 1190 nm 1192 nm 1194 nm 1196 nm 1198 nm 1200 nm 1202 nm
## 1.4053538 1.7190034 1.9802975 1.8517752 1.2902405 1.1725342 1.9719580 2.7105622
## 1204 nm 1206 nm 1208 nm 1210 nm 1212 nm 1214 nm 1216 nm 1218 nm
## 3.1481006 3.3203835 3.2997525 3.2212928 3.0743068 2.9301012 2.6847044 2.4089456
## 1220 nm 1222 nm 1224 nm 1226 nm 1228 nm 1230 nm 1232 nm 1234 nm
## 2.1187125 1.8357523 1.6288461 1.4053006 1.2351583 1.0604709 0.9464608 0.8555358
## 1236 nm 1238 nm 1240 nm 1242 nm 1244 nm 1246 nm 1248 nm 1250 nm
## 0.7519073 0.6693124 0.6190502 0.5426565 0.4912462 0.4355729 0.3896690 0.3535458
## 1252 nm 1254 nm 1256 nm 1258 nm 1260 nm 1262 nm 1264 nm 1266 nm
## 0.3262476 0.3104856 0.2958596 0.2834389 0.2706242 0.2631681 0.2558319 0.2509942
## 1268 nm 1270 nm 1272 nm 1274 nm 1276 nm 1278 nm 1280 nm 1282 nm
## 0.2456936 0.2396629 0.2377048 0.2285042 0.2246010 0.2114436 0.2068035 0.2042099
## 1284 nm 1286 nm 1288 nm 1290 nm 1292 nm 1294 nm 1296 nm 1298 nm
## 0.1965025 0.1949848 0.1895823 0.1804870 0.1746245 0.1689881 0.1674122 0.1705444
## 1300 nm 1302 nm 1304 nm 1306 nm 1308 nm 1310 nm 1312 nm 1314 nm
## 0.1686110 0.1711039 0.1665947 0.1758912 0.1718508 0.1786189 0.1733876 0.1808341
## 1316 nm 1318 nm 1320 nm 1322 nm 1324 nm 1326 nm 1328 nm 1330 nm
## 0.1800080 0.1821854 0.1754082 0.1667783 0.1552978 0.1469117 0.1436292 0.1425013
## 1332 nm 1334 nm 1336 nm 1338 nm 1340 nm 1342 nm 1344 nm 1346 nm
## 0.1393666 0.1358602 0.1294355 0.1332246 0.1519338 0.1716226 0.2194534 0.2903217
## 1348 nm 1350 nm 1352 nm 1354 nm 1356 nm 1358 nm 1360 nm 1362 nm
## 0.3945911 0.5488020 0.7490226 0.9859314 1.2263247 1.4697384 1.6741200 1.8069878
## 1364 nm 1366 nm 1368 nm 1370 nm 1372 nm 1374 nm 1376 nm 1378 nm
## 1.8728065 1.9027164 1.9419905 1.9946274 2.0841138 2.1318083 2.0769073 2.0224221
## 1380 nm 1382 nm 1384 nm 1386 nm 1388 nm 1390 nm 1392 nm 1394 nm
## 1.8779534 1.7663843 1.6495505 1.5399997 1.4425053 1.3805055 1.3827607 1.4252602
## 1396 nm 1398 nm 1400 nm 1402 nm 1404 nm 1406 nm 1408 nm 1410 nm
## 1.4099603 1.3087945 1.1321591 0.9278557 0.7538290 0.6315712 0.6218756 0.6808989
## 1412 nm 1414 nm 1416 nm 1418 nm 1420 nm 1422 nm 1424 nm 1426 nm
## 0.7421757 0.7769044 0.8179101 0.8585295 0.8836795 0.9043375 0.8739804 0.8105936
## 1428 nm 1430 nm 1432 nm 1434 nm 1436 nm 1438 nm 1440 nm 1442 nm
## 0.7142704 0.6143313 0.5440115 0.5224638 0.5261790 0.5272218 0.5193421 0.4957662
## 1444 nm 1446 nm 1448 nm 1450 nm 1452 nm 1454 nm 1456 nm 1458 nm
## 0.4689173 0.4447327 0.4290226 0.4158464 0.4001543 0.3971178 0.3939996 0.3934582
## 1460 nm 1462 nm 1464 nm 1466 nm 1468 nm 1470 nm 1472 nm 1474 nm
## 0.3934421 0.3952749 0.3973631 0.3955701 0.3882707 0.3606490 0.3565064 0.3437247
## 1476 nm 1478 nm 1480 nm 1482 nm 1484 nm 1486 nm 1488 nm 1490 nm
## 0.3285122 0.3125687 0.3162551 0.3136561 0.3172010 0.3368234 0.3544416 0.3901589
## 1492 nm 1494 nm 1496 nm 1498 nm 1500 nm 1502 nm 1504 nm 1506 nm
## 0.4163884 0.4525728 0.4886232 0.4930304 0.4634403 0.4099207 0.3427633 0.3039732
## 1508 nm 1510 nm 1512 nm 1514 nm 1516 nm 1518 nm 1520 nm 1522 nm
## 0.2771592 0.2593871 0.2560647 0.2535287 0.2541353 0.2702627 0.2784899 0.2904819
## 1524 nm 1526 nm 1528 nm 1530 nm 1532 nm 1534 nm 1536 nm 1538 nm
## 0.2900829 0.2941141 0.3095008 0.2936192 0.3159499 0.3102032 0.3107711 0.2922935
## 1540 nm 1542 nm 1544 nm 1546 nm 1548 nm 1550 nm 1552 nm 1554 nm
## 0.2652497 0.2625761 0.2367097 0.2236757 0.2303524 0.2216735 0.2177601 0.2082711
## 1556 nm 1558 nm 1560 nm 1562 nm 1564 nm 1566 nm 1568 nm 1570 nm
## 0.2122006 0.2119241 0.2042510 0.2046693 0.2115436 0.2237222 0.2314283 0.2120457
## 1572 nm 1574 nm 1576 nm 1578 nm 1580 nm 1582 nm 1584 nm 1586 nm
## 0.2010191 0.2088539 0.2142081 0.2175453 0.2187649 0.2340833 0.2375675 0.2510524
## 1588 nm 1590 nm 1592 nm 1594 nm 1596 nm 1598 nm 1600 nm 1602 nm
## 0.2626914 0.2565681 0.2649631 0.2650900 0.2677792 0.2615896 0.2943406 0.2855536
## 1604 nm 1606 nm 1608 nm 1610 nm 1612 nm 1614 nm 1616 nm 1618 nm
## 0.3074166 0.3051334 0.3370570 0.3638774 0.3851392 0.4166682 0.4526956 0.4622451
## 1620 nm 1622 nm 1624 nm 1626 nm 1628 nm 1630 nm 1632 nm 1634 nm
## 0.4523225 0.3297532 0.1964628 0.3689999 0.7857965 1.4903056 1.9234985 2.3958122
## 1636 nm 1638 nm 1640 nm 1642 nm 1644 nm 1646 nm 1648 nm 1650 nm
## 2.5952856 2.6432551 2.3672764 2.0936969 1.8114199 1.6312617 1.5371581 1.5486927
## 1652 nm 1654 nm 1656 nm 1658 nm 1660 nm 1662 nm 1664 nm 1666 nm
## 1.6256073 1.6889816 1.8278971 2.0879135 2.2646793 2.4444846 2.7937172 3.0306645
## 1668 nm 1670 nm 1672 nm 1674 nm 1676 nm 1678 nm 1680 nm 1682 nm
## 3.2130487 3.2641719 3.1206490 3.1032153 2.8949577 2.6172434 2.2305012 1.8245969
## 1684 nm 1686 nm 1688 nm 1690 nm 1692 nm 1694 nm 1696 nm 1698 nm
## 1.4775513 1.3951413 1.0097740 1.6268817 1.7226562 1.9024337 1.4982248 1.0659153
## 1700 nm
## 1.2536950
plot(vip, model$coefficients[, , ncomp.onesigma])
#デフォルトではvipが1.0以上を使うような形で絞り込む
#vipの絞り込みはbve_pls()関数を用いる
y <- gasoline$octane
X <- gasoline$NIR
vip.selected <- bve_pls(y, X, ncomp = ncomp)
vip.selected
## $bve.selection
## [1] 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162
## [20] 163 164 165 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243
## [39] 244 245 246 247 248 249 250 251 366 367 368 369 370 371 372 373 374 375 376
## [58] 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395
## [77] 396 397 398 399 400 401
#絞り込んだ変数で再度plsにかける
new_data <- data.frame(octane = gasoline$octane,
NIR = gasoline$NIR[, vip.selected$bve.selection])
new_model <- plsr(octane~., data = new_data, validation = "CV")
#以下上と同じ処理
summary(new_model)
## Data: X dimension: 60 82
## Y dimension: 60 1
## Fit method: kernelpls
## Number of components considered: 53
##
## VALIDATION: RMSEP
## Cross-validated using 10 random segments.
## (Intercept) 1 comps 2 comps 3 comps 4 comps 5 comps 6 comps
## CV 1.543 1.281 0.4205 0.3135 0.2730 0.2562 0.2542
## adjCV 1.543 1.281 0.4149 0.3120 0.2712 0.2528 0.2501
## 7 comps 8 comps 9 comps 10 comps 11 comps 12 comps 13 comps
## CV 0.2581 0.2642 0.2798 0.2843 0.2954 0.3018 0.3005
## adjCV 0.2529 0.2583 0.2725 0.2769 0.2869 0.2926 0.2912
## 14 comps 15 comps 16 comps 17 comps 18 comps 19 comps 20 comps
## CV 0.3074 0.3172 0.3401 0.3686 0.4072 0.4302 0.4421
## adjCV 0.2973 0.3064 0.3274 0.3547 0.3900 0.4115 0.4228
## 21 comps 22 comps 23 comps 24 comps 25 comps 26 comps 27 comps
## CV 0.4392 0.4680 0.4583 0.4614 0.4714 0.4783 0.4728
## adjCV 0.4193 0.4462 0.4368 0.4393 0.4491 0.4554 0.4503
## 28 comps 29 comps 30 comps 31 comps 32 comps 33 comps 34 comps
## CV 0.4610 0.4729 0.4900 0.4942 0.4826 0.4740 0.4771
## adjCV 0.4389 0.4500 0.4659 0.4696 0.4585 0.4503 0.4530
## 35 comps 36 comps 37 comps 38 comps 39 comps 40 comps 41 comps
## CV 0.4721 0.4664 0.4654 0.4659 0.4671 0.4684 0.4688
## adjCV 0.4482 0.4428 0.4416 0.4420 0.4432 0.4444 0.4447
## 42 comps 43 comps 44 comps 45 comps 46 comps 47 comps 48 comps
## CV 0.4696 0.4701 0.4706 0.4710 0.4712 0.4713 0.4714
## adjCV 0.4455 0.4460 0.4464 0.4468 0.4470 0.4471 0.4472
## 49 comps 50 comps 51 comps 52 comps 53 comps
## CV 0.4714 0.4714 0.4714 0.4714 0.4714
## adjCV 0.4472 0.4472 0.4472 0.4472 0.4472
##
## TRAINING: % variance explained
## 1 comps 2 comps 3 comps 4 comps 5 comps 6 comps 7 comps 8 comps
## X 73.27 82.38 93.03 95.26 96.28 97.13 97.94 98.28
## octane 31.88 93.60 96.36 97.76 98.43 98.73 98.88 98.96
## 9 comps 10 comps 11 comps 12 comps 13 comps 14 comps 15 comps
## X 98.58 98.96 99.11 99.37 99.49 99.56 99.62
## octane 99.02 99.06 99.15 99.20 99.25 99.30 99.34
## 16 comps 17 comps 18 comps 19 comps 20 comps 21 comps 22 comps
## X 99.65 99.73 99.77 99.82 99.88 99.91 99.93
## octane 99.38 99.41 99.50 99.55 99.58 99.66 99.70
## 23 comps 24 comps 25 comps 26 comps 27 comps 28 comps 29 comps
## X 99.94 99.94 99.96 99.96 99.97 99.98 99.98
## octane 99.74 99.78 99.80 99.82 99.84 99.87 99.89
## 30 comps 31 comps 32 comps 33 comps 34 comps 35 comps 36 comps
## X 99.98 99.99 99.99 99.99 99.99 99.99 99.99
## octane 99.92 99.94 99.95 99.96 99.98 99.98 99.99
## 37 comps 38 comps 39 comps 40 comps 41 comps 42 comps 43 comps
## X 100.00 100 100 100 100 100 100
## octane 99.99 100 100 100 100 100 100
## 44 comps 45 comps 46 comps 47 comps 48 comps 49 comps 50 comps
## X 100 100 100 100 100 100 100
## octane 100 100 100 100 100 100 100
## 51 comps 52 comps 53 comps
## X 100 100 100
## octane 100 100 100
plot(new_model, "validation")
new_ncomp.onesigma <- selectNcomp(new_model, method = "onesigma", plot = TRUE, ylim = c(.1, 1))
new_ncomp.onesigma
## [1] 4
#new_model$coefficients[, , new_ncomp.onesigma] #長いのでコメントアウト
library(gaselect)
ctrl <- genAlgControl(populationSize = 100, numGenerations = 30, minVariables = 5,
maxVariables = 20, verbosity = 1)
#Single repeated CV のやり方
#outerSegments に変数を与えないとSRCVになる。
evaluatorSRCV <- evaluatorPLS(numReplications = 2, innerSegments = 7, testSetSize = 0.4,
numThreads = 1)
#Repeated Double CV のやり方
#outerSegments に1以上を入れるとRDになる。
evaluatorRDCV <- evaluatorPLS(numReplications = 2, innerSegments = 5, outerSegments = 3,
numThreads = 1)
#Xはmatrix型にする必要がある。
#よくわからないが、ちょっと回りくどいことをしないとASISになってしまう。
y <- gasoline$octane
X <- as.matrix.data.frame(gasoline$NIR)
resultSRCV <- genAlg(y, X, control = ctrl, evaluator = evaluatorSRCV, seed = 0)
## Generating initial population
## Generating generation 1
## Generating generation 2
## Generating generation 3
## Generating generation 4
## Generating generation 5
## Generating generation 6
## Generating generation 7
## Generating generation 8
## Generating generation 9
## Generating generation 10
## Generating generation 11
## Generating generation 12
## Generating generation 13
## Generating generation 14
## Generating generation 15
## Generating generation 16
## Generating generation 17
## Generating generation 18
## Generating generation 19
## Generating generation 20
## Generating generation 21
## Generating generation 22
## Generating generation 23
## Generating generation 24
## Generating generation 25
## Generating generation 26
## Generating generation 27
## Generating generation 28
## Generating generation 29
## Generating generation 30
resultRDCV <- genAlg(y, X, control = ctrl, evaluator = evaluatorRDCV, seed = 0)
## Generating initial population
## Generating generation 1
## Generating generation 2
## Generating generation 3
## Generating generation 4
## Generating generation 5
## Generating generation 6
## Generating generation 7
## Generating generation 8
## Generating generation 9
## Generating generation 10
## Generating generation 11
## Generating generation 12
## Generating generation 13
## Generating generation 14
## Generating generation 15
## Generating generation 16
## Generating generation 17
## Generating generation 18
## Generating generation 19
## Generating generation 20
## Generating generation 21
## Generating generation 22
## Generating generation 23
## Generating generation 24
## Generating generation 25
## Generating generation 26
## Generating generation 27
## Generating generation 28
## Generating generation 29
## Generating generation 30
#変数選択のベスト5...ってコト?
subsets(resultSRCV, 1:5)
## $`1`
## [1] 17 21 49 62 96 120 124 168 188 234 265 286 378
##
## $`2`
## [1] 13 97 120 124 168 186 234 265 286 378 395
##
## $`3`
## [1] 17 41 69 97 120 124 168 186 234 265 286 378
##
## $`4`
## [1] 17 41 97 120 124 168 186 234 265 286 378
##
## $`5`
## [1] 13 69 97 120 124 168 186 234 265 286 378 395
subsets(resultRDCV, 1:5)
## $`1`
## [1] 64 119 154 163 224 237 284 370
##
## $`2`
## [1] 64 119 154 163 210 237 284 370
##
## $`3`
## [1] 64 119 154 163 210 224 237 284 370
##
## $`4`
## [1] 64 119 154 163 198 210 224 237 284 370
##
## $`5`
## [1] 74 119 154 163 198 237 284 370
#VIPで変数選択->plsした結果のmseを計算
VIP.est.Y <- predict(new_model, newx = gasoline$NIR, ncomp = new_ncomp.onesigma)
VIP.mse <- sum((gasoline$octane - VIP.est.Y)^2) / length(VIP.est.Y)
VIP.mse
## [1] 0.05153734
plot(gasoline$octane, VIP.est.Y)
#GAでで変数選択->plsした結果のmseを計算
RDCV.list <- subsets(resultRDCV, 1)$`1`
RDCV.new_data <- data.frame(octane = gasoline$octane,
NIR = gasoline$NIR[, RDCV.list])
RDCV.model <- plsr(octane~., data = RDCV.new_data, validation = "CV")
RDCV.ncomp.onesigma <- selectNcomp(RDCV.model, method = "onesigma", plot = FALSE)
GA.est.Y <- predict(RDCV.model, newx = gasoline$NIR, ncomp = RDCV.ncomp.onesigma)
GA.mse <- sum((gasoline$octane - GA.est.Y)^2) / length(GA.est.Y)
GA.mse
## [1] 0.02993481
plot(gasoline$octane, GA.est.Y)