Multiplicative Scatter Correction (MSC)

In the late 1970s and early 1980s, near-infrared (NIR) spectroscopy was rapidly gaining ground as a fast, non-destructive way to analyse agricultural and food products. You could point an NIR instrument at ground grain, flour or minced meat and, in seconds, obtain a spectrum related to its protein, moisture or fat content. But there was a stubborn problem: the spectra were dominated by light scattering. When NIR light enters a powder or a granular sample, much of it bounces between particles before it reaches the detector. How much it bounces depends on physical properties (particle size, surface texture, how tightly the sample is packed) that have nothing to do with chemistry. Two samples with the same composition but a different grind could give visibly different spectra.

Harald Martens and his colleagues in Norway attacked this problem with a strikingly simple idea. If scattering mostly shifts a spectrum up or down and stretches it, then each spectrum should look like a shifted and stretched copy of some "typical" spectrum. Estimate the shift and the stretch with a straight-line regression, undo them, and what remains should be mostly chemistry. The method was first presented by Martens, Jensen and Geladi in 1983 [1] and became widely known through the 1985 paper by Geladi, MacDougall and Martens on NIR reflectance spectra of meat [2]. They called it Multiplicative Scatter Correction (MSC); you will also see it called multiplicative signal correction, with the same acronym, the reading used by Martens and Stark when they extended the method in 1991 [10].

MSC was then given a central place in Martens and Næs' book Multivariate Calibration (1989) [3], one of the foundational textbooks of chemometrics. Its simplicity (one small regression per spectrum) and its clear physical meaning made it one of those rare methods that practitioners adopted almost universally. Decades later, it is still one of the most common preprocessing steps for NIR spectra of solids [4].

In diffuse reflectance, light that comes back from the sample has travelled through it along many different paths. Fine particles scatter more and send light back after a short path; coarse particles let it travel further. The chemistry is the same, but the effective path length is not, and the spectrum changes accordingly.

In practice, scattering produces two main effects:

• an additive effect: the whole spectrum moves up or down (a baseline offset),
• a multiplicative effect: the whole spectrum is stretched or compressed, so tall bands change more than small ones.

The problem is that these physical effects are often larger than the chemical differences we want to measure. A PCA of raw NIR spectra of powders very often spends its first component describing particle size. A PLS model can still cope, but it needs extra components just to model scatter, and it becomes harder to interpret and less robust.

MSC tries to remove both effects before modelling.

MSC assumes that every measured spectrum ${\mathbf{x}}_i$ is approximately a shifted and stretched copy of a reference spectrum ${\mathbf{x}}_{\text{ref}}$:

Here:

• $a_i$ is the additive scatter term (the offset),
• $b_i$ is the multiplicative scatter term (the stretch),
• $1$ is a vector of ones, one per wavelength,
• ${\mathbf{e}}_i$ is whatever the line cannot explain. Ideally, that is the chemical information plus noise.

This is just a simple linear regression, but with an unusual twist: the "x" of the regression is the reference spectrum and the "y" is the sample spectrum, and the data points are the wavelengths. With a spectrum of a few hundred wavelengths, we fit a straight line through hundreds of points to estimate only two numbers.

Once $a_i$ and $b_i$ are known, we undo the distortion:

Subtracting $a_i$ removes the offset, and dividing by $b_i$ removes the stretch. After the correction, every spectrum has been brought onto the same "scatter level" as the reference.

In practice $b_i$ should be close to 1. A value near zero, or a negative one, means that the spectrum does not resemble the reference at all (a failed measurement, a different kind of sample), and it should be inspected rather than corrected: dividing by a tiny $b_i$ would blow the spectrum up.

The demo below shows what that regression looks like for one spectrum. Every dot is one wavelength. If the sample differed from the reference only by scatter, all the dots would fall exactly on a straight line. They nearly do: the line is the scatter. The few dots that drift away from the line, more visible for the samples at the extremes of the colour scale, are where the chemistry of this sample differs from the average.

Pick different samples and compare the fitted $a$ and $b$ with the values used to simulate the data. They are close, but not identical, because the chemical differences also nudge the regression a little. We come back to that in the limitations. Keep in mind that the fitted values are always relative to the reference: the mean spectrum itself carries the average scatter of the set, so in general the fitted $b_i$ estimates the true stretch divided by the average stretch. Here the average simulated stretch is almost exactly 1, so the two look alike.

The third plot shows what the straight line leaves behind at each wavelength. The residual peaks at the analyte bands: that is the chemistry MSC is supposed to keep.

Step 1: choose a reference spectrum. Almost always, this is the mean spectrum of the calibration (training) set:

Step 2: regress each spectrum on the reference. For every sample, fit ${\mathbf{x}}_i\approx a_i+b_i{\mathbf{x}}_{\text{ref}}$ by ordinary least squares. In matrix form, with $\mathbf{Z}=[1,{\mathbf{x}}_{\text{ref}}]$, this is exactly the problem solved by the normal equations:

Step 3: correct.

Step 4: store the reference. New samples (test set, future measurements) are corrected against the same stored reference. They are never used to recompute it.

That is all. MSC is one of the cheapest preprocessing methods there is: two numbers per spectrum.

Now let us look at a whole set of spectra. The eight simulated samples differ in their content of one analyte, which has a band at 1730 nm (dashed line). Their colours go from blue (low content) to amber (high content). On top of that chemistry, each spectrum has received a random offset and stretch.

Things to try:

1. Set the scatter to zero. Raw and corrected spectra look the same, and the R² of the band is already high.
2. Raise the scatter. The raw spectra fan out, the colour order at 1730 nm is scrambled, and the R² of the raw band collapses. After MSC, the spectra lie on top of each other except at the chemical bands, and the R² comes back.
3. Now look at analyte signal kept. It compares the spread of the corrected band heights with the true, scatter-free ones, and it stays around 84 % whatever the scatter, even at zero. Part of the chemical difference looks like a "stretch" to the regression, so MSC removes about 15 % of the analyte band together with the scatter. The R² stays high because the loss is proportional (every sample loses the same fraction), so a calibration model can still work, but band heights after MSC are biased towards the mean and should not be read as absolute amounts.
4. Move the chemical variation slider. The fraction kept barely moves in this simulation, because only the analyte changes from sample to sample. The leak is worst when the varying chemistry looks like the reference spectrum itself. Set the chemical variation to zero and the R² values show "n/a": with identical chemistry there is nothing to measure.

MSC is a good candidate when:

• you measure diffuse reflectance (or transmittance) of solids: powders, grains, tablets, granulates, minced or ground materials,
• the spectra show a fan shape: offset and spread that grow together across the spectrum,
• the first principal component of the raw spectra mostly describes particle size or packing rather than chemistry,
• the samples are chemically similar enough that a common reference makes sense.

MSC is a poor choice when:

• the samples are chemically very different (the "same spectrum, different scatter" assumption breaks down),
• scatter is negligible, for example clear liquids in transmission,
• the scatter effect is strongly wavelength dependent, so that a single offset and a single stretch cannot describe it. That is the situation EMSC was designed for.

import numpy as np
import matplotlib.pyplot as plt

def msc(spectra, reference=None, mask=None):
    """
    Multiplicative Scatter Correction.

    spectra   : array (n_samples, n_wavelengths)
    reference : array (n_wavelengths,), optional. Defaults to the mean spectrum.
    mask      : boolean array (n_wavelengths,), optional. If given, a and b are
                estimated only from these wavelengths (regions with little
                chemical variation); the correction is applied to all of them.
    Returns (corrected spectra, reference). Keep the reference to correct new samples.
    """
    spectra = np.asarray(spectra, dtype=float)
    if reference is None:
        reference = spectra.mean(axis=0)
    if mask is None:
        mask = np.ones(reference.size, dtype=bool)

    Z = np.column_stack([np.ones_like(reference), reference])   # [1, x_ref]
    corrected = np.empty_like(spectra)
    for i, x in enumerate(spectra):
        a, b = np.linalg.lstsq(Z[mask], x[mask], rcond=None)[0]   # x ≈ a + b * x_ref
        corrected[i] = (x - a) / b
    return corrected, reference

# Simulated NIR spectra: same chemistry, different scatter
rng = np.random.default_rng(42)
wl = np.linspace(1100, 2500, 300)
band = 0.5 * np.exp(-((wl - 1500) ** 2) / 5000) + 0.3 * np.exp(-((wl - 2000) ** 2) / 3000)
chem = band + 0.2                    # identical for every sample

raw = []
for _ in range(20):
    a = rng.uniform(-0.1, 0.1)       # additive scatter
    b = rng.uniform(0.8, 1.2)        # multiplicative scatter
    raw.append(a + b * chem + rng.normal(0, 0.002, wl.size))
raw = np.array(raw)

train, test = raw[:15], raw[15:]
train_msc, ref = msc(train)          # reference from the training set only
test_msc, _ = msc(test, reference=ref)

print('spread between raw spectra      :', np.ptp(train, axis=0).max())
print('spread between corrected spectra:', np.ptp(train_msc, axis=0).max())   # only noise left

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4), sharey=False)
ax1.plot(wl, train.T, lw=1); ax1.set_title('Raw')
ax2.plot(wl, train_msc.T, lw=1); ax2.set_title('After MSC')
for ax in (ax1, ax2):
    ax.set_xlabel('Wavelength (nm)')
plt.tight_layout()
plt.show()

% Simulated NIR spectra: same chemistry, different scatter
rng(42);
wl = linspace(1100, 2500, 300);
band = 0.5 * exp(-((wl - 1500).^2) / 5000) + 0.3 * exp(-((wl - 2000).^2) / 3000);
chem = band + 0.2;                   % identical for every sample

raw = zeros(20, numel(wl));
for i = 1:20
    a = -0.1 + 0.2 * rand;           % additive scatter
    b = 0.8 + 0.4 * rand;            % multiplicative scatter
    raw(i, :) = a + b * chem + 0.002 * randn(size(wl));
end

train = raw(1:15, :);
test  = raw(16:20, :);
[trainMSC, ref] = msc(train);        % reference from the training set only
testMSC = msc(test, ref);

subplot(1, 2, 1); plot(wl, train'); title('Raw'); xlabel('Wavelength (nm)');
subplot(1, 2, 2); plot(wl, trainMSC'); title('After MSC'); xlabel('Wavelength (nm)');

function [corrected, reference] = msc(spectra, reference, mask)
    % Multiplicative Scatter Correction. Rows are samples.
    % mask (optional): logical vector; a and b are estimated only from these
    % wavelengths, and the correction is applied to all of them.
    if nargin < 2 || isempty(reference)
        reference = mean(spectra, 1);
    end
    if nargin < 3 || isempty(mask)
        mask = true(1, numel(reference));
    end
    Z = [ones(numel(reference), 1), reference(:)];
    corrected = zeros(size(spectra));
    for i = 1:size(spectra, 1)
        coef = Z(mask, :) \ spectra(i, mask)';    % [a; b]
        corrected(i, :) = (spectra(i, :) - coef(1)) / coef(2);
    end
end

msc <- function(spectra, reference = NULL, mask = NULL) {
  # Multiplicative Scatter Correction. Rows are samples.
  # mask (optional): logical vector; a and b are estimated only from these
  # wavelengths, and the correction is applied to all of them.
  if (is.null(reference)) reference <- colMeans(spectra)
  if (is.null(mask)) mask <- rep(TRUE, length(reference))
  Z <- cbind(1, reference)
  corrected <- spectra
  for (i in seq_len(nrow(spectra))) {
    coef <- lm.fit(Z[mask, , drop = FALSE], spectra[i, mask])$coefficients   # a, b
    corrected[i, ] <- (spectra[i, ] - coef[1]) / coef[2]
  }
  list(corrected = corrected, reference = reference)
}

# Simulated NIR spectra: same chemistry, different scatter
set.seed(42)
wl <- seq(1100, 2500, length.out = 300)
band <- 0.5 * exp(-((wl - 1500)^2) / 5000) + 0.3 * exp(-((wl - 2000)^2) / 3000)

chem <- band + 0.2                                   # identical for every sample

raw <- t(sapply(1:20, function(i) {
  a <- runif(1, -0.1, 0.1)
  b <- runif(1, 0.8, 1.2)
  a + b * chem + rnorm(length(wl), 0, 0.002)
}))

train <- raw[1:15, ]
test  <- raw[16:20, ]
fit <- msc(train)                                    # reference from training only
test_msc <- msc(test, reference = fit$reference)$corrected

par(mfrow = c(1, 2))
matplot(wl, t(train), type = "l", lty = 1, main = "Raw", xlab = "Wavelength (nm)")
matplot(wl, t(fit$corrected), type = "l", lty = 1, main = "After MSC", xlab = "Wavelength (nm)")

The simulated spectra share exactly the same chemistry, so after MSC the only differences left between them are noise. Beware of a tempting variation: if you let the height of the single band vary from sample to sample, MSC removes that variation too. With only one chemical shape, a change in its intensity is indistinguishable from a stretch.

Strengths

• Simple and fast. One straight-line regression per spectrum.
• Physically meaningful. The two coefficients have a clear interpretation as offset and stretch, and they can be plotted: a plot of $b_i$ can often be related to particle size or packing.
• Often improves models. Removing scatter typically gives PCA components and PLS models that are easier to interpret and need fewer components [4].

Limitations

• It assumes the samples are chemically similar. The regression cannot tell a chemical difference from a stretch. When the chemical bands are large, part of them leaks into $b_i$ and gets "corrected" away. Some practitioners therefore estimate $a_i$ and $b_i$ only from wavelength regions with little chemical variation, and then apply the correction to the whole spectrum [3]. That is what the mask argument of the code above does. A related idea, piecewise MSC, fits separate corrections in moving windows along the spectrum [8].
• The reference matters. The correction is relative to the reference spectrum. A reference dominated by a few unusual samples gives a poor correction.
• Only one offset and one stretch. Real scatter often depends on wavelength in a more complicated way (for example, it can produce sloping or curved baselines). MSC cannot remove that; EMSC can.

Choosing the reference.

• The mean spectrum of the training set is the standard choice.
• A median spectrum is more robust if there are a few outlying samples.
• Avoid using a single, randomly chosen sample as the reference: its own noise and peculiarities are copied into every corrected spectrum.

Signs that MSC will help.

• The spectra fan out, with offsets and spread that grow together.
• The first principal component of the raw spectra looks like a smooth, positive "average spectrum" shape, and its scores follow particle size or packing.
• A PLS model needs more components than you would expect from the chemistry.

Combining with other preprocessing.

• MSC followed by mean centering is a standard combination before PCA or PLS.
• MSC can be followed by a derivative if residual baseline effects remain.
• MSC and SNV correct the same effects, so you normally choose one of them, not both.

SNV (Standard Normal Variate) was proposed by Barnes, Dhanoa and Lister in 1989 [5]. It corrects each spectrum on its own: subtract the spectrum's own mean ${\bar{x}}_i$ (the average over its wavelengths) and divide by its own standard deviation $s_i$:

def snv(spectra):
    spectra = np.asarray(spectra, dtype=float)
    mean = spectra.mean(axis=1, keepdims=True)
    std = spectra.std(axis=1, keepdims=True)
    return (spectra - mean) / std

After SNV every spectrum has mean 0 and standard deviation 1, so an offset and a stretch disappear, just as with MSC. The difference is that SNV needs no reference, so there is nothing to store and nothing that can leak from a test set.

• Recomputing the reference for new data. The reference is part of the model. Compute it from the training set, store it, and correct test and future spectra against it.
• Letting test samples into the reference. If the mean spectrum includes the test set, information from the test set leaks into the training data and the validation looks better than it should.
• Thinking MSC cannot handle a single new spectrum. It can, as long as you have the stored reference. What you cannot do is build a meaningful reference from one spectrum.
• Stacking MSC and SNV. They do almost the same job; applying both rarely helps.
• Ignoring outliers in the reference. A few extreme spectra can distort the mean. Look at your data first, or use a median reference.

MSC treats scatter as a straight-line distortion: every spectrum is a shifted ($a_i$) and stretched ($b_i$) copy of a reference spectrum. A least squares fit across wavelengths estimates both numbers, and the correction $({\mathbf{x}}_i-a_i)/b_i$ undoes them.

It works best when samples are chemically similar and scatter is a single offset and stretch. Keep the reference from the training set, and if the scatter has a more complicated shape, move on to EMSC.