An overdetermined system can take on many forms. The basic idea behind over determined systems is that there is more data than there are variables to be solved. Many times, it's possible to use a least squares method to assist with solving an over determined system. (Click here for more information on overdetermined systems and the least squares methodology.) Mapack uses QR decomposition to solve the over determined system.
I've been working on an overdetermined system of the form:
(a1⋅x + b1)⋅(a2⋅y + b2) = f(x,y)
which resolves to the following:
a1⋅a2⋅x⋅y + a1⋅b2⋅x + a2⋅b1⋅y + b1⋅b2 = f(x,y)
the coefficients resolve to the following for my purposes:
a'1⋅x⋅y + a'2⋅x + a'3⋅y + a'4 = f(x,y).
(All I've done is combine the coefficients, since I am not interested in solving for the original equation. I am only interested in a function that fits my data smoothly.)
The system of equations I want to solve for the variables a'1..a'4 is as follows:
a'1 ⋅ 65775 + a'2 ⋅ 75 + a'3 ⋅ 877 + a'4 = 2949
a'1 ⋅ 75750 + a'2 ⋅ 75 + a'3 ⋅ 1010 + a'4 = 4086
a'1 ⋅ 78675 + a'2 ⋅ 75 + a'3 ⋅ 1049 + a'4 = 4420
a'1 ⋅ 80775 + a'2 ⋅ 75 + a'3 ⋅ 1077 + a'4 = 4659
a'1 ⋅ 51743 + a'2 ⋅ 50 + a'3 ⋅ 877 + a'4 = 3046
a'1 ⋅ 59590 + a'2 ⋅ 50 + a'3 ⋅ 1010 + a'4 = 4111
a'1 ⋅ 61891 + a'2 ⋅ 50 + a'3 ⋅ 1049 + a'4 = 4419
a'1 ⋅ 63543 + a'2 ⋅ 50 + a'3 ⋅ 1077 + a'4 = 4643
a'1 ⋅ 21925 + a'2 ⋅ 25 + a'3 ⋅ 877 + a'4 = 3120
a'1 ⋅ 25250 + a'2 ⋅ 25 + a'3 ⋅ 1010 + a'4 = 4126
a'1 ⋅ 26225 + a'2 ⋅ 25 + a'3 ⋅ 1049 + a'4 = 4420
a'1 ⋅ 26925 + a'2 ⋅ 25 + a'3 ⋅ 1077 + a'4 = 4632
I have four unknowns and 12 equations. This problem is well suited for the QR decomp algorithim. Here's the c# code for the problem:
using System;
using Mapack;
class Example
{
public static void Main(String[] args)
{
// The extra '1' at the end of the array is a constant placeholder
double[,] Data = new double[12, 4]
{
// Coefficient array
// {x, y, x ⋅ y, 1}
{75, 877, 65775, 1},
{75, 1010, 75750, 1},
{75, 1049, 78675, 1},
{75, 1077, 80775, 1},
{50, 877, 51743, 1},
{50, 1010, 59590, 1},
{50, 1049, 61891, 1},
{50, 1077, 63543, 1},
{25, 877, 21925, 1},
{25, 1010, 25250, 1},
{25, 1049, 26225, 1},
{25, 1077, 26925, 1}
};
Matrix A = new Matrix(Data);
// f(x,y) values
double[] Y = new double[12]
{
2949,
4086,
4420,
4659,
3046,
4111,
4419,
4643,
3120,
4126,
4420,
4632
};
Matrix B = new Matrix(Y);
Matrix X = A.Solve(B);
Console.WriteLine("A.Solve(B)");
Console.WriteLine(X.ToString());
Console.WriteLine("B = ");
Console.WriteLine(B.ToString());
Matrix T = A * X;
Console.WriteLine("A * A.Solve(B) = ");
Console.WriteLine(T.ToString());
Matrix R = B - T;
Console.WriteLine("Residuals = ");
Console.WriteLine(R.ToString());
int Row = 0;
double SumSq = 0.0;
double MaxRes = 0.0;
for (int i = 0; i < R.Rows; i++)
{
SumSq += R[i, 0] * R[i, 0];
if (Math.Abs(R[i, 0]) > Math.Abs(MaxRes))
{
MaxRes = R[i, 0];
Row = i;
}
}
SumSq = Math.Sqrt(SumSq / R.Rows);
Console.Write("RMS Residuals = ");
Console.WriteLine(SumSq);
Console.WriteLine("Max Res = {0} Row = {1}", MaxRes, Row);
}
}
The output of this program is shown below. In this case, the fit is not very good, but the output is a good example of what to expect. (I actually had to expand the terms to a much higher order of x and a quadratic in y to achieve an acceptable fit.)
A.Solve(B)
-3.55322662491747
7.89265922972302
0.0026246963617418
-3827.62671876104
B =
2949
4086
4420
4659
3046
4111
4419
4643
3120
4126
4420
4632
A * A.Solve(B) = B =
3000.38283203081
4076.28785579234
4391.77880260963
4618.28512340154
3052.38375830578
4122.70342820953
4436.5565644971
4661.88702131894
3062.9512278143
4121.40202077025
4431.77480968215
4654.60655556761
Residuals =
-51.3828320308057
9.71214420765864
28.2211973903659
40.7148765984648
-6.38375830578025
-11.7034282095301
-17.5565644970966
-18.8870213189375
57.0487721856989
4.59797922974576
-11.7748096821497
-22.6065555676132
RMS Residuals = 28.8050050381453
Max Res = 57.0487721856989 Row = 8
The final form of the equation is:
-3.55322662491747 x +
7.89265922972302 ⋅ y + 0.0026246963617418 ⋅ x ⋅ y - 3827.62671876104 = f(x,y)
You should now be able to approach your own solutions using mapack and QR decomposition.