The idea of smoothing data by averaging neighboring values is older than the name we give it today. In 1901, the English statistician R. H. Hooker published a study correlating marriage rates with trade prosperity, where he computed what he called “instantaneous averages” — replacing each year’s figure with the mean of itself and its neighbors to reveal long-term trends hidden under year-to-year fluctuations. Eight years later, G. U. Yule described Hooker’s technique as “moving-averages” in the Journal of the Royal Statistical Society (1909), and the term entered general use after W. I. King adopted it in his textbook Elements of Statistical Method (1912).
But the practice predates the name by at least a century. Actuaries smoothing mortality tables in the 18th and 19th centuries routinely averaged adjacent age-group death rates to iron out irregularities caused by small sample sizes — a running mean in all but name. Astronomers did the same when reducing series of observations, averaging consecutive readings to suppress random measurement error long before anyone formalized the procedure.
What makes the moving average so enduring is its simplicity. When digital signal processing emerged in the mid-20th century and engineers needed filters they could implement on the earliest computers, the moving average — or “boxcar filter,” as it became known in that context — was the obvious first choice. It requires only addition and division, no trigonometric functions, no matrix algebra, no iterative fitting. That simplicity is why it remains the first smoothing method most scientists learn, and a natural starting point before moving to more sophisticated techniques like Savitzky-Golay or Whittaker smoothing.
Real spectroscopic data is noisy. When you measure a fluorescence emission spectrum, an IR absorption scan, or a Raman spectrum, the signal doesn’t come out as a smooth curve. Instead, it wiggles and jumps (sometimes going higher than expected, sometimes lower) even when measuring the exact same sample under identical conditions.
This isn’t a flaw in your instrument. It’s fundamental physics. Spectroscopic measurements involve counting photons (or electrons, or other discrete events), and these arrivals are inherently random. Even the best instruments produce noisy data. (We discussed the nature of measurement uncertainty and where noise comes from in detail in Least Squares: your first step into chemometrics.)
Your brain can easily see past the noise to identify the underlying smooth pattern: a peak here, a baseline there. But computers just see numbers: 0.421, 0.398, 0.447, 0.389… And this creates real problems for analysis. Random fluctuations can create false peaks or hide real ones, making accurate peak identification difficult. Noise adds uncertainty to measurements of peak heights and areas, which are critical for quantitative work. When building calibration models, regression algorithms may fit noise instead of signal, reducing prediction accuracy. Noise also obscures real differences between samples during spectral comparison and can completely hide subtle spectral features.
To solve this, we need smoothing - methods that remove random fluctuations while preserving the underlying chemical information. This makes downstream analysis more reliable and reproducible.
The simplest smoothing method is based on one elegant observation:
Noise is random. Signal, sometimes, is continuous.
Many spectroscopic signals arise from physical phenomena (molecular vibrations, electronic transitions, fluorescence emission) that use to change smoothly with wavelength or frequency. A molecule’s absorption doesn’t jump wildly from one wavelength to the next - it follows continuous, smooth shapes. If the intensity at 500 nm is 0.40, the intensity at 501 nm will be something close to 0.40, not 0.95 or 0.02. This is continuity: neighboring points in a spectrum carry similar information.
Noise, on the other hand, has no such behavior. It’s random at every point, independent of its neighbors. If noise pushes the value at 500 nm upward, it’s just as likely to push the value at 501 nm downward.
This difference is what we exploit. Imagine the true signal at wavelengths 499, 500, and 501 nm is 0.39, 0.40, and 0.41 - nearly the same, because the signal is continuous. But noise makes your measurements come out as 0.35, 0.46, 0.38 - scattered around those true values. If you average these three measured points, you get 0.397, which is very close to the true value at 500 nm (0.40). The noise, being random (some positive, some negative), tends to cancel itself out. The signal, being nearly the same across neighbors, survives the averaging almost unchanged.
Not all signals are continuous, though. Some spectroscopic techniques (like mass spectrometry or atomic emission) produce sharp, discrete features where neighboring points are not similar. Smoothing methods that assume continuity can distort or destroy these features, so keep your signal’s nature in mind before applying any smoothing.
Moving average smoothing (also called running mean or boxcar filter) applies this principle by replacing each point with the average of itself and its neighbors:
Where:
That’s it! Since noise pushes values randomly higher and lower, averaging nearby points tends to cancel out the randomness, leaving the smooth underlying trend.
Before diving into the math, let’s see how moving average actually works:
Try this:
The purple star shows the smoothed value being calculated, while the window (shaded blue) highlights which neighbors are being averaged.
Now that you’ve seen the process, let’s break down the algorithm step by step.
For a fluorescence spectrum with intensity values at different wavelengths:
Choose a window size
This determines how many neighbors each point consults for averaging.
For each point i, gather neighboring values
Collect values from to (the point itself plus neighbors on both sides).
Calculate the average
Sum all values in the window and divide by the window size.
Replace the original value
The smoothed value at point i becomes this average.
Suppose you have these intensity values:
Original: [10.2, 15.1, 12.3, 19.8, 18.2, 15.9, 14.1, 18.7, 17.2, 15.3]With a window size of 3 (each point averaged with its two immediate neighbors):
Smoothed: [12.5, 15.7, 16.8, 18.0, 16.1, 16.1, 16.2, 16.7, ...]Notice how the wild jump from 12.3 to 19.8 becomes a gentler transition. The noise is being averaged away.
The window size determines how aggressively you smooth. Larger windows give you smoother data, but they also blur your features. There’s always a tradeoff:
Smaller windows (relative to your peak width):
Moderate windows (balance point):
Larger windows (approaching peak width):
So how do you choose? The key is to think about your data characteristics:
For a spectrum with 100 points and peaks spanning ~10-20 points:
For a high-resolution spectrum with 1000 points and peaks spanning ~100 points:
For ultra-high-resolution data with 50,000 points:
Window sizes are typically odd (3, 5, 7, 9, 11) for symmetry. With an odd window, you get the same number of neighbors on both sides of the central point:
With an even window (like 4 or 8), you’d have to offset the window left or right, which subtly shifts your data. Odd windows keep everything centered and symmetric.
At the start and end of your spectrum, there aren’t enough neighbors on one side to fill the window. Different approaches handle this differently:
Truncation (most common):
Edge padding:
convolve with mode='same' does thisReflection:
Reduced window:
For most spectroscopic applications, the interesting features are in the middle of the spectrum, so edge effects aren’t critical. But if you need accurate edge values (e.g., for integration), pay attention to how your software handles edges.
Everything we’ve discussed so far assumes you have your complete dataset available. But what if you’re collecting data in real time, point by point?
In offline spectroscopy, you collect the entire spectrum first, then apply preprocessing. The centered moving average we’ve described works perfectly: you can use points both before and after the current position.
But in real-time applications, you can’t use future data that hasn’t been measured yet. This is critical in:
In these scenarios, you need a causal (or one-sided) moving average that only uses the current point and past measurements:
Instead of centering the window around point i, the window extends backward in time from the current point.
Notice the difference from the centered version above: the window (red dots) only appears to the left of the current point. There are no “future” measurements to use. Try increasing the window size and watch how the smoothed signal lags further behind the true signal - this is the fundamental cost of real-time smoothing.
Causal moving average introduces a time lag (or phase shift). The smoothed signal is delayed compared to the true signal by approximately time points.
Example: With w = 9 in a PAT application measuring NIR spectra every 30 seconds:
Practical implications:
Choosing window size for real-time:
The optimal window size depends on your process dynamics, noise level, and acceptable lag time. Consider:
Balance noise reduction against the acceptable delay for your specific application.
def causal_moving_average(y, window_size): """ One-sided moving average for real-time applications. Only uses current and past data points. """ smoothed = [] for i in range(len(y)): # Window extends backward from current point start = max(0, i - window_size + 1) window_values = y[start:i+1] smoothed.append(sum(window_values) / len(window_values)) return smoothed
# Real-time PAT example: NIR spectra acquired every 30 secondsspectra = acquire_nir_spectrum() # Latest measurementsmoothed_spectrum = causal_moving_average(spectra, window_size=5)You can, but it’s not recommended. Each pass blurs features more. If one pass with w = 7 isn’t enough, use a single pass with w = 11 instead - you’ll get better results.
Yes, absolutely! Consistency is critical. If you’re building a calibration model or comparing spectra, use the same window size for all data (calibration, validation, and test sets). Inconsistent smoothing creates inconsistent features.
Check the residuals:
Good smoothing: Residuals look like random noise (no patterns, symmetric around zero)
Over-smoothing: Residuals show structured patterns (like half of a peak) - you’re removing signal. Use a smaller window.
Averaging a sharp peak with its neighboring baseline values pulls the peak value down. The peak spreads out to become shorter and broader, though the area stays roughly the same.
This is why moving average isn’t ideal for quantitative work on sharp peaks. Savitzky-Golay smoothing preserves peak heights much better.
Here’s how to implement moving average in Python, MATLAB, and R:
import numpy as npimport matplotlib.pyplot as plt
def moving_average(y, window_size): """ Apply moving average smoothing to a 1D signal.
Parameters: ----------- y : array-like Input signal to be smoothed (e.g., spectrum intensities) window_size : int Number of points in the averaging window (use odd numbers!)
Returns: -------- y_smooth : array Smoothed signal (same length as input) """ # Create the averaging kernel (array of 1/w values) window = np.ones(window_size) / window_size
# Apply convolution (NumPy handles edge padding automatically) y_smooth = np.convolve(y, window, mode='same')
return y_smooth
# Example: noisy fluorescence emission spectrumnp.random.seed(42)
# Wavelength range typical for fluorescence (nm)wavelength = np.linspace(400, 600, 100)
# Simulate a realistic fluorescence spectrum with vibronic structuretrue_signal = (0.85 * np.exp(-((wavelength - 480)**2) / 400) + 0.35 * np.exp(-((wavelength - 510)**2) / 300) + 0.15 * np.exp(-((wavelength - 545)**2) / 250))
# Add realistic shot noise (proportional to sqrt of signal)noise = np.random.normal(0, 0.08 * np.sqrt(true_signal))noisy_spectrum = true_signal + noise
# Try different window sizessmooth_w5 = moving_average(noisy_spectrum, window_size=5)smooth_w9 = moving_average(noisy_spectrum, window_size=9)smooth_w15 = moving_average(noisy_spectrum, window_size=15)
# Plot and compareplt.figure(figsize=(12, 6))plt.plot(wavelength, noisy_spectrum, 'o-', alpha=0.4, label='Noisy data', color='gray', markersize=3)plt.plot(wavelength, smooth_w5, linewidth=2.5, label='Window = 5')plt.plot(wavelength, smooth_w9, linewidth=2.5, label='Window = 9')plt.plot(wavelength, smooth_w15, linewidth=2.5, label='Window = 15')plt.xlabel('Wavelength (nm)', fontsize=12)plt.ylabel('Fluorescence Intensity (a.u.)', fontsize=12)plt.title('Moving Average Smoothing: Window Size Comparison', fontsize=14)plt.legend(fontsize=10)plt.grid(True, alpha=0.3)plt.tight_layout()plt.show()
# Check residualsresiduals = noisy_spectrum - smooth_w9plt.figure(figsize=(12, 4))plt.plot(wavelength, residuals, 'o-', markersize=3)plt.axhline(y=0, color='red', linestyle='--', alpha=0.5)plt.xlabel('Wavelength (nm)')plt.ylabel('Residual (Original - Smoothed)')plt.title('Residuals: Should look like random noise with no patterns')plt.grid(True, alpha=0.3)plt.tight_layout()plt.show()function y_smooth = moving_average(y, window_size) % Apply moving average smoothing to a 1D signal % % Parameters: % y - Input signal to smooth (column or row vector) % window_size - Number of points in averaging window (use odd!) % % Returns: % y_smooth - Smoothed signal (same size as input)
% Create the averaging kernel window = ones(1, window_size) / window_size;
% Apply convolution (MATLAB handles edge padding) y_smooth = conv(y, window, 'same');end
% Example: fluorescence emission spectrumrng(42);
% Typical fluorescence wavelength range (nm)wavelength = linspace(400, 600, 100);
% Simulate realistic fluorescence with vibronic bandstrue_signal = 0.85 * exp(-((wavelength - 480).^2) / 400) + ... 0.35 * exp(-((wavelength - 510).^2) / 300) + ... 0.15 * exp(-((wavelength - 545).^2) / 250);
% Add realistic shot noisenoise = 0.08 * sqrt(true_signal) .* randn(size(wavelength));noisy_spectrum = true_signal + noise;
% Try different window sizessmooth_w5 = moving_average(noisy_spectrum, 5);smooth_w9 = moving_average(noisy_spectrum, 9);smooth_w15 = moving_average(noisy_spectrum, 15);
% Visualize the comparisonfigure('Position', [100 100 1000 500]);plot(wavelength, noisy_spectrum, 'o-', 'Color', [0.7 0.7 0.7], ... 'MarkerSize', 3, 'DisplayName', 'Noisy data');hold on;plot(wavelength, smooth_w5, 'LineWidth', 2.5, 'DisplayName', 'Window = 5');plot(wavelength, smooth_w9, 'LineWidth', 2.5, 'DisplayName', 'Window = 9');plot(wavelength, smooth_w15, 'LineWidth', 2.5, 'DisplayName', 'Window = 15');xlabel('Wavelength (nm)', 'FontSize', 12);ylabel('Fluorescence Intensity (a.u.)', 'FontSize', 12);title('Moving Average: Window Size Comparison', 'FontSize', 14);legend('Location', 'northeast', 'FontSize', 10);grid on;
% Check residualsfigure('Position', [100 100 1000 300]);residuals = noisy_spectrum - smooth_w9;plot(wavelength, residuals, 'o-', 'MarkerSize', 3);yline(0, '--r', 'LineWidth', 1.5, 'Alpha', 0.5);xlabel('Wavelength (nm)');ylabel('Residual');title('Residuals: Should be random noise, no structure');grid on;moving_average <- function(y, window_size) { # Apply moving average smoothing to a 1D signal # # Parameters: # y - Input signal to smooth (numeric vector) # window_size - Number of points in window (use odd numbers!) # # Returns: # y_smooth - Smoothed signal
# Create averaging kernel window <- rep(1 / window_size, window_size)
# Apply convolution with centered window y_smooth <- stats::filter(y, window, sides = 2)
# Handle edge NAs by using original edge values # (R's filter leaves NAs at edges by default) y_smooth[is.na(y_smooth)] <- y[is.na(y_smooth)]
return(as.vector(y_smooth))}
# Example: fluorescence emission spectrumset.seed(42)
# Wavelength range for fluorescence (nm)wavelength <- seq(400, 600, length.out = 100)
# Simulate realistic fluorescence with vibronic structuretrue_signal <- 0.85 * exp(-((wavelength - 480)^2) / 400) + 0.35 * exp(-((wavelength - 510)^2) / 300) + 0.15 * exp(-((wavelength - 545)^2) / 250)
# Add realistic shot noisenoise <- 0.08 * sqrt(true_signal) * rnorm(length(wavelength))noisy_spectrum <- true_signal + noise
# Compare different window sizessmooth_w5 <- moving_average(noisy_spectrum, 5)smooth_w9 <- moving_average(noisy_spectrum, 9)smooth_w15 <- moving_average(noisy_spectrum, 15)
# Visualizepar(mfrow = c(1, 1), mar = c(4, 4, 3, 1))plot(wavelength, noisy_spectrum, type = 'o', col = 'gray60', pch = 16, cex = 0.6, xlab = 'Wavelength (nm)', ylab = 'Fluorescence Intensity (a.u.)', main = 'Moving Average: Window Size Comparison')lines(wavelength, smooth_w5, col = '#2563eb', lwd = 2.5)lines(wavelength, smooth_w9, col = '#059669', lwd = 2.5)lines(wavelength, smooth_w15, col = '#dc2626', lwd = 2.5)legend('topright', legend = c('Noisy', 'w=5', 'w=9', 'w=15'), col = c('gray60', '#2563eb', '#059669', '#dc2626'), lwd = c(1, 2.5, 2.5, 2.5), pch = c(16, NA, NA, NA), pt.cex = 0.6, cex = 0.9)grid()
# Check residualspar(mar = c(4, 4, 3, 1))residuals <- noisy_spectrum - smooth_w9plot(wavelength, residuals, type = 'o', pch = 16, cex = 0.6, xlab = 'Wavelength (nm)', ylab = 'Residual', main = 'Residuals: No patterns = good smoothing')abline(h = 0, col = 'red', lty = 2, lwd = 2)grid()✅ Works well for
Quick exploratory analysis Fast and simple, perfect for initial data inspection
Broad, gentle features Wide, smooth peaks (like fluorescence emission) are handled well
Teaching and communication Intuitive concept (“averaging with neighbors”) that’s easy to explain
Baseline estimation Use a very large window (50+ points) to capture slowly-varying baseline drift
Speed-critical applications Extremely fast computation, ideal for batch processing thousands of spectra
❌ Use something better for
Sharp Raman peaks Moving average significantly reduces peak heights and increases widths. Use Savitzky-Golay instead.
Quantitative analysis Peak broadening and height reduction introduce bias in calibration curves
High-resolution work Resolution matters (vibronic structure, closely-spaced peaks) and moving average’s blurring is unacceptable
Publication-quality figures Moving average can leave visible artifacts (slight stair-stepping effects)
Advantages
✅ Extremely simple Easy to understand and implement
✅ Very fast Minimal computational cost
✅ No parameters beyond window size Just choose w and you’re done
✅ Predictable behavior Always smooths, never does anything surprising
Limitations
❌ Reduces peak heights Averaging pulls sharp peaks down toward their neighbors
❌ Broadens features Peaks get wider, closely-spaced features may merge
❌ Equal weighting All neighbors are treated equally, even those far from the center point
❌ Not optimal for most applications Other methods (Savitzky-Golay, Whittaker) preserve features better
Window size selection:
Quality control:
Computational efficiency:
convolve is highly optimizedMoving average teaches the fundamental concepts of smoothing: noise reduction vs. feature preservation, window size tradeoffs, and edge effects. These principles apply to all smoothing methods.
For most serious spectroscopic work, consider these more advanced methods:
Savitzky-Golay smoothing - Fits local polynomials instead of averaging. Much better peak preservation. The most common choice for Raman and IR spectroscopy.
Whittaker smoothing - Optimization-based approach that balances smoothness against data fidelity. Excellent general-purpose method.
Gaussian smoothing - Weighted averaging where close neighbors matter more than distant ones. Smoother than moving average with similar simplicity.
Moving average is the simplest smoothing method, making it perfect for quick exploratory analysis and teaching. But for quantitative work or publication-quality results, the more sophisticated methods above are usually worth the extra complexity.