Hat man mehr als n Wertepaare, so ist das Gleichungssystem überbestimmt. Falls es sich bei den Wertepaaren um Meßwerte handelt, sind sie oftmals mit Fehlern behaftet, so daß man eine Funktion, die alle Werte genau interpoliert, überhaupt nicht finden kann. In diesen Fällen kann man sich von der Interpolationsforderung lösen und fordern, daß die Abweichung der Funktionswerte der zu findenden Funktion von den gemessenen y-Werten so klein wie möglich ist.

Die "Abweichung" kann hierbei im Grunde recht unterschiedlich definiert werden. Plausibel erscheint eine Definition, in der ein größerer Fehler überproportional mehr wiegt als ein kleiner; außerdem darf das Fehlermaß nicht vorzeichenbehaftet sein, denn sonst würde ein negativer Fehler einen betragsgleichen positiven Fehler ausgleichen. Die einfachen Quadrate der Differenzen zwischen Meßwert und Funktionswert erfüllen diese beiden Forderungen. (f(x) - y)² ist also ein geeignetes Maß für die Abweichung des bei x gemessenen Wertes y vom Wert, den die Funktion f für x ausgibt.

Man sucht also eine Funktion für die die Summe der Quadrate der einzelnen Differenzen zwischen Funktions- und Meßwerten minimal wird. Für diese Summe schreiben wir S:

S = (f(x₁)-y₁)² + (f(x₂)-y₂)² + ... + (f(x_n)-y_n)²
bzw.
S = ∑((f(x_i) - y_i)²)

Die notwendige Bedingung dafür, daß diese Summe minimal ist, ist daß ihre erste Ableitung 0 wird. Nun hängt ja diese Summe nicht von x, sondern von der Wahl der Koeffizienten innerhalb der Funktion f ab. Zur Erinnerung: f ist eine Linearkombination a₁φ₁(x) + a₂φ₂(x) + ... mit den zunächst unbekannten Koeffizienten a_i. S ist also eine Funktion von a₁, a₂, ..., a_m. Wir suchen also Koeffizienten a_i mit S'(a₁, a₂, ..., a_m) = 0.

Ohne in die Tiefe zu gehen, wie die Ableitung einer mehrdimensionalen Funktion bestimmt wird, und ohne auf das Gaußsche Fehlerfortpflanzungsgesetz einzugehen, das hier eigentlich dahinter steckt, gehen wir einmal "naiv" davon aus, daß wir nur den bestimmten Koeffizienten a_j suchen, während die übrigen a_i festliegen. Dann hängt S von a_j ab, ist also quasi eine Funktion von a_j, und die Ableitung von S nach a_j ist mit der Anwendung von Kettenregel und Summenregel leicht anzugeben: S' = ∂∑/∂a_j = ∑( 2(f(x_i) - y_i)·f'(x_i) )

Auch die Ableitung f' nach a_j ist leicht zu bestimmen, denn f ist ja eine Linearkombination, in der a_j nur in einem Summand auftritt, und das linear. Die Ableitung ∂f'/∂a_j ist also gerade der zweite Faktor von a_j·φ_j(x), also φ_j(x). Alle anderen Summanden von f sind konstant, weil sie nicht von a_j abhängen, und fallen in der Ableitung weg. Somit ergibt sich für S':

S' = ∑( 2(f(x_i) - y_i)·φ_j(x_i) )

Dies läßt sich vereinfachen, indem 2 ausgeklammert wird: S' = 2·∑((f(x_i) - y_i)·φ_j(x_i))

Aus S'=0 ergibt sich somit die Gleichung 2·∑((f(x_i) - y_i)φ_j(x_i)) = 0

Wir dividieren durch 2, multiplizieren die innere Klammer aus, teilen die Summe auf und schreiben die y nach rechts:
∑((f(x_i) - y_i)φ_j(x_i)) = 0
∑((f(x_i)·φ_j(x_i) - y_i·φ_j(x_i)) = 0
∑f(x_i)·φ_j(x_i) - ∑y_i·φ_j(x_i) = 0
∑f(x_i)·φ_j(x_i) = ∑y_i·φ_j(x_i)

Wegen f(x) = a₁φ₁(x) + a₂φ₂(x) + ... + a_mφ_m(x), kann man schreiben
∑(a₁φ₁(x_i) + a₂φ₂(x_i) + ... + a_mφ_m(x_i))·φ_j(x_i) = ∑y_i·φ_j(x_i)
∑(a₁φ₁(x_i)·φ_j(x_i) + a₂φ₂(x_i)·φ_j(x_i) + ... + a_mφ_m(x_i)·φ_j(x_i)) = ∑y_i·φ_j(x_i)
∑a₁φ₁(x_i)·φ_j(x_i) + ∑a₂φ₂(x_i)·φ_j(x_i) + ... + ∑a_mφ_m(x_i)·φ_j(x_i) = ∑y_i·φ_j(x_i)
a₁·∑φ₁(x_i)·φ_j(x_i) + a₂·∑φ₂(x_i)·φ_j(x_i) + ... + a_m·∑φ_m(x_i)·φ_j(x_i) = ∑y_i·φ_j(x_i)

Gauß führte zur Vereinfachung die Schreibweise [φ₁·φ₂] für die Summe ∑φ₁(x_i)·φ₂(x_i) über alle i ein.
Die obige Gleichung läßt sich damit einfacher schreiben: a₁[φ₁·φ_j] + a₂[φ₂·φ_j] + ... + a_m[φ_m·φ_j] = [y·φ_j(x)]. Dies ist eine lineare Gleichung über die unbekannten Koeffizienten a.

Führt man dieses Verfahren für jedes a durch, ergeben sich m lineare Gleichungen, mit denen die Koeffizienten a₁ bis a_m eindeutig bestimmt werden können:

a₁[φ₁·φ₁] + a₂[φ₂·φ₁] + ... + a_m[φ_m·φ₁] = [y·φ₁(x)]
a₁[φ₁·φ₂] + a₂[φ₂·φ₂] + ... + a_m[φ_m·φ₂] = [y·φ₂(x)]
...
a₁[φ₁·φ_m] + a₂[φ₂·φ_m] + ... + a_m[φ_m·φ_m] = [y·φ_m(x)]

Zum Nachvollziehen rechne ich das nochmal anhand eines konkreten Falles durch. Gesucht sei eine kubische Funktion f(x) = ax³ + bx² + cx + d, die eine Reihe von n Meßwertpaaren (x₁|y₁), (x₂|y₂), ... (x_n|y_n) optimal, d.h. mit minimaler Summe der Fehlerquadrate, approximiert. (Da vier Koeffizienten gesucht sind, sei n>4.)

Die Summe S der Fehlerquadrate ist S = ∑(ax_i³ + bx_i² + cx_i + d - y_i)².

Das wird minimal bei S'=0. Da S von jedem der vier Koeffizienten abhängt, suchen wir Koeffizienten, für die die partielle Ableitung S' nach jedem Koeffizient 0 wird:
Die Ableitung von S nach a ist: ∂S/∂a = ∑(2(ax_i³ + bx_i² + cx_i + d - y_i)·x_i³)
Die Ableitung von S nach b ist: ∂S/∂a = ∑(2(ax_i³ + bx_i² + cx_i + d - y_i)·x_i²)
Die Ableitung von S nach c ist: ∂S/∂a = ∑(2(ax_i³ + bx_i² + cx_i + d - y_i)·x_i)
Die Ableitung von S nach d ist: ∂S/∂a = ∑(2(ax_i³ + bx_i² + cx_i + d - y_i))

Setzt man die Ableitung nach a in der ersten Zeile gleich 0, dividiert durch 2 und formt um, gewinnt man nach und nach:
∑(ax_i³·x_i³ + bx_i²·x_i³ + cx_i·x_i³ + d·x_i³ - y_i·x_i³) = 0
∑ax_i³·x_i³ + ∑bx_i²·x_i³ + ∑cx_i·x_i³ + ∑d·x_i³ - ∑y_i·x_i³ = 0
∑ax_i³·x_i³ + ∑bx_i²·x_i³ + ∑cx_i·x_i³ + ∑d·x_i³ = ∑y_i·x_i³
a·∑x_i³·x_i³ + b·∑x_i²·x_i³ + c·∑x_i·x_i³ + d·∑x_i³ = ∑y_i·x_i³

Natürlich lassen sich die Potenzen zusammenfassen, so daß sich ergibt:
a·∑x_i⁶ + b·∑x_i⁵ + c·∑x_i⁴ + d·∑x_i³ = ∑y_i·x_i³
bzw. in der vereinfachten Schreibweise mit Gaußschen Summenklammern:
a·[x⁶] + b·[x⁵] + c·[x⁴] + d·[x³] = [y·x³]

Für die übrigen partiellen Ableitungen ergibt sich analog:
a·[x⁵] + b·[x⁴] + c·[x³] + d·[x²] = [y·x²]
a·[x⁴] + b·[x³] + c·[x²] + d·[x] = [y·x]
a·[x³] + b·[x²] + c·[x] + d·[1] = [y]

Der Ausdruck [1] in der letzten Gleichung bedeutet ∑(1) für i=1 bis n, ist also gleich n, der Anzahl der Meßwerte bzw. der gegebenen Wertepaare.

(unfertige Version!)