2

I have an overdetermined system of linear equations A*X=B. Where X is a column matrix with N unknowns. The least square solution of this system can be obtained in Matlab as X = mldivide(A,B). However, I want a solution X, such that the sum of second order differences between different elements in X is minimum. That is, to find a solution X, such that the following error term is minimum:

$E = ||\mathbf{A}\mathbf{X}-\mathbf{B}||^2+\sum_{j=2}^{N-1}[X_{j-1}-2*X_j+X_{j+1}]^2$

DSingh
  • 41
  • 6

2 Answers2

3

What you're describing sounds exactly like what is sometimes called a "graph ridge"/"graphridge" regression, which is like a ridge regression (if that's something you're familiar with) but with a penalty not on the (squared) regression coefficients themselves, but on the (squared) differences between coefficients that are connected in an assumed graph structure. In your case, you'd let $X_j$ and $X_k$ be connected in this graph if $|k-j|=1$.

This type of regression is convenient because like the "regular" ridge regression, it has a closed-form solution. For more info see e.g. this paper which makes use of this method (esp. equation 14).

Ruben van Bergen
  • 6,511
  • 1
  • 20
  • 38
2

After getting hints from the comments and lot of literature survey, I have solved the above problem by myself. The above said error can be minimized by means of Tikhonov regularization of second order. Error term in general Tikhonov regularization is given as

$$E = {\|\textbf{A}\textbf{X}-\textbf{B}\|^2_2+\lambda\|LX\|^2_2}$$

Where $\lambda$ is damping factor and $L$ is regularization matrix. $L$ can have many forms. If $L$ is replaced by an identity matrix $I$, it is known as ridge regresion. Another common practice is to replace $L$ by finite order difference matrix. For our concerned error term, $L$ will replaced by second order difference matrix, i.e. $$L= \left[ {\begin{array}{cccc} -1& 2 & -1 & & & 0\\ & -1 & 2 & -1 & & \\ & & ... & ... & ... & \\ & & & -1 & 2 & -1\\ \end{array} } \right]$$

The solution by Tikhonov regularization can be then obtained by solving the following linear system:

$$\left[ \begin{array}{cc} A \\ \lambda L\end{array}\right]X = \left[ \begin{array}{cc} B \\ 0 \end{array}\right]$$ You will have to decide on value of $\lambda$ (see L Curve in references).

References:

  1. Noschese, Silvia, and Lothar Reichel. "Inverse problems for regularization matrices." Numerical Algorithms (2012): 1-14.
  2. Calvetti, Daniela, Lothar Reichel, and Abdallah Shuibi. "Invertible smoothing preconditioners for linear discrete ill-posed problems." Applied numerical mathematics 54.2 (2005): 135-149.
  3. Hansen, Per Christian. The L-curve and its use in the numerical treatment of inverse problems. IMM, Department of Mathematical Modelling, Technical Universityof Denmark, 1999.
DSingh
  • 41
  • 6
  • FYI this would typically be called [regularized](https://stats.stackexchange.com/questions/250722/the-origin-of-the-term-regularization) regression, rather than constrained. (Although Tikhonov can be expressed as an inequality constraint.) – GeoMatt22 May 04 '17 at 04:48
  • @GeoMatt22 Thanks. I have edited the question title to fit better into the context. – DSingh May 04 '17 at 05:38
  • Glad you solved your problem, but in case it wasn't clear: this is exactly the method I was talking about in my answer, and that is used in the paper I referenced. – Ruben van Bergen May 04 '17 at 09:51