User:Hikenstuff/sandbox

Polynomial least squares

Limited aspects of the general subject of polynomial least squares are also addressed under many other titles: polynomial regression, curve fitting, linear regression, least squares, ordinary least squares, simple linear regression, linear least squares, approximation theory and method of moments. Polynomial least squares has application in radar trackers, estimation theory, signal processing, statistics, and econometrics.

There are two fundamental applications of polynomial least squares:

(1) To approximate a complicated function or set of observations with a simple low degree polynomial. This is commonly used in statistics and econometrics to fit a scatter plot (often called a scatter gram) with a straight line in the form of a first degree polynomial. This is addressed under afore mentioned titles.

(2) To estimate an assumed underlying deterministic polynomial that is corrupted with statistically described additive errors (generally called noise in engineering) from observations or measurements. This is commonly used in target tracking in the form of the Kalman filter, which is effectively a recursive implementation of polynomial least squares. Estimating an assumed underlying deterministic polynomial can be used econometrics as well. Processing noisy measurements is uniquely addressed here.

The term "estimate" is derived from statistical estimation theory and is perhaps better suited when assuming that a polynomial is corrupted with statistical measurement or observation errors. The term "approximate" is perhaps better suited when no statistical measurement or observation errors are assumed, such as when conventionally fitting a scatter plot or complicated function.

In effect, both applications produce average curves as generalizations of the common average of a set of numbers, which is equivalent to zero degree polynomial regression and least squares.

Polynomial least squares estimate of a deterministic first degree polynomial corrupted with observation errors

Assume the deterministic function y with unknown coefficients 𝜶 and 𝜷 as follows:

$y=\alpha +\beta t$

which is corrupted with an additive stochastic process $\epsilon$, described as an error (noise in tracking) written as

$z=y=\alpha +\beta t+\epsilon$

Given samples $z_{n}$ where the subscript $n$ is the sample index, the problem is to apply polynomial least squares to estimate y(t), and to determine its variance along with its expected value.

Assumptions and definitions

(1) The error $\epsilon$ is modeled as a zero mean stochastic process, samples of which are random variables that are uncorrelated and assumed to have identical probability distributions (specifically same mean and variance), but not necessarily Gaussian, treated as inputs to polynomial least squares. Stochastic processes and random variables are described only by probability distributions.

(2) Polynomial least squares is modeled as a linear signal processing "system" which processes statistical inputs deterministically, the output being the linearly processed empirically determined statistical estimate, variance, and expected value.

(3) Polynomial least squares processing produces deterministic moments (analogous to mechanical moments), which may be considered as moments of sample statistics, but not of statistical moments.

Polynomial least squares and the orthogonality principle

Approximating a function z(t) with a polynomial

${\hat {z}}(t)=\sum _{j=1}^{J}a_{j}t^{j-1}$

where hat (^) denotes the estimate and (J-1) is the polynomial degree, can be performed by applying the orthogonality principle. The error e in the sum of the squared errors can be written as

$e=\sum _{n=1}^{N}(z_{n}-{\hat {z}}_{n})^{2}$

According to the orthogonality principle , e is minimum when the error ($z$-${\hat {z}}$) is orthogonal to the estimate ${\hat {z}}$, that is

$\sum _{n=1}^{N}(z_{n}-{\hat {z}}_{n}){\hat {z}}_{n}=0$

This can be described as the orthogonal projection of the data $z_{n}$ onto a solution in the form of the polynomial ${\hat {z}}(t)$. For N > J-1, orthogonal projection yields the standard overdetermined system of equations (often called normal equations) used to compute the coefficients in the polynomial approximation. The minimum e is then

$e_{min}=\sum _{n=1}^{N}(z_{n}-{\hat {z}}_{n})z_{n}$

The advantage of using orthogonal projection is that $e_{min}$ can be determined for use in the polynomial least squares processed statistical variance of the estimate.

The empirically determined polynomial least squares output of a first degree polynomial corrupted with a observation errors

To fully determine the output of polynomial least squares, a weighting function describing the processing must first be structured and then the statistical moments can be computed.

The weighting function describing the linear polynomial least squares "system"

Given estimates of the coefficients 𝜶 and 𝜷 from polynomial least squares, the weighting function $w_{n}(\tau )$ can be formulated to estimate the unknown y as follows:

${\hat {y}}(\tau )={\frac {1}{N}}\sum _{n=1}^{N}z_{n}w_{n}(\tau )={\frac {1}{N}}\sum _{n=1}^{N}(\alpha +\beta t_{n}+\epsilon _{n})w_{n}(\tau )$

where N is the number of samples, $z_{n}$ are random variables as samples of the stochastic $z$ (noisy signal), and the first degree polynomial data weights are

$w_{n}(\tau )\equiv {\frac {[{\bar {t^{2}}}-{\bar {t}}t_{n}+(t_{n}-{\bar {t}})\tau ]}{({\bar {t^{2}}}-{\bar {t}}^{2})}}$

which represent the linear polynomial least squares "system" and describe its processing. The Greek letter 𝜏 is the independent variable t when estimating the dependent variable y after data fitting has been performed. (The letter 𝜏 is used to avoid confusion with t before and sampling during polynomial least squares processing.) The overbar ( ¯ ) defines the deterministic centroid of $u_{n}$ as processed by polynomial least squares – i.e., it defines the deterministic first order moment, which may be considered a sample average, but does not here approximate a first order statistical moment:

${\bar {u}}{\overset {\underset {\mathrm {def} }{}}{=}}{\frac {1}{N}}\sum _{n=1}^{N}u_{n}$

Empirically determined statistical moments

Applying $w_{n}(\tau )$ yields

${\hat {y}}(\tau )={\hat {\alpha }}+{\hat {\beta }}\tau$

where

${\hat {\alpha }}={\frac {({\bar {z}}{\bar {t^{2}}}-{\bar {zt}}{\bar {t}})}{({\bar {t^{2}}}-{\bar {t}}^{2})}}=\alpha +{\frac {({\bar {\epsilon }}{\bar {t^{2}}}-{\bar {{\epsilon }t}}{\bar {t}})}{({\bar {t^{2}}}-{\bar {t}}^{2})}}$

and

${\hat {\beta }}={\frac {({\bar {zt}}-{\bar {z}}{\bar {t}})}{({\bar {t^{2}}}-{\bar {t}}^{2})}}=\beta +{\frac {({\bar {\epsilon t}}-{\bar {\epsilon }}{\bar {t}})}{({\bar {t^{2}}}-{\bar {t}}^{2})}}$

As linear functions of the random variables $\epsilon _{n}$, both coefficient estimates ${\hat {\alpha }}$ and ${\hat {\beta }}$ are random variables. In the absence of the errors $\epsilon _{n}$, ${\hat {\alpha }}=\alpha$ and ${\hat {\beta }}=\beta$, as they should to meet that boundary condition.

Because the statistical expectation operator E[•] is a linear function and the sampled stochastic process errors $\epsilon _{n}$ are zero mean, the expected value of the estimate ${\hat {y}}$ is the first order statistical moment as follows:

$E[{\hat {y}}(\tau )]=\alpha +\beta \tau +{\frac {1}{N}}\sum _{n=1}^{N}E[\epsilon _{n}]w_{n}(\tau )=\alpha +\beta \tau =\alpha +\beta t$

The statistical variance in ${\hat {y}}$ is given by the second order statistical central moment as follows:

$\sigma _{\hat {y}}^{2}=E[({\hat {y}})-E[{\hat {y}}])^{2}]={\frac {1}{N}}{\frac {1}{N}}\sum _{n=1}^{N}\sum _{i=1}^{N}w_{n}(\tau )E[\epsilon _{n}\epsilon _{i}]w_{i}(\tau )$ $=\sigma _{\epsilon }^{2}{\frac {1}{N}}{\frac {1}{N}}\sum _{n=1}^{N}\sum _{i=1}^{N}w_{n}^{2}(\tau )$

because

$\sum _{i=1}^{N}E[\epsilon _{n}\epsilon _{i}]w_{i}(\tau )=\sigma _{\epsilon }^{2}w_{n}(\tau )$

where $\sigma _{\epsilon }^{2}$ is the statistical variance of random variables $\epsilon _{n}$; i.e., $E[\epsilon _{n}\epsilon _{i}]=\sigma _{\epsilon }^{2}$ for i = n and (because $\epsilon _{n}$ are uncorrelated) $\sigma _{\epsilon }^{2}=0$ for $i\neq n$

Carrying out the multiplications and summations in $\sigma _{\hat {y}}^{2}$ yields

$\sigma _{\hat {y}}^{2}=\sigma _{\epsilon }^{2}{\frac {({\bar {t^{2}}}-2{\bar {t}}\tau +\tau ^{2})}{N({\bar {t^{2}}}-{\bar {t}}^{2})}}$

Measuring or approximating the statistical variance of the random errors

In a hardware system, such as a tracking radar, the measurement noise variance $\sigma _{\epsilon }^{2}$ can be determined from measurements when there is no target return – i.e., by just taking measurements of the noise alone.

However, if polynomial least squares is used when the variance $\sigma _{\epsilon }^{2}$ is not measureable (such as in econometrics or statistics), it can be estimated with observations in $e_{min}$ from orthogonal projection as follows:

$\sigma _{\epsilon }^{2}\approx {\hat {\sigma _{\epsilon }^{2}}}=({\bar {z^{2}}}-{\hat {\alpha }}{\bar {z}}-{\hat {\beta }}{\bar {zt}})$

As a result, to the first order approximation from the estimates ${\hat {\alpha }}$ and ${\hat {\beta }}$ as functions of sampled $z$ and $t$

$\sigma _{\hat {y}}^{2}\approx {\bigg [}{\frac {({\bar {z^{2}}}-{\bar {z}}^{2})}{({\bar {t^{2}}}-{\bar {t}}^{2})}}-{\Biggl (}{\frac {({\bar {zt}}-{\bar {z}}{\bar {t}})}{({\bar {t^{2}}}-{\bar {t}})}}{\Biggl )}^{2}{\bigg ]}{\frac {({\bar {t^{2}}}-2{\bar {t}}\tau +\tau ^{2})}{N}}$

which goes to zero in the absence of the errors $\epsilon _{n}$, as it should to meet that boundary condition.

As a result, the samples $z_{n}$ (noisy signal) are considered to be the input to the linear polynomial least squares "system" which transforms the samples into the empirically determined statistical estimate ${\hat {y}}(\tau )$, the expected value $E[{\hat {y}}]$, and the variance $\sigma _{\hat {y}}^{2}$.

Properties of polynomial least squares modeled as a linear "system"

(1) The empirical statistical variance $\sigma _{\hat {y}}^{2}$ is a function of $\sigma _{\epsilon }^{2}$, N and $\tau$. Setting the derivative of $\sigma _{\hat {y}}^{2}$ with respect to $\tau$ equal to zero shows the minimum to occur at $\tau ={\bar {t}}$; i.e., at the centroid (sample average) of the samples $t_{n}$. The minimum statistical variance thus becomes ${\frac {\sigma _{\epsilon }^{2}}{N}}$. This is equivalent to the statistical variance from polynomial least squares of a zero degree polynomial – i.e., of the centroid (sample average) of $\alpha$.

(2) The empirical statistical variance $\sigma _{\hat {y}}^{2}$ is a function of the quadratic $\tau ^{2}$ . Moreover, the further $\tau$ deviates from ${\bar {t}}$ (even within the data window), the larger is the variance $\sigma _{\hat {y}}^{2}$ due to the random variable errors $\epsilon _{n}$ . The independent variable $\tau$ can take any value on the $t$ axis. It is not limited to the data window. It can extend beyond the data window – and likely will at times depending on the application. If it is within the data window, estimation is described as interpolation. If it is outside the data window, estimation is described as extrapolation. It is both intuitive and well known that the further is extrapolation, the larger is the error.

(3) The empirical statistical variance $\sigma _{\hat {y}}^{2}$ due to the random variable errors $\epsilon _{n}$ is inversely proportional to N. As N increases, the statistical variance decreases. This is well known and what filtering out the errors $\epsilon _{n}$ is all about. The underlying purpose of polynomial least squares is to filter out the errors to improve estimation accuracy by reducing the empirical statistical estimation variance. In reality, only two data points are required to estimate $\alpha$ and $\beta$; albeit the more data points with zero mean statistical errors included, the smaller is the empirical statistical estimation variance as established by N samples.

(4) There is an additional issue to be considered when the noise variance is not measureable: Independent of the polynomial least squares estimation, any new observations would be described by the variance $\sigma _{\epsilon }^{2}\approx {\hat {\sigma _{\epsilon }^{2}}}=({\bar {z^{2}}}-{\hat {\alpha }}{\bar {z}}-{\hat {\beta }}{\bar {zt}})$.

Thus, the polynomial least squares statistical estimation variance $\sigma _{\hat {y}}^{2}$ and the statistical variance of any new sample in $\sigma _{\epsilon }^{2}$ would both contribute to the uncertainty of any future observation. Both variances are clearly determined by polynomial least squares in advance.

(5) This concept also applies to higher degree polynomials. However, the weighting function $w_{n}(\tau )$ is obviously more complicated. In addition, the estimation variances increase exponentially as polynomial degrees increase linearly (i.e., in unit steps). However, there are ways of dealing with this as described in .

The synergy of integrating polynomial least squares with statistical estimation theory

Modeling polynomial least squares as a linear signal processing "system" creates the synergy of integrating polynomial least squares with statistical estimation theory to deterministically process corrupted samples of an assumed polynomial. In the absence of the error ε, statistical estimation theory is irrelevant and polynomial least squares reverts back to the conventional approximation of complicated functions and scatter plots.