mcrals: Multivariate curve resolution using Alternating Least Squares

View source: R/mcrals.R

mcralsR Documentation

Multivariate curve resolution using Alternating Least Squares

Description

mcralls allows to resolve spectroscopic data to linear combination of individual spectra and contributions using the alternating least squares (ALS) algorithm with constraints.

Usage

mcrals(
  x,
  ncomp,
  cont.constraints = list(),
  spec.constraints = list(),
  spec.ini = matrix(runif(ncol(x) * ncomp), ncol(x), ncomp),
  cont.forced = matrix(NA, nrow(x), ncomp),
  spec.forced = matrix(NA, ncol(x), ncomp),
  cont.solver = mcrals.nnls,
  spec.solver = mcrals.nnls,
  exclrows = NULL,
  exclcols = NULL,
  verbose = FALSE,
  max.niter = 100,
  tol = 10^-6,
  info = ""
)

Arguments

x

spectra of mixtures (matrix or data frame).

ncomp

number of components to calculate.

cont.constraints

a list with constraints to be applied to contributions (see details).

spec.constraints

a list with constraints to be applied to spectra (see details).

spec.ini

a matrix with initial estimation of the pure components spectra.

cont.forced

a matrix which allows to force some of the concentration values (see details).

spec.forced

a matrix which allows to force some of the spectra values (see details).

cont.solver

which function to use as a solver for resolving of pure components contributions (see detials).

spec.solver

which function to use as a solver for resolving of pure components spectra (see detials).

exclrows

rows to be excluded from calculations (numbers, names or vector with logical values).

exclcols

columns to be excluded from calculations (numbers, names or vector with logical values).

verbose

logical, if TRUE information about every iteration will be shown.

max.niter

maximum number of iterations.

tol

tolerance, when explained variance change is smaller than this value, iterations stop.

info

a short text with description of the case (optional).

Details

The method implements the iterative ALS algorithm, where, at each iteration, spectra and contributions of each chemical component are estimated and then a set of constraints is applied to each. The method is well described in [1, 2].

The method assumes that the spectra (D) is a linear combination of pure components spectra (S) and pure component concentrations (C):

D = CS' + E

So the task is to get C and S by knowing D. In order to do that you need to provide:

1. Constraints for spectra and contributions. The constraints should be provided as a list with name of the constraint and all necessary parameters. You can see which constraints and parameters are currently supported by running constraintList(). See the code examples below or a Bookdown tutorial for more details.

2. Initial estimation of the pure components spectra, S. By default method uses a matrix with random numbers but you can provide a better guess (for example by running mcrpure) as a first step.

3. Which solver to use for resolving spectra and concentrations. There are two built in solvers: mcrals.nnls (default) and mcrals.ols. The first implements non-negative least squares method which gives non-negative (thus physically meaningful) solutions. The second is ordinary least squares and if you want to get non-negative spectra and/or contributions in this case you need to provide a non-negativity constraint.

The algorithm iteratively resolves C and S and checks how well CS' is to D. The iterations stop either when number exceeds value in max.niter or when improvements (difference between explained variance on current and previous steps) is smaller than tol value.

Parameters cont.force and spec.force allows you to force some parts of the contributions or the spectra to be equal to particular pre-defined values. In this case you need to provide the parameters (or just one of them) in form of a matrix. For example cont.force should have as many rows as many you have in the original spectral data x and as many columns as many pure components you want to resolve. Feel all values of this matrix with NA and the values you want to force with real numbers. For example if you know that in the first measurement concentration of 2 and 3 components was zero, set the corresponding values of cont.force to zero. See also the last case in the examples section.

Value

Returns an object of mcrpure class with the following fields:

resspec

matrix with resolved spectra.

rescont

matrix with resolved contributions.

cont.constraints

list with contribution constraints provided by user.

spec.constraints

list with spectra constraints provided by user.

expvar

vector with explained variance for each component (in percent).

cumexpvar

vector with cumulative explained variance for each component (in percent).

ncomp

number of resolved components

max.niter

maximum number of iterations

info

information about the model, provided by user when build the model.

More details and examples can be found in the Bookdown tutorial.

Author(s)

Sergey Kucheryavskiy (svkucheryavski@gmail.com)

References

1. J. Jaumot, R. Gargallo, A. de Juan, and R. Tauler, "A graphical user-friendly interface for MCR-ALS: a new tool for multivariate curve resolution in MATLAB", Chemometrics and Intelligent #' Laboratory Systems 76, 101-110 (2005).

See Also

Methods for mcrals objects:

summary.mcrals shows some statistics for the case.
predict.mcrals computes contributions by projection of new spectra to the resolved ones.

Plotting methods for mcrals objects:

plotSpectra.mcr shows plot with resolved spectra.
plotContributions.mcr shows plot with resolved contributions.
plotVariance.mcr shows plot with explained variance.
plotCumVariance.mcr shows plot with cumulative explained variance.

Examples




library(mdatools)

# resolve mixture of carbonhydrates Raman spectra

data(carbs)

# define constraints for contributions
cc <- list(
   constraint("nonneg")
)

# define constraints for spectra
cs <- list(
   constraint("nonneg"),
   constraint("norm", params = list(type = "area"))
)

# because by default initial approximation is made by using random numbers
# we need to seed the generator in order to get reproducable results
set.seed(6)

# run ALS
m <- mcrals(carbs$D, ncomp = 3, cont.constraints = cc, spec.constraints = cs)
summary(m)

# plot cumulative and individual explained variance
par(mfrow = c(1, 2))
plotVariance(m)
plotCumVariance(m)

# plot resolved spectra (all of them or individually)
par(mfrow = c(2, 1))
plotSpectra(m)
plotSpectra(m, comp = 2:3)

# plot resolved contributions (all of them or individually)
par(mfrow = c(2, 1))
plotContributions(m)
plotContributions(m, comp = 2:3)

# of course you can do this manually as well, e.g. show original
# and resolved spectra
par(mfrow = c(1, 1))
mdaplotg(
   list(
      "original" = prep.norm(carbs$D, "area"),
      "resolved" = prep.norm(mda.subset(mda.t(m$resspec), 1), "area")
   ), col = c("gray", "red"), type = "l"
)

# in case if you have reference spectra of components you can compare them with
# the resolved ones:
par(mfrow = c(3, 1))
for (i in 1:3) {
   mdaplotg(
      list(
         "pure" = prep.norm(mda.subset(mda.t(carbs$S), 1), "area"),
         "resolved" = prep.norm(mda.subset(mda.t(m$resspec), 1), "area")
      ), col = c("gray", "red"), type = "l", lwd = c(3, 1)
   )
}

# This example shows how to force some of the contribution values
# First of all we combine the matrix with mixtures and the pure spectra, so the pure
# spectra are on top of the combined matrix
Dplus <- mda.rbind(mda.t(carbs$S), carbs$D)

# since we know that concentration of C2 and C3 is zero in the first row (it is a pure
# spectrum of first component), we can force them to be zero in the optimization procedure.
# Similarly we can do this for second and third rows.

cont.forced <- matrix(NA, nrow(Dplus), 3)
cont.forced[1, ] <- c(NA, 0, 0)
cont.forced[2, ] <- c(0, NA, 0)
cont.forced[3, ] <- c(0, 0, NA)

m <- mcrals(Dplus, 3, cont.forced = cont.forced, cont.constraints = cc, spec.constraints = cs)
plot(m)

# See bookdown tutorial for more details.




svkucheryavski/mdatools documentation built on Aug. 25, 2023, 12:27 p.m.