Search code examples
rplotlatticerglpls

Plot 3D in R: Classical spectral gasoline in function of octane


Load a data set comprising spectral intensities of 60 samples of gasoline at 401 wavelengths, and their octane ratings in pls package. In my example:

library(pls)

#Data set

data(gasoline)

'data.frame':   60 obs. of  2 variables:
 $ octane: num  85.3 85.2 88.5 83.4 87.9 ...
 $ NIR   : 'AsIs' num [1:60, 1:401] -0.0502 -0.0442 -0.0469 -0.0467 -0.0509 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr  "1" "2" "3" "4" ...
  .. ..$ : chr  "900 nm" "902 nm" "904 nm" "906 nm"

I've like to plot the 3D representation of NIR in function of octane, my holy grail plot is:

https://it.mathworks.com/help/stats/examples/partial-least-squares-regression-and-principal-components-regression.html

But this output plot was create in matlab and not in R. The original code is (https://it.mathworks.com/help/stats/examples/partial-least-squares-regression-and-principal-components-regression.html):

load spectra
whos NIR octane

[dummy,h] = sort(octane);
oldorder = get(gcf,'DefaultAxesColorOrder');
set(gcf,'DefaultAxesColorOrder',jet(60));
plot3(repmat(1:401,60,1)',repmat(octane(h),1,401)',NIR(h,:)');
set(gcf,'DefaultAxesColorOrder',oldorder);
xlabel('Wavelength Index'); ylabel('Octane'); axis('tight');
grid on

Please any ideas or similar function/packages in R? Thanks in advance!


Solution

  • This gives something like the Matlab plot:

    library(rgl) 
    library(pls) 
    data(gasoline) 
    str(gasoline) 
    gasoline$ID <- seq(1:length(gasoline[,1])) 
    open3d() 
    x <- gasoline$octane 
    y <- gasoline$NIR 
    z <- gasoline$ID 
    plot3d(x, 1:ncol(y), y, type = "n", xlab = "", ylab = "", zlab = "")
    cols <- rainbow(1000)[(x - min(x))/(max(x) - min(x))*999 + 1]
    for (i in seq_along(x))
      lines3d(x[i], 1:ncol(y), y[i,], col = cols[i])
    

    screenshot

    There are several differences:

    • The handedness of the coordinate system. (Octane increases towards the viewer in the rgl plot.) You can change it in rgl using par3d("userMatrix" = par3d("userMatrix") %*% diag(c(-1, 1,1,1))) if you really want the Matlab style.
    • I didn't draw axis labels. The default ones in rgl are kind of ugly, and I was too lazy to use mtext3d to draw nicer ones.
    • I didn't bother with background grids. I don't like them, but if you really want them, you can get them with grid3d.
    • Matlab doesn't show perspective. If that's what you want, use par3d(FOV = 0).
    • The aspect ratios are different. I think Matlab is using something like aspect3d(1, 1, 0.75).
    • The colors. rgl uses the R color system, so if you can find a palette you like better than rainbow(1000), use that.
    • The price.