Interpolation

Numerical analysis technique to fit a function to data points.

by 50kgTorus

Author Avatar

Introduction Today, I will show you mathematically, how one would fit a function type to known data, then I will give a code.

Prerequisites: Knowledge of algebra, some linear algebra(though I will go over the linear algebra given that this is a game meant for children)

Functions If you recall from algebra, there are these two kinds of functions: linear or nonlinear, which you can express as a general polynomial function. For instance:

(Eq 1) f(x) = (mx+c)^n -- where ^ indicates exponent and m and c are real numbers and n is a positive integer

From Eq 1, you can see that if you set n = 1, you will get the commonly known linear function:

(Eq 2) f(x) = mx + c

And if n = 2, then:

(Eq 3) f(x) = (mx)^2 + 2cmx + c^2 -- Resembles a quadratic function format!

Knowing these common functions formats will be vital to later instruction. I didn't give much application, though, this entire thread will be using this idea, so I think it suffices.

System of Equations At times, you may come across what I like to think of as a set of knowns. For instance:

(Eqs 4) a + b = 3 2a+4b = 5

If you conduct 3-5, this is the equivalence of a+b - (2a+4b), and you can also multiply single equations by scalar values(i.e. 2*3 - 5 = 2a+2b - (2a + 4b)). In this instance, we are interested in conducting an operation among our system of equations such that we create a new equation dropping the number of unknowns by 1.

In our example of 2*3 - 5 = 2a+2b - (2a +4b), we can reduce this to 1 = -2b, and we get that b = -1/2. We now know a value of b that satisfies our conditions. Let's find a. We know b, so we can use one of our equations from (Eqs 4) to find a:

a + (-1/2) = 3 --plugging -1/2 in for b a = 3 + 1/2 = 7/2 --addings +1/2 to both sides to isolate a

To test with our other equation from (Eqs 4): 2(7/2) + 4(-1/2) = 5 -- Using the second equation to check our a and b 7 - 2 = 5 5 = 5 -- check

Next, we will be discussing linear algebra, which basically uses this same idea.

Linear Algebra In this section, I will show how we can represent a system of equations as a matrix. I will be using LUA Vector2 and Vector3 objects to better explain this process because this object has desired properties. These vector objects will be used to indicate rows of a matrix. For instance:

a + b = 3 => A*{a,b} = B -- Where A = {Vector2.new(1,1),Vector2.new(2,4)} and B = {3,5} 2a + 4b = 5

A is the matrix of coefficients for this system of equations, and B is the conditions of rows 1 and 2. For instance, given the coefficients of row 1, those coefficients shall give 3. For our application, understanding this form is not so important, but if the reader's interest is piqued, then he or she should look up matrix multiplication. We shall strip away A and B from this form and create a new matrix C.

C = {Vector3.new(1,1,3),Vector3.new(2,4,5)}

And we wish to get matrix C into the following form called Row Echelon Form or REF:

{Vector3.new(1,g,f),Vector3.new(0,1,b)} -- Indicating the answer to b will be the Z component of row 2.

To get this desired result, we must manipulate the C matrix. Row 1 already has our desired conditions, so we will leave it alone. Though, if the X component of row 1 wasn't 1, we would divide the first row by its X component.

C[0] = C[0]/C[0].X --The convenience of using vectors.

After doing so, we need to get a 0 in row 2's X component, so similar to our system of equations section, we will do the following:

C[1] = C[1] - C[1].X*C[0] -- Subtracting the X component times the first row from the second row C = {Vector3.new(1,1,3),Vector3.new(0,2,-1)} -- Result of this operation

To get a 1 in the Y component of the second row, we shall divide the second row by that value.

C[1] = C[1]/C[1].Y C = {Vector3.new(1,1,3),Vector3.new(0,1,-1/2)}

Equivalent system of equations: a + b = 3 -> a = 3 - b = 3 + 1/2 = 7/2 0*a + b = -1/2 -> b = -1/2

The reason we use this method of solving this system of equations is because this process provides a systematic approach to finding the solutions to our unknown variables when there are more than 2 unknowns.

Interpolation We have all the necessary components to teach interpolation. To begin, we will use a more generalized version of (Eq 2) for what is known as quadratic interpolation:

(Eq 5) f(x) = (mx)^2 + 2cmx + c^2 => ax^2 + bx + c

Suppose we know 3 points to which we wish to fit this function: (x,y) belonging to {(1,1),(3,2),(4,3)}

(Eqs 6) f(1) = a(1^2) + b(1) + c = 1 <br/> f(3) = a(3^2 )+ b(3) + c = 2 <br/> f(4) = a(4^2) + b(4) + c = 3

=> A = {Vector3.new(1^2,1,1),Vector3.new(3^2,3,2),Vector3.new(4^2,4,3)} B = Vector3.new(1,2,3)

Like the linear algebra section, we wish to manipulate matrix A such that we have:

A = {Vector3.new(1,g,f),Vector3.new(0,1,j),Vector3.new(0,0,1)} B = Vector3.new(k,m,c) -- Indicating that once we get A and B into REF form, B will have the answer to c in B's Z component and then a and b can be solved by

a = k - gb - fc b = m - jc

and you will have a new coefficient vector V = Vector3.new(a,b,c)

Once you have your coefficients known, you can approximate new points on the interval of the 3 points you interpolated with. For instance, with (Eqs 6), your domain would be points from [1,4]. You will approximate this point with this form:

f(x') = V:Dot(Vector3.new((x')^2,x',1)) or f(x') = V.X(x')^2 + V.Y(x') + V.Z -- The dot product version looks nicer for me.

Exercise for the reader: Write an algorithm that performs quadratic interpolation. Below, I have a solution; however, try writing the algorithm yourself before looking at my solution.

Solution:

local A = {Vector3.new(math.pow(x0,2),x0,1), 
Vector3.new(math.pow(x1,2),x1,1),
Vector3.new(math.pow(x2,2),x2,1)}

        local B = {y0,y1,y2}

        local c = A[3].X
        A[3] = A[3]/c
        B[3] = B[3]/c
        c = A[2].X
        A[2] = A[2] - c*A[3]
        B[2] = B[2] - c*B[3]
        c = A[1].X
        A[1] = A[1] - c*A[3]
        B[1] = B[1] - c*B[3]
        c = A[2].Y
        A[2] = A[2]/c
        B[2] = B[2]/c
        c = A[1].Y
        A[1] = A[1] - c*A[2]
        B[1] = B[1] - c*B[2]
        c = A[1].Z
        A[1] = A[1]/c
        B[1] = B[1]/c

        local coeff = Vector3.new(B[3] - A[3].Y*(B[2]-A[2].Z*B[1]) - A[3].Z*B[1],
                 B[2]-A[2].Z*B[1],
                 B[1])

 local f  = coeff:Dot(Vector3.new(math.pow(t,2),t,1))

Example of Use

img|80x45

View in-game to comment, award, and more!