How to interpret a specific data transformation? I came across this specific data transformation in the context of a physics application, which by itself is rather complex and hence out of the scope of this question. However since this transformation is applicable to general data sets, I wondered if there's a way to describe general properties of or to interpret this specific data transformation.
Suppose a set of $N$ data points $\{x_1, \ldots, x_N\}$ where $N = m\cdot n$ for two integer numbers $m, n > 1$. Then the transformation from $x_i$ to $v_i$ is given by the following system of linear equations (with $m+1$ rows and $m\cdot n$ columns):
$$
\begin{pmatrix}
1 & \cdots & 1 & 0 & \cdots & 0 & \cdots\cdots & 0 & \cdots & 0 \\
0 & \cdots & 0 & 1 & \cdots & 1 & \cdots\cdots & 0 & \cdots & 1 \\
\vdots & & \vdots & \vdots & & \vdots & & \vdots & & \vdots \\
0 & \cdots & 0 & 0 & \cdots & 0 & \cdots\cdots & 1 & \cdots & 1 \\
x_{1} & \cdots & x_{n} & x_{n+1} & \cdots & x_{2n} & \cdots\cdots & x_{(m-1)\cdot n+1} & \cdots & x_{m\cdot n}
\end{pmatrix}
\cdot
\begin{pmatrix}
v_1 \\
\vdots \\
\vdots \\
\vdots \\
v_{m\cdot n}
\end{pmatrix}
=
\begin{pmatrix}
1 \\
0 \\
\vdots \\
0 \\
0
\end{pmatrix}
$$
The first $m$ rows of the above coefficient matrix $A$ can be described as follows:
$$
A_{i,j} = \begin{cases}
1 \quad \textrm{if} \;\; (i-1)\cdot n < j \leq i\cdot n \\
0 \quad \textrm{otherwise}
\end{cases}
$$
for $i, j \in \{1, \ldots N\}$.
The last row of $A$ simply holds all the data points $x_i$.
The r.h.s. vector $b$ is all zeros except for the first row which equals $1$ (in fact any other row except the last one could be chosen to hold the $1$; it's important that there is only one non-zero entry though).
This represents an underdetermined system of linear equations for which a solution can be found by minimizing the L2-norm. The solution $V^{\ast}$ is then given by the following $N\times N$ system of equations:
$$
A^T A V^{\ast} = A^T b
$$
where $A^T$ is the transpose of $A$ and $b$ is the r.h.s. vector of the original system of equations. The solution $V^{\ast}$ can be found as the least-squares solution by pseudo-inversion of the matrix $A^T A$ via singular value decomposition.
Now what I'm interested in are general properties of this transformation and how the transformed data $v_i$ reflects the original data $x_i$. Even when describing all of the involved steps, I'm having difficulties interpreting the results.
Example application
This is some sample Python code which applies the above transformation to different data sets (Gaussian, sin, cos, exp).
import matplotlib.pyplot as plt
import numpy as np

m = 10
n = 100
N = m*n

x = np.exp(-np.linspace(-2, 2, N)**2)
# x = np.sin(np.linspace(0, 2*np.pi, N))
# x = np.cos(np.linspace(0, 2*np.pi, N))
# x = np.exp(-np.linspace(0, 4, N))

A = np.zeros((m+1, N))
for i in range(m):
    A[i, i*n:(i+1)*n] = 1
A[-1, :] = x

b = np.zeros(m+1)
b[0] = 1

v, *info = np.linalg.lstsq(A, b, rcond=None)

fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(10, 4))
ax1.set_title('Data X')
ax1.plot(x)
ax2.set_title('Transformed V')
ax2.plot(v, '--o', lw=0.7, ms=1)

plt.show()

Output for Normal Distribution

Output for Sine Function

Output for Cosine Function

Output for Exponential Function

 A: This transformation results in something very close to least squares residuals.
At first step you can regress vextor $\pmb x$ on the set of $m$ dummies. In such case all other than $\pmb x$ columns of the $\pmb A^T$matrix would be OLS matrix $\pmb X$. When performing such action the plot of negative residuals will be very similar to your transformation results, but not exactly the same.
However, if we want to achieve really close result, all residual values require multiplying by a constant dependent both on $m$, and values included in the data. Also, the dummy variable which corresponds to the nonzero value in your $b$ vector is moved upwards by $1/n$. Here, it seems to be always this value.
If this reasoning is true (see numerical example below for doubts), the interpretation would be equal to OLS/ANOVA residuals.

To see numerical example, after executing code from the question please try this one:
import statsmodels.formula.api as sm

# Creating OLS data:
data = pd.DataFrame(A.transpose(), columns = ["v" + str(id) for id in range(m+1)])

# Estimating OLS:
f = "v" + str(m) + "~"
for i in range(m-1):
    f = f + "+v" + str(i+1)
result = sm.ols(formula=f, data=data).fit()

# Calculating unknown constant:
constant = v[m*n-1]/result.resid[m*n-1]

# Calculating values of OLS transformation:
reg_v = constant * result.resid + data["v0"]/n

# Differences between Original and OLS:
v/reg_v
v-reg_v

# Upgreaded plot:
fig, (ax1, ax2, ax3) = plt.subplots(ncols=3, figsize=(10, 4))
ax1.set_title('Data X')
ax1.plot(x)
ax2.set_title('Transformed V')
ax2.plot(v, '--o', lw=0.7, ms=1)
ax3.set_title('OLS, const = ' + str(constant))
ax3.plot(reg_v, '--o', lw=0.7, ms=1)
plt.show()

It recreates the transformation using simple OLS estimation resulting in very similar result:
The key idea here is not to try to derive the constant manually, but to numerically gather it from the data. This however results in that the OLS idea and original transformation results differ at small values ($10^{-15}$ to $10^{-20}$). This could be numerical problem with inversion of different matrices or... the constant could not be so constant.

Another solution I can think of would be forcing the set of equations into OLS notation structure. However I am not sure at what point this could be done. Probably somewhere closer to optimisation, what was not included in your description. If my reasoning is true, and visible differences are in fact numerical in nature, such step should be possible to perform.
