Gegeben: Kennfeld mit Eingangsmatrix X und einem Ergebnisvektor y
x1 x2 y
----------------
1,2 3,1 5,7
2,5 3,1 8,2
3,5 4,5 5,0
4,0 4,5 8,2
6,0 5,0 9,5
7,5 6,5 10,5
Gesucht: Lösungsvektor b der Flächengleichung y = f(x1,x2) mit beliebiger Näherung
y = b0 + b1x1 + b2x1^2 + b3x2 + b4x2^2
Lösungsprinzip
Lösungsgleichung über Matrizenrechnung (Ordinary Least Squares, OLS)
b = (X^T * X)^-1 * X^T * y
Erweiterte Tabelle mit Eingangswerten, die fehlenden Werte für x1² und x2² werden ergänzt
x1 x1² x2 x2²
x1 x2 x3 x4 y
---------------------------------
1,2 1,44 3,1 9,61 5,7
2,5 6,25 3,1 9,61 8,2
3,5 12,25 4,5 20,25 5,0
4,0 16 4,5 20,25 8,2
6,0 36 5,0 25 9,5
7,5 56,25 6,5 42,25 10,5
SciLab-Code
// Eingabedaten
x1 = [1.2 2.5 3.5 4.0 6.0 7.5]
x2 = [1.44 6.25 12.25 16.0 36.0 56.25]
x3 = [3.1 3.1 4.5 4.5 5.0 6.5]
x4 = [9.61 9.61 20.25 20.25 25.0 42.25]
y = [5.7 8.2 5.0 8.2 9.5 10.5]
// Auffüllen der ersten Spalte der X-Matrix mit 1 ist wichtig (Startwerte)!
X = [ones(6,1) x1' x2' x3' x4']
// Koeffizientenvektor ausrechnen (y muss auch transponiert werden, da Eingabe als Zeilenvektor)
b = inv(X'*X)*X'*y'
// Koeffizientenvektor testen
erg = (X*b)'
Ergebnis
b = [19.6184110 2.3267959 -0.0386958 -6.9383043 0.4893346]
erg = [5.548606304 8.387314176 5.9748294 6.9931181 9.7279811 10.4681504]
Grafische Anzeige
xdel
clf()
// Laufvariable für x1 von 1.2 bis 7.5, Schrittweite 0.1
x = [1.2:0.1:7.5];
// Laufvariable für x1 von 1.2 bis 7.5, Schrittweite 0.1
y = [3.1:0.1:6.5];
// Funktionsdefinition
function z=f(x,y)
// ausgerechnete Werte für b verwenden
z = 19.6184110 + 2.3267959 * x + -0.0386958 * x^2 + -6.9383043 * y + 0.4893346 * y^2;
endfunction
// 3D-Plot anzeigen
fplot3d(x, y, f);
Ausblick
Das Verfahren kann auch für Gleichungen mit Wurzelfunktionen oder trigonometrischen Termen benutzt werden. Je nach Beschaffenheit des “Gebirges” und gewünschtem Rechenaufwand können bel. Terme hinzugefügt werden. Es wäre z.B. also auch vorstellbar die Koeffizienten folgender Flächengleichung zu approximieren:
y = b0 + b1x1 + b2sin(x1) + b3x2 + b4sin(x2) + b5x1x2
Entsprechend ist die Tabelle mit Eingangswerten anzupassen. Je mehr Terme in der Näherungsgleichung addiert werden, um so genauer wird das Ergebnis in Abhängigkeit der verwendeten Glieder. Gleichzeitig steigt der Rechenaufwand bei der Matrizeninversion.
Links
deutsche Wikipedia
englische Wikipedia
Egwald Statistics: Multiple Regression
Michael Palmer: MULTIPLE REGRESSION