# Probability Density Optimization

I am working on an optimization problem which I am stuck on towards the end.

Essentially, I have two probability density functions in $\mathbb{R}^2$, call them $q(x,y)$ and $p(x,y)$, now I define the objective functional to be:

$C(p,q) = \int p(x,y) \ln\left[\frac{p(x,y)}{q(x,y)}\right] dxdy \cdots (1)$

Now assume we already know what $q(x,y)$ is, ie, it is pre-determined, so my aim is to minimise $(1)$ with respect to the unknown density $p(x,y)$ subject to the following constraints:

$\int_{-\infty}^{\infty} \int_{X^x_d}^{\infty} p(x,y) dx dy = PoD^x$

$\int_{-\infty}^{\infty} \int_{X^y_d}^{\infty} p(x,y) dy dx = PoD^y$

$\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} p(x,y) dx dy = 1$

$p(x,y) \ge 0 \ \forall (x,y) \in \mathbb{R}^2$

Now if we define:

$I_{[X^x_d, \infty)} = \begin{cases} 1 & \text{if} x \ge X^x_d \\ 0 & \text{if} x < X^x_d \end{cases}$

$I_{[X^y_d, \infty)} = \begin{cases} 1 & \text{if} y \ge X^y_d \\ 0 & \text{if} y < X^y_d \end{cases}$

Then our constraints become:

$\int \int p(x,y) I_{[X^x_d, \infty)} dx dy = PoD^x \cdots (2)$

$\int \int p(x,y) I_{[X^y_d, \infty)} dy dx = PoD^y \cdots (3)$

$\int \int p(x,y) dx dy = 1 \cdots (4)$

So I set up the Lagrangian as follows:

$L(p,q) = \int \int p(x,y) \ln[p(x,y)] dxdy - \int \int p(x,y) \ln[q(x,y)] dxdy + \lambda_1 \left[ \int \int p(x,y) I_{[X^x_d, \infty)} dx dy - PoD^x\right] + \lambda_2 \left[\int \int p(x,y) I_{[X^y_d, \infty)} dydx - PoD^y \right] + \mu \left[\int \int p(x,y) dxdy -1\right]$

which can further be simplified:

$L(p,q) = \int \int p(x,y)\left[\ln[p(x,y)] - \ln[q(x,y)]\right] dxdy + \int \int p(x,y) \left[\lambda_1 I_{[X^x_d, \infty)} + \lambda_2 I_{[X^y_d, \infty)} + \mu \right] dxdy - \lambda_1 PoD^x - \lambda_2 PoD^y - \mu$

Now to solve this for $p(x,y)$ I was told to use the calculus of variations and the optimal solution would be:

$\widehat{p(x,y)} = q(x,y) \exp\{-\left[1+\hat{\mu} + \left(\hat{\lambda_1} I_{[X^x_d, \infty)} \right) + \left(\hat{\lambda_2} I_{[X^y_d, \infty)} \right) \right] \}$

Now I am not very familiar with calculus of variations, I only really know how to optimize functionals with respect to algebraic constraints, however I am not sure what to do when the constraints are (multiple) integrals, any help would be appreciated :)

-
I'm not an expert in probability theory and statistics, so forgive me if my question is misguided. Already in your first equation, what happens if $q = 0$ but $p\neq 0$ on a set of positive measure? –  Theo Johnson-Freyd Jul 23 '13 at 3:12
@TheoJohnson-Freyd: In such a case $C(p,q)$ is defined to be $\infty$. –  Ashok Nov 11 '14 at 11:29

Just differentiate $L(p,q)$ with respect to $p(x,y)$, assuming regularities, to get
$$\frac{\partial}{\partial p(x,y)}L(p,q)=1+\ln p(x,y)-\ln q(x,y)+\lambda_1 I_{[X^x_d, \infty)} + \lambda_2 I_{[X^y_d, \infty)}+\mu$$ and equate to zero (the necessary condition in Lagrange multiplier method).
Now using the constraint that $\int \int p(x,y) dx dy = 1$, you see the optimising $p$ is given by $$p(x,y) = c\cdot q(x,y) \exp{\left(-(\lambda_1 I_{[X^x_d, \infty)}+\lambda_2 I_{[X^y_d, \infty)})\right)}$$
where $c$ is the normalising constant.
$C(p,q)$ is a convex functional with respect to function $p$, so the stationary solution can be shown to be the unique minimal point. –  Hans Aug 27 '14 at 16:20