1
$\begingroup$

I'm trying to estimate $X$ using MAP Estimator. If $Y = aX + N$, where $a$ is a constant, $X \sim \mathcal{N}(0,\sigma_{x}^{2})$ and $N \sim \mathcal{N}(0,\sigma_{n}^{2})$.

The likelihood function for $Y = aX + N$ can be evaluated as \begin{equation}\label{llf} f_{Y|X}(y|x)=\frac{1}{\sqrt{2\pi}\sigma^2}\exp\Bigg\{-\frac{1}{2}\frac{(y-\mu)^2}{\sigma^2}\Bigg\} \end{equation} where $\mu$ and $\sigma^2$ are mean and variance of $f_{Y|X}(y|x)$, which can be evaluated as \begin{align}\label{mean} \mu&=\mathbb{E}[Y\vert X] \notag\\ &= \mathbb{E}[aX + N\vert X] \notag\\ &= a\mathbb{E}[X\vert X]+\mathbb{E}[N\vert X] \notag\\ &=ax \end{align} and the variance can be evaluated as \begin{align}\label{var} \sigma^2 &= \mathbb{E}[(Y - \mu)^2\vert X]\notag\\ &= \mathbb{E}[Y^2\vert X] - \mu^2 \notag\\ &= \mathbb{E}[(aX + N)^2\vert X] - a^2x^2\notag\\ &= a^2\mathbb{E}[X^2\vert X]+\mathbb{E}[N^2\vert X] + 2a\mathbb{E}[X N\vert X] - a^2x^2 \notag\\ &= a^2x^2 + \sigma_{n}^2 + 2a\mathbb{E}[X\vert X]\mathbb{E}[N\vert X]-a^2x^2\notag\\ &= \sigma_{n}^2 \end{align} And a priori density function is \begin{equation}\label{ppdf} f_{X}(x)=\frac{1}{\sqrt{2\pi}\sigma_{x}^2}\exp\Bigg\{-\frac{x^2}{2\sigma_{x}^2}\Bigg\} \end{equation} Therefore the Bayesian log-likelihood function will be \begin{align} \mathcal{L}(X,Y) &= \ln f_{Y\vert X}(y\vert x)f_{X}(x)\notag\\ &=\ln\frac{1}{2\pi\sigma_{n}^2\sigma_{x}^2}\exp\Bigg\{-\frac{1}{2}\Bigg[\frac{(y-ax)^2}{\sigma_{n}^2} + \frac{x^2}{\sigma_{x}^2}\Bigg]\Bigg\}\notag\\ &= \ln\frac{1}{2\pi\sigma_{n}^2\sigma_{x}^2} -\frac{1}{2}\Bigg[\frac{(y-ax)^2}{\sigma_{n}^2} + \frac{x^2}{\sigma_{x}^2}\Bigg] \end{align} On differentiating above equation with respect to $x$, we get \begin{align} \frac{\partial}{\partial x}\mathcal{L}(X,Y)=\frac{\partial}{\partial x}\ln\frac{1}{2\pi\sigma_{n}^2\sigma_{x}^2} -\frac{1}{2}\frac{\partial}{\partial x}\Bigg[\frac{(y-ax)^2}{\sigma_{n}^2} + \frac{x^2}{\sigma_{x}^2}\Bigg] \end{align} And now equating the above equation for $x$ equal to zero. \begin{align} \frac{1}{2\sigma_{n}^2}\frac{\partial}{\partial x}(y-ax)^2 - \frac{1}{2\sigma_{x}^2}\frac{\partial}{\partial x}x^2 &= 0\notag\\ -\frac{(a^2x - ay)}{\sigma_{n}^2} - \frac{x}{\sigma_{x}^2}&=0\notag\\ x\Bigg(\frac{a^2\sigma_{x}^2 + \sigma_{n}^2}{\sigma_{x}^2\sigma_{n}^2}\Bigg) &= \frac{ay}{\sigma_{n}^2}\notag\\ x &= \frac{a\sigma_{x}^2y}{a^2\sigma_{x}^2 + \sigma_{n}^2} \end{align} Hence \begin{equation} \hat{X}_{MAP} = \frac{a\sigma_{x}^2y}{a^2\sigma_{x}^2 + \sigma_{n}^2} \end{equation}

Below is the MATLAB code which I'm using to simulate and validate the theoretical result obtained.

%% Release memory and clear screen
clear;clc;

%% Theoretical
varX = 3;                   % Variance of X
varN = 1;                   % Variance of N
a = 2;                      % Constant
L = 1e3;                    % Number of points
x = linspace(-10,10,L);

%% Generate Gaussianly Distributed Random Number
X = varX*randn;
N = varN*randn;
Y = a*X + N;
% MAP Estimator of X
Xhatsim = a*varX*Y/(a*a*varX + varN);
mu = a*X;
fx = normpdf(x,0,sqrt(varX));     % Density of X
fn = normpdf(x,0,sqrt(varN));     % Density of N
fyx = normpdf(x,mu,sqrt(varN));   % Density of Y|X
fxy = fyx.*fx;                    % Density of X|Y
[Xhattheo,loc] = max(fxy);

%% Plot Results
plot(x,fx,'r',x,fn,'b',x,fyx,'k',x,fxy,'m',...
     x(loc),Xhattheo,'o',x(loc),Xhatsim,'*');
axis([min(x),max(x),...
      min(min(min(min(fx,fn),fyx),fxy)),...
      max(max(max(max(fx,fn),fyx),fxy))]);
str = {'$f_X(x)$','$f_N(n)$','$f_{Y\vert X}(y\vert x)$',...
       '$f_{X\vert Y}(x\vert y)$','Theoretical','Simulated'};
legend(str,'Interpreter', 'latex')

But as it can be seen from the result shown in the figure doesn't match. Where am I doing a mistake? Is it in the derivation or in Simulation code? enter image description here

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.