Signal smoothing algorithm (Matlab's moving average)

Question

I have written a simple code that performs a 3-point moving average smoothing algorithm. It is meant to follow the same basic algorithm as Matlab's smooth(...) function as described here.

However, the result of my code is very different from that of Matlab's. Matlab's 3-point filter appears to perform a much more aggressive smoothing.

Here is a comparison of a noisy data smoothed using my code (red) and Matlab's function (blue):

Here is my code written in the form of a function:

function [NewSignal] = smoothing(signal)
NewSignal = signal;
for i = 2 : length(signal)-1
    NewSignal(i,:) = (NewSignal(i,:)+NewSignal(i-1,:)+NewSignal(i+1,:))./3;
end
end

Matlab's function is used as follows:

signal = smooth(time, signal, 3, 'moving');

As far as I understand Matlab's function works the same way; it averages 3 adjacent bins to a single bin. So I expected both algorithms to produce the same results.

So, what is the reason for the discrepancy? And how could I tweak my code to produce the same results?

Edit:

My sample data can be found here. It can be accessed using:

M = csvread('DS0009.csv');
time = M(:,1);
signal = M(:,2);

Here is the new result (red plot) using rinkert's correction:

Answer 1

To answer your question in the comments: the Function filter and smooth perform arithmetically the same (in the case that they are applied for moving average). however, there are the special cases at the beginning and endpoints which are handled differently. This is also stated in the documentation of smooth "Because of the way smooth handles endpoints, the result differs from the result returned by the filter function."

Here you see it in an example:

%generate randonm data
signal=rand(1,50);

%plot data
plot(signal,'LineWidth',2)
hold on
%plot filtered data
plot(filter(ones(3,1)/3,1,signal),'r-','LineWidth',2)
%plot smoothed data
plot( smooth(signal,3,'moving'),'m--','LineWidth',2)
%plot smoothed and delayed
plot([zeros(1,1); smooth(signal,3,'moving')],'k--','LineWidth',2)
hold off
legend({'Data','Filter','Smooth','Smooth-Delay'})

As you can see the filtered data (in red) is just a delayed version of the smoothed data (in magenta). Additionally, they differ in the beginning. Delaying the smoothed data results in an identical waveform as the filtered data (besides the beginning). As rinkert pointed out, your approach overwrites the data points which you are accessing in the next step. This is a different issue.

In the next example you will see that rinkerts implementation (smooth-rinkert) is identical to matlabs smooth, and that your approach differs from both due to overwriting the values:

So it is your function which low passes the input stronger. (as pointed out by Cris)

Answer 2

One reason for the difference could be that you are partially using your smoothed signal during smoothing. In your loop, you store the smoothed value in NewSignal(i,:), and for the next sample to smooth this value will be called by NewSignal(i-1,:).

Let NewSignal be determined by the original signal only:

function [NewSignal] = smoothing(signal)
    NewSignal = signal;
    for i = 2 : length(signal)-1
        NewSignal(i,:) = (signal(i,:)+signal(i-1,:)+signal(i+1,:))./3;
    end
end

---

Update: To show that the function above in fact does the same as Matlab's smooth function, let's consider this MVCE:

t = (0:0.01:10).';      % time vector 
y = sin(t) + 0.5*randn(size(t));

y_smooth1 = smooth(t,y,3,'moving');
y_smooth2 = smoothing(y);

difference_methods = abs(y_smooth1-y_smooth2);

So creating a sine wave, add some noise, and determine the absolute difference between the two methods. If you take the sum of all the differences you will see that this adds up to something like 7.5137e-14, which cannot explain the differences you see.

Plotting the smooth signal (blue original, red smoothed):

figure(1); clf; hold on
plot(t,y)
plot(t,y_smooth2)

And then plotting the difference between the two methods:

figure(2); clf; hold on;
plot(t,y_smooth1-y_smooth2)

As you can see, the difference is of the order 1e-16, so influenced by the Floating-point relative accuracy (see eps).