# Constrained Least Squares Approximation of a Plane

I have a situation in which I wish to fit a plane $ax+by+cz=d$ to a set of points in 3D space using least squares approximation. However I also know for a fact that one of these points $p=[p_1,p_2,p_3]^T$ lies exactly on the plane. I would therefore like to constrain the plane to fit this constraint.

My idea to solve this would be to define the plane in normal form, so that $[x-p]* n=0$, where $p$ is my point on the plane and $n$ is the normal vector. Then for the cartesian form of the plane $a=n_1, b=n_2, c=n_3, d=n\cdot p=n_1p_1+n_2p_2+n_3p_3$. Substituting all that back into the plane equation gives me $n_1(x-p)+n_2(y-p)+n_3(z-p)=0$, which as I add $x,y,z$ data points will provide me with the overdetermined system I want.

Is this approach correct? Are there better alternatives? Many thanks!

This is a linear least squares problem with linear equality constraints which you can formulate as: find $\mathbf{x}^*$ such that $$\tag{1} \|\mathbf{b}-\mathbf{A}\mathbf{x}^*\|_2=\min\{\|\mathbf{b}-\mathbf{A}\mathbf{x}\|_2:\mathbf{C}\mathbf{x}=\mathbf{d}\}.$$
In this problem, I have $m$ points supposedly in a plane and let's say that the point $k$ is where I want an equality. Note that if $d\neq 0$, we can scale the plane equation such that $ax+by+cz=1$. Then $$\mathbf{A}=\begin{bmatrix} x_1 & y_1 & z_1 \\ \vdots & \vdots & \vdots \\ x_{k-1} & y_{k-1} & z_{k-1} \\ x_{k+1} & y_{k+1} & z_{k+1} \\ \vdots & \vdots & \vdots \\ x_m & y_m & z_m \end{bmatrix}, \quad \mathbf{x}=\begin{bmatrix}a\\b\\c\end{bmatrix}, \quad \mathbf{b}=\begin{bmatrix}1 \\ \vdots \\ 1\end{bmatrix},$$ $$\mathbf{C}=\begin{bmatrix}x_k & y_k & z_k\end{bmatrix}, \quad \mathbf{d}=[1].$$
The most direct way how to solve this is using the Lagrangian $$f(\mathbf{x})=\frac{1}{2}(\mathbf{b}-\mathbf{A}\mathbf{x})^T(\mathbf{b}-\mathbf{A}\mathbf{x})+\boldsymbol{\lambda}^T(\mathbf{C}\mathbf{x}-\mathbf{d}).$$ Differentiating w.r.t. $\mathbf{x}$ and $\boldsymbol{\lambda}$ gives a linear system $$\begin{bmatrix} \mathbf{A}^T\mathbf{A} & \mathbf{C}^T \\ \mathbf{C} & 0 \end{bmatrix} \begin{bmatrix} \mathbf{x}\\\boldsymbol{\lambda} \end{bmatrix} = \begin{bmatrix} \mathbf{A}^T\mathbf{b} \\ \mathbf{d} \end{bmatrix}.$$
Note that we need to eliminate one parameter from the plane equation. If we included $d$ into the vector of unknowns $\vec{x}$, the optimal solution would be zero ($a=b=c=d=0$).