# Using least squares for computing gradients

I am developing a code where I am using the least squares method to compute gradients.

Generally, we use least squares to obtain some model based on a set of data (${q_1 \cdots q_N}$) at locations (${x_1 \cdots x_N}$), and use this model to predict $q_c$ at $x_c$.

In my case, I am looking for $\nabla{q}$ at $x_c$ using least squares. Should the $q$ value at $x_c$ be included in a least squares model for its gradient at the location, or only the neighboring values?

I would like use the weighted least squares method. It seems the weighting matrix is defined as a diagonal matrix with the inverse of the variance squared along the main diagonal. The popular choice for the variance seems to be the distance between the neighboring data point and the point in question. So if I include the $q_c$ value at $x_c$ in my least squares model, then the variance would be zero for this data point, which would be problematic unless I add some $\epsilon$ value to the variance before taking it's inverse.

The least squares equation is

\begin{align} \nabla q = q^TZM^T\nabla\phi(x) \\ \\ M=(Z^TZ)^{-1}\\ \\ \phi = \begin{bmatrix} 1 & x & y & z & xy & xz & yz & 0.5x^2 & 0.5y^2 & 0.5z^2 \\ \end{bmatrix}\\ \\ Z_{unweight} = \begin{bmatrix} 1 & x_1 & y_1 & z_1 & x_1y_1 & x_1z_1 & y_1z_1 & 0.5x_1^2 & 0.5y_1^2 & 0.5z_1^2 \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots &\vdots & \vdots \\ 1 & x_N & y_N & z_N & x_Ny_N & x_Nz_N & y_Nz_N & 0.5x_N^2 & 0.5y_N^2 & 0.5z_N^2 \\ \end{bmatrix} \\ \\ q = \begin{bmatrix} q_1 & \cdots & q_N \end{bmatrix} \end{align}

$q_c$ is known, but I am not sure if this should be included in the above model.

• Just to be clear, you're finding the gradient of $y(x_{c))$ with respect to what independent variables? – Brian Borchers Mar 25 '18 at 21:44
• Yes, I am finding the gradient of y at $x_c$. The independent variable are the x's, which represents spatial coordinates in my case. – David Mar 25 '18 at 23:57
• So you want the gradient of $y(x_{c})$ with respect to $x_{1}$, $x_{2}$, $\ldots$, $x_{n}$, keeping $y_{1}$, $y_{2}$, $\ldots$, $y_{n}$ constant? What's the equation that relates $y(x_{c})$ to $x_{1}$, $x_{2}$, $\ldots$, $x_{n}$? – Brian Borchers Mar 26 '18 at 0:08
• That is correct, I will update the original post, with the least squares equation. – David Mar 26 '18 at 0:22
• Updated it. Note that I changed 'y' to 'q' because I am using 'y' to represent a spatial coordinate in the above eqn. – David Mar 26 '18 at 0:32

You are asking the following question. Given $(q_0,q_1,q_2,\ldots,q_N)$ find $\nabla q$ at $\vec{r}=\vec{r}_0$. You can do this two ways.

I. Assume

$$q(x,y,z) = a + \vec{b} \cdot (\vec{r} - \vec{r}_0)$$ and find $a$ and $\vec{b}$ by $$\min_{a,\vec{b}}\sum_{i=0}^N(q(x_i,y_i,z_i) - q_i)^2$$

II. Assume

$$q(x,y,z) = q_0 + \vec{b} \cdot (\vec{r} - \vec{r}_0)$$ and find $\vec{b}$ by $$\min_{\vec{b}}\sum_{i=1}^N(q(x_i,y_i,z_i) - q_i)^2$$

In either case, $\vec{b}$ will give an approximation to $\nabla q$ at $\vec{r}_0$.

In (II), you have $q(x_0,y_0,z_0)=q_0$ but in (I) you only have $q(x_0,y_0,z_0) \approx q_0$.

If you want to put in weights based on some inverse distance in above norm, then even in (I) you are forcing $q(x_0,y_0,z_0)=q_0$ to hold. In this case, you can just use approach (II) after putting in the weights appropriately in the norm.