Simple Fit

lua-users home
wiki

Showing revision 8
Here is a way to fit curves--e.g. straight lines, parabolas, and exponential functions--using LuaMatrix.

Files:wiki_insecure/users/chill/fit-0.1.zip

The code and method is very simple.

First you get the x values into a table. Then you concatenate it with the y-values. Use the Gauss-Jordan Method to get the results for your variables.

The only thing one had to think about with exponential functions was how to make them linear. That is easy; for example:

y = a * x^b | ln ==> ln(y) = ln( a ) + b * ln( x )

Then one can use fit:linear() again to get the variables a and b.

--///////////////--
-- Curve Fitting --
--///////////////--

-- v 0.1

-- Lua 5.1 compatible

-- chillcode, http://lua-users.org/wiki/SimpleFit
-- Licensed under the same terms as Lua itself.


-- requires matrix module
local matrix = require "matrix"

-- The Fit Table
local fit = {}

-- Note all these Algos use the Gauß-Jordan Method to caculate equation systems

-- Fit a straight Line
-- model (  y = a + b * x  )
function fit:linear( x_values,y_values )
	-- x_values = { x1,x2,x3,...,xn }
	-- y_values = { y1,y2,y3,...,yn }
	
	-- values for A matrix
	local a_vals = {}
	-- values for Y vector
	local y_vals = {}

	for i,v in ipairs( x_values ) do
		a_vals[i] = { 1, v }
		y_vals[i] = { y_values[i] }
	end

	-- create both Matrixes
	local A = matrix:new( a_vals )
	local Y = matrix:new( y_vals )

	local ATA = matrix.mul( matrix.transpose(A), A )
	local ATY = matrix.mul( matrix.transpose(A), Y )

	local ATAATY = matrix.concath(ATA,ATY)

	--print( "=>    y = ( "..a.." )  +  ( "..b.." ) * x")
	-- return a,b
	return self:getresults( ATAATY )
end

-- Fit a Parabola
-- model (  y = a + b * x + c * x² )
function fit:parabolic( x_values,y_values )
	-- x_values = { x1,x2,x3,...,xn }
	-- y_values = { y1,y2,y3,...,yn }

	-- values for A matrix
	local a_vals = {}
	-- values for Y vector
	local y_vals = {}

	for i,v in ipairs( x_values ) do
		a_vals[i] = { 1, v, v*v }
		y_vals[i] = { y_values[i] }
	end

	-- create both Matrixes
	-- create both Matrixes
	local A = matrix:new( a_vals )
	local Y = matrix:new( y_vals )

	local ATA = matrix.mul( matrix.transpose(A), A )
	local ATY = matrix.mul( matrix.transpose(A), Y )

	local ATAATY = matrix.concath(ATA,ATY)

	--print( "=>    y = ( "..a.." )  +  ( "..b.." ) * x  +  ( "..c.." ) * x²")
	-- return a,b,c
	return self:getresults( ATAATY )
end

-- Fit Exponential
-- model (  y = a * x^b )
function fit:exponential( x_values,y_values )
	-- convert to linear problem
	-- ln(y) = ln(a) + b * ln(x)
	for i,v in ipairs( x_values ) do
		x_values[i] = math.log( v )
		y_values[i] = math.log( y_values[i] )
	end

	local a,b = self:linear( x_values,y_values )

	--print( "=>    y = ( "..a.." )  *  x^( "..b.." )")
	return math.exp(a), b
end

-- get Results
function fit:getresults( mtx )
	assert( #mtx+1 == #mtx[1], "Cannot calculate Results" )
	mtx:dogauss()
	-- tresults
	local cols = #mtx[1]
	local tres = {}
	for i = 1,#mtx do
		tres[i] = mtx[i][cols]
	end
	return unpack( tres )
end
return fit

Testcode in the *.zip file


RecentChanges · preferences
edit · history · current revision
Edited August 11, 2007 4:54 pm GMT (diff)