In this post I will talk about Minimum Mean Squared Error a little bit. This is an estimation method in which the mean of squared error between the actual and the estimated error is minimized.

Assume that we want to estimate some parameter \(\theta\) from the available set of observations, \(x = [x[0], x[1], ... , x[N-1]]\). For example suppose the observed data comes from the following process,

\(x[n] = A + w[n]\), w[n] is additive noise with a normal distribution. In this example we would be trying to estimate the parameter A(in real scenarios the actual value of A will not be available) from N observations. The most simple method for estimating A (this is the \(theta\) for this problem), as you have probably guessed it, is the sample mean of the observed data.

\[\hat{A} = frac{x[0] + x[1] + ... + x[N-1]}{N}\]

But then the question comes is this the best estimate of A we can have. Well, that depends on how are we judging a estimator function,

\(\bar{A} = f(x[0], x[1], ... , x[N-1])\).

Now to decide on whether a estimation is the best one we need a Objective function. There are quite a few of objective functions that are used in different estimation problems.

  • Mean Squared Error
  • Root Mean Squared Logarithmic Error
  • Hit or Miss
  • Sum of Absolute Difference and etc. etc.

Here we would be talking about the first one, MSE or Mean Squared Error.

\[MSE(\hat{A}) = E[(A - \hat{A})^2]\]

Thus a Minimum Mean Squared Error or MMSE estimate would be that $latex \hat{A}$, which minimizes the value of $latex MSE(\hat{A})$.

Now we verify the performance of MMSE estimator using the following MATLAB Code.

p = zeros(3, 1);
q = zeros(3, 1);
r = zeros(3, 1);
A = 3;
for j = 1:5
    N = 10;
    sg_n = sqrt(j);
    samp = 10000;

    arr_mse = zeros(samp, 1);
    arr_smean = zeros(samp, 1);

    mse = @(u, v) mean((u - v).^2);
    for i = 1:samp
        mu_a = 1;
        sg_a = sqrt(0.5);

        A = normrnd(mu_a, sg_a);
        x = normrnd(A, sg_n, N, 1);

        M = sum(x)/sg_n^2 + mu_a/sg_a^2;
        C = N/sg_n^2 + 1/sg_a^2;

        A_e = M/C;
        arr_mse(i) = mse(A, A_e);
        arr_smean(i) = mse(A, mean(x));

    disp(['Var: ' num2str(j)])
    mse_exp = mean(arr_mse)
    mse_smean = mean(arr_smean)

    % Theoretical
    mse_t = 1/(C)
    p(j) = mse_t;
    q(j) = mse_exp;
    r(j) = mse_smean;
hold on
plot(p, 'o')
legend('Theoretical', 'Post. Mean', 'Sample Mean');
hold off