0.16 Appendix: optimization theory

Optimization theory is the branch of applied mathematics whose purpose is to consider a mathematical expression in order to find a set of parameters that either maximize or minimize it. Being an applied discipline, problems usually arise from real-life situations including areas like science, engineering and finance (among many other). This section presents some basic concepts for completeness and is not meant to replace a treaty on the subject. The reader is encouraged to consult further references for more information.

Solution of linear weighted least squares problems

Consider the quadratic problem

$$\underset{h}{\text{min}}{\parallel d-\mathbf{C}h\parallel }_2$$

which can be written as

$$\underset{h}{\text{min}}{\left(d,-,\mathbf{C},h\right)}^T\left(d,-,\mathbf{C},h\right)$$

omitting the square root since this problem is a strictly convex one. Therefore its unique (and thus global) solution is found at the point where the partial derivatives with respect to the optimization variable are equal to zero. That is,

$$\begin{array}{cc}\hfill \frac{\partial }{\partial h}\left\{{\left(d,-,\mathbf{C},h\right)}^T,\left(d,-,\mathbf{C},h\right)\right\}& =\frac{\partial }{\partial h}\left\{d^T,d,-,2,d^T,\mathbf{C},h,+,{\left(\mathbf{C},h\right)}^T,\mathbf{C},h\right\}\hfill \\ & =-2{\mathbf{C}}^{T}d+2{\mathbf{C}}^T\mathbf{C}h=0\hfill \\ & \Rightarrow {\mathbf{C}}^T\mathbf{C}h={\mathbf{C}}^{T}d\hfill \end{array}$$

The solution of [link] is given by

$$h={\left({\mathbf{C}}^T,\mathbf{C}\right)}^{-1}{\mathbf{C}}^{T}d$$

where the inverted term is referred [link] , [link] as the Moore-Pentrose pseudoinverse of ${\mathbf{C}}^T\mathbf{C}$ .

In the case of a weighted version of [link] ,

$$\underset{h}{\text{min}}{\parallel \sqrt{w}\left(d,-,\mathbf{C},h\right)\parallel }_2^2=\sum _k{w}_k{|d_k-C_{k}h|}^2$$

where $C_k$ is the $k$ -th row of $\mathbf{C}$ , one can write [link] as

$$\underset{h}{\text{min}}\left(\mathbf{W},(,d,-,\mathbf{C},h\right){)}^T\left(\mathbf{W},(,d,-,\mathbf{C},h\right))$$

where $\mathbf{W}=\text{diag}\left(\sqrt{w}\right)$ contains the weighting vector $w$ . The solution is therefore given by

$$h={\left({\mathbf{C}}^T,{\mathbf{W}}^T,\mathbf{W},\mathbf{C}\right)}^{-1}{\mathbf{C}}^T{\mathbf{W}}^T\mathbf{W}d$$

Newton's method and the approximation of linear systems in an $l_p$ Sense

Newton's method and $l_p$ Linear phase systems

Consider the problem

$$\underset{a}{\text{min}}g\left(a\right)={\parallel A(\omega ;a)-D\left(\omega \right)\parallel }_p$$

for $a\in {\mathbb{R}}^{M+1}$ . Problem [link] is equivalent to the better posed problem

$$\begin{array}{ccc}\hfill \underset{a}{\text{min}}f\left(a\right)=g{\left(a\right)}^p& =& {\parallel A(\omega ;a)-D\left(\omega \right)\parallel }_p^p\hfill \\ & =& \sum _{i=0}^L\mid C_{i}a-D_i{\mid }^p\hfill \end{array}$$

where $D_i=D\left({\omega }_i\right)$ , ${\omega }_i\in [0,\pi ]$ , $C_i=[C_{i,0},...,C_{i,M}]$ , and

$$\mathbf{C}=\left[\begin{array}{c}{C}_0\\ ⋮\\ C_L\end{array}\right]$$

The $ij$ -th element of $\mathbf{C}$ is given by $C_{i,j}=\text{cos}{\omega }_i(M-j)$ , where $0\le i\le L$ and $0\le j\le M$ . From [link] we have that

$$\nabla f\left(a\right)=\left[\begin{array}{c}\frac{\partial }{\partial a_0}f\left(a\right)\\ ⋮\\ \frac{\partial }{\partial a_M}f\left(a\right)\end{array}\right]$$

where $a_j$ is the $j$ -th element of $a\in {\mathbb{R}}^{M+1}$ and

$$\begin{array}{ccc}\hfill \frac{\partial }{\partial a_j}f\left(a\right)& =& \frac{\partial }{\partial a_j}\sum _{i=0}^L\mid C_{i}a-D_i{\mid }^p\hfill \\ & =& \sum _{i=0}^L\frac{\partial }{\partial a_j}\mid C_{i}a-D_i{\mid }^p\hfill \\ & =& p\sum _{i=0}^L\mid C_{i}a-D_i{\mid }^{p-1}·\frac{\partial }{\partial a_j}\mid C_{i}a-D_i\mid \hfill \end{array}$$

Now,

$$\frac{\partial }{\partial a_j}\mid C_{i}a-D_i\mid =\text{sign}(C_{i}a-D_i)·\frac{\partial }{\partial a_j}(C_{i}a-D_i)=C_{i,j}\text{sign}(C_{i}a-D_i)$$

where Note that

$$\underset{u\left(a\right)\to 0^{+}}{lim}\frac{\partial }{\partial a_j}\mid u\left(a\right){\mid }^p=\underset{u\left(a\right)\to 0^{-}}{lim}\frac{\partial }{\partial a_j}\mid u\left(a\right){\mid }^p=0$$

Therefore the Jacobian of $f\left(a\right)$ is given by

The Hessian of $f\left(a\right)$ is the matrix ${\nabla }^{2}f\left(a\right)$ whose $jm$ -th element ( $0\le j,m\le M$ ) is given by

where adequate substitutions have been made for the sake of simplicity. We have

Note that the partial derivative of $d\left(a\right)$ at $D_i-C_{i}a=0$ is not defined. Therefore

Note that ${\text{sign}}^2(C_{i}a-D_i)=1$ for all $D_i-C_{i}a\ne 0$ where it is not defined. Then

except at $D_i-C_{i}a=0$ where it is not defined.

Based on [link] and [link] , one can apply Newton's method to problem [link] as follows,

• Given $a_0\in {\mathbb{R}}^{M+1}$ , $D\in {\mathbb{R}}^{L+1}$ , $\mathbf{C}\in {\mathbb{R}}^{L+1\times M+1}$
• For $i=0,1,...$
  1. Find $\nabla f\left(a_i\right)$ .
  2. Find ${\nabla }^{2}f\left(a_i\right)$ .
  3. Solve ${\nabla }^{2}f\left(a_i\right)s=-\nabla f\left(a_i\right)$ for $s$ .
  4. Let $a_{+}=a_i+s$ .
  5. Check for convergence and iterate if necessary.

Note that for problem [link] the Jacobian of $f\left(a\right)$ can be written as

where

Also,

where

and

Therefore

From [link] , the Hessian ${\nabla }^{2}f\left(a\right)$ can be expressed as

where

The matrix $\mathbf{C}\in {\mathbb{R}}^{(L+1)\times (M+1)}$ is given by

The matrix $\mathbf{H}={\nabla }^{2}f\left(a\right)$ is positive definite (for $p>1$ ). To see this, consider $\mathbf{H}={\mathbf{K}}^T\mathbf{K}$ where $\mathbf{K}=\mathbf{W}\mathbf{C}$ . Let $z\in {\mathbb{R}}^{M+1}$ , $z\ne 0$ . Then

unless $z\in N\left(\mathbf{K}\right)$ . But since $\mathbf{W}$ is diagonal and $\mathbf{C}$ is full column rank, $N\left(\mathbf{K}\right)=0$ . Thus $z^T\mathbf{H}z\ge 0$ (identity only if $z=0$ ) and so $\mathbf{H}$ is positive definite.

Newton's method and $l_p$ Complex linear systems

Consider the problem

where $\mathbf{A}\in {\mathbb{C}}^{m\times n}$ , $x\in {\mathbb{R}}^n$ and $b\in {\mathbb{C}}^m$ . One can write [link] in terms of the real and imaginary parts of $\mathbf{A}$ and $b$ ,

where $\mathbf{A}=\mathbf{R}+j\mathbf{Z}$ and $b=\alpha +j\gamma$ . The gradient $\nabla e\left(x\right)$ is the vector whose $k$ -th element is given by

where $q_k$ is the row vector whose $i$ -th element is

Therefore one can express the gradient of $e\left(x\right)$ by $\nabla e\left(x\right)=\frac{p}{2}\mathbf{Q}\widehat{g}$ , where $\mathbf{Q}=\left[q_{k,i}\right]$ as above. Note that one can also write the gradient in vector form as follows

The Hessian $\mathbf{H}\left(x\right)$ is the matrix of second derivatives whose $kl$ -th entry is given by

Now,

Substituting [link] and [link] into [link] we obtain

Note that $\mathbf{H}\left(x\right)$ can be written in matrix form as

Therefore to solve [link] one can use Newton's method as follows: given an initial point $x_0$ , each iteration gives a new estimate $x^{+}$ according to the formulas

where $\mathbf{H}\left(x^c\right)$ and $\nabla e\left(x^c\right)$ correspond to the Hessian and gradient of $e\left(x\right)$ as defined previously, evaluated at the current point $x^c$ . Since the $p$ -norm is convex for $1<p<\infty$ , problem [link] is convex. Therefore Newton's method will converge to the global minimizer $x^{☆}$ as long as $\mathbf{H}\left(x^c\right)$ is not ill-conditioned.