In my experiment, I diluted a chemical in water and recorded its Raman spectra. Then I repeated this with different concentrations. As I can clearly see some peaks decreasing in intensity when I dilute my solution, I hope that through Principal Component Regression, I would be able to automatically figure out which wavenumbers are useful for determining the concentration of my solution.

If I understand PCA correctly, I am looking for the loadings of the first few components.

I have the following data that I prepared:

> concentration.hcl
0.001 0.001 0.001 1e-04 1e-04 1e-04 1e-04 1e-09 1e-09 1e-09 1e-09 1e-11 1e-11
1e-03 1e-03 1e-03 1e-04 1e-04 1e-04 1e-04 1e-09 1e-09 1e-09 1e-09 1e-11 1e-11
1e-11 1e-11 1e-11
1e-11 1e-11 1e-11


and spectra.hcl with each row being a reading corresponding to an entry in concentration.hcl and the columns being a specific wavenumber.

> dim(spectra.hcl)
[1]   16 2048

> colnames(spectra.hcl)
[1] "-818.41" "-813.78" "-809.16" "-804.53" "-799.92" "-795.31" "-790.7"
[8] "-786.1"  "-781.5"  "-776.91" "-772.32" "-767.74" "-763.16" "-758.59"
[15] "-754.02" "-749.45" "-744.9"  "-740.34" "-735.79" "-731.25" "-726.71" ...


I ran the following code:

pmodel <- pcr(concentration.hcl ~ spectra.hcl)

> summary(pmodel)
Data:   X dimension: 16 2048
Y dimension: 16 1
Fit method: svdpc
Number of components considered: 15
TRAINING: % variance explained
1 comps  2 comps  3 comps  4 comps  5 comps  6 comps
X                     85.7    96.91    99.29    99.64    99.81    99.92
concentration.hcl     45.4    90.14    98.35    98.72    98.73    99.27
...


> pmodel\$loadings
... (truncated)
3291.07
3291.68
3292.29
3292.89
3293.49
3294.1
3294.7
3295.3
3295.9
3296.5
3297.09
3297.69
3298.29
3298.88
3299.48
3300.07
3300.66
3301.26
3301.85
3302.44
3303.03
3303.62
3304.2
3304.79
3305.38
3305.96
3306.55
3307.14

Comp 1 Comp 2 Comp 3 Comp 4 Comp 5 Comp 6 Comp 7 Comp 8 Comp 9
Proportion Var      0  0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
Cumulative Var      0  0.001  0.001  0.002  0.002  0.003  0.003  0.004  0.004
Comp 10 Comp 11 Comp 12 Comp 13 Comp 14 Comp 15
Proportion Var   0.000   0.000   0.000   0.000   0.000   0.000
Cumulative Var   0.005   0.005   0.006   0.006   0.007   0.007

2. How can I find out which wavenumbers are contributing the most to the varying concentrations in my experiments?

NOTE: The data is available at https://gist.githubusercontent.com/anonymous/1d994685c3b3a5133b29/raw/696662a3aa199890958c162483c0cb7aa95e18ec/file1.txt You just need to download this file and then type load('file1.txt') and the variables will show up in your R environment

• What do you mean "why can't I see my loadings"? What exactly did you expect to see and why the output that you copy-pasted does not satisfy you? – amoeba Jul 8 '15 at 13:12
• You should clarify which library is loaded for pcr() function? Is it qualityTools package? – rnso Jul 8 '15 at 13:16
• If you are using the pls package, you can obtain the loadings with the command loadings(pmodel). As @amoeba suggests, you might want to expand on your objective to obtain more appropriate help. The rationale for principal components regression has always escaped me because one is doing a manipulation of the predictor variables in complete isolation of the dependent variable and somehow expecting to get better predictions especially if one only uses the principal components that "explain" most of the variation in the predictor variables. (And this should work for all possible dep. vars?) – JimB Jul 8 '15 at 14:27
• @Jim, as an aside, principal components regression is intimately related to ridge regression, and can be seen as a particular form of regularization. Nobody complains that ridge regression works "in isolation" from the dependent variable, but this is often said about PCR; in fact, they are not very much different. See e.g. The Elements of Statistical Learning for details, or search this site for relevant discussions. – amoeba Jul 8 '15 at 15:02
• @amoeba. Thanks. I will look at that. But one can certainly have a case where the principal component that explains the least amount of variation is the best predictor. For whatever it's worth, my other issue (which might also be unfounded) is that principal components (with or without a regression aspect) is many times applied completely ignoring the experimental design or how the data was collected. – JimB Jul 8 '15 at 15:14

There are 2 things going on here:

• the print method for loadings objects in package pls tidies the printed values so that close-to-zero values are zapped, and
• you look only at the last lines of the output.

E.g. if you look at the PCR example with the yarn data, the first few lines of the output are

Loadings:
Comp 1 Comp 2 Comp 3 Comp 4 Comp 5 Comp 6
NIR1                 -0.115               -0.243
NIR2                 -0.134         0.170 -0.325
NIR3                 -0.113         0.206 -0.369
NIR4                                0.154 -0.352
NIR5                                      -0.269
NIR6                                      -0.180
NIR7    0.108                             -0.114


matplot (loadings (model), type = "l")