[Numpy-discussion] gtkLeastSquares.py

Adrian Feiguin feiguin at magnet.fsu.edu
Mon Jul 2 11:03:52 CDT 2001


Hi everybody,
	I'm the author of SciGraphica
(http://scigraphica.sourceforge.net), an application for scientific
graphics and data analysis that uses Python and Numpy. Among other things,
it features spreadheets that can contain Python data, like formulas.
	The Python stuff was mainly done by Conrad Steenberg, because I
don't have much experience, specially embedding Python into C.
	However, I decided to get started this weekend because we needed
badly a GUI for non-linear least squares fits. The result is
gtkLeastSquares.py, which is attached below. It's going to be distributed
with SG, but it works independently and needs gtk.py and Scientific.
Basically it is a gtk GUI for Scientific.Functions.LeastSquares. I'm quite
satisfied with the result, considering that I'm a newbie ;-)
	I'm writting because I'd like to have some feedback, and people
interested in contributing. The code is still immature, and needs
enhancements, like user defined functions and parameters (either editing
the existing in the GUI, or adding new ones interactively). I still have
to learn.
	It is very simple, and easy to use. You'll find the code below,
with an example at the bottom of the file.
	I hope you like it. All comments, and suggestions are very welcome
(as well as contributors ;-)
	Thank you,
<ADRIAN>
-------------- Cut here ----------------------
#
# gtkLeastSquares.py : 	GUI for the module Scientific.Functions.LeastSquares
#		Implementation of the Levenberg-Marquardt algorithm for general
#		non-linear least-squares fits.
#
# Copyright (C) 2001 Adrian E. Feiguin <feiguin at ifir.edu.ar>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
#


from string import *
from gtk import *
from Numeric import *
from Scientific.Functions.LeastSquares import *
from Scientific.Functions.FirstDerivatives import *

def fit_linear(parameters, values):
	a, b = parameters
	x = values
	return (a + b * x)

def fit_cuadratic(parameters, values):
	a, b, c, x0 = parameters
	x = values
	return(a + b * (x - x0) + c * (x - x0)**2)
	return

def fit_gauss(parameters, values):
	y0, x0, a, w = parameters
	x = values
	return(y0 + a * exp(-2*(x-x0)**2/w**2))
	return

def fit_lorentz(parameters, values):
	x0, y0, a, b = parameters
	x = values
	return(y0 + 2*a/pi * w / (4 * (x - x0)**2 + w**2))
	return

def fit_boltzman(parameters, values):
	x0, a1, a2, dx = parameters
	x = values
	return((a1 - a2)/(1 + exp((x - x0)/dx)) + a2)

def fit_logistic(parameters, values):
	x0, a1, a2, p = parameters
	x = values
	return((a1 - a2)/(1 + (x/x0)**p) + a2)

def fit_expdecay(parameters, values):
	x0, y0, a, t = parameters
	x = values
	return(y0 + a * exp(-(x - x0)/t))

def fit_expgrow(parameters, values):
	x0, y0, a, t = parameters
	x = values
	return(y0 + a * exp((x - x0)/t))

def fit_expassoc(parameters, values):
	y0, a1, t1, a2, t2 = parameters
	x = values
	return(y0 + a1 * (1 + exp(-x/t1)) + a2 * (1 + exp(-x/t2)))

def fit_hyperbl(parameters, values):
	p1, p2 = parameters
	x = values
	return(p1 * x/ (p2 + x))

def fit_pulse(parameters, values):
	x0, y0, a, t1, t2 = parameters
	x = values
	return(y0 + a * (1 + exp(-(x - x0)/t1)) * exp(-(x - x0)/t2))

def fit_rational0(parameters, values):
	a, b, c = parameters
	x = values
	return((b + c*x)/(1 + a*x))

def fit_sine(parameters, values):
	x0, a, w = parameters
	x = values
	return(a * sin(pi*(x - x0)/w))

def fit_gaussamp(parameters, values):
	x0, y0, a, w = parameters
	x = values
	return(y0 + a * exp(-(x - x0)**2/(2*w**2)))

def fit_allometric(parameters, values):
	a, b = parameters
	x = values
	return(a * x**b)



fit_linear_dic = {
	"Doc" : "Linear Function",
	"Exp" : "y = a + b * x",
	"Par" :	("a", "b"),
	"NumPar" : 2,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_linear
}
fit_cuadratic_dic = {
	"Doc" : "Cuadratic Function",
	"Exp" : "y = a + b * (x - x0) + c * (x - x0)**2",
	"Par" :	("a", "b", "c", "x0"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_cuadratic
}
fit_gauss_dic = {
	"Doc" : "Amplitude version of Gaussian Function",
	"Exp" : "y = y0 + a * exp(-2*(x-x0)**2/w**2)",
	"Par" :	("x0", "y0", "a", "w"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_gauss
}
fit_lorentz_dic = {
	"Doc" : "Lorentzian Peak Function",
	"Exp" : "y = y0 + 2*a/pi * w / (4 * (x - x0)**2 + w**2)",
	"Par" :	("x0", "y0", "a", "w"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_lorentz
}
fit_boltzman_dic = {
	"Doc" : "Boltzman Function: sigmoidal curve",
	"Exp" : "y = (a1 - a2)/(1 + exp((x - x0)/dx)) + a2",
	"Par" :	("x0", "a1", "a2", "dx"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_boltzman
}
fit_logistic_dic = {
	"Doc" : "Logistic dose/response",
	"Exp" : "y = (a1 - a2)/(1 + (x/x0)**p) + a2",
	"Par" :	("x0", "a1", "a2", "p"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_logistic
}
fit_expdecay_dic = {
	"Doc" : "Exponential Decay",
	"Exp" : "y = y0 + a * exp(-(x - x0)/t)",
	"Par" :	("x0", "y0", "a", "t"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_expdecay
}
fit_expgrow_dic = {
	"Doc" : "Exponential Growth",
	"Exp" : "y = y0 + a * exp((x - x0)/t)",
	"Par" :	("x0", "y0", "a", "t"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_expgrow
}
fit_expassoc_dic = {
	"Doc" : "Exponential Associate",
	"Exp" : "y = y0 + a1 * (1 + exp(-x/t1)) + a2 * (1 + exp(-x/t2))",
	"Par" :	("y0", "a1", "t1", "a2", "t2"),
	"NumPar" : 5,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_expassoc
}
fit_hyperbl_dic = {
	"Doc" : "Hyperbola Function",
	"Exp" : "y = p1 * x/ (p2 + x)",
	"Par" :	("p1", "p2"),
	"NumPar" : 2,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_hyperbl
}
fit_pulse_dic = {
	"Doc" : "Pulse Function",
	"Exp" : "y = y0 + a * (1 + exp(-(x - x0)/t1)) * exp(-(x - x0)/t2)",
	"Par" :	("y0", "a", "t1", "t2"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_pulse
}
fit_rational0_dic = {
	"Doc" : "Rational Function , type 0",
	"Exp" : "y = (b + c*x)/(1 + a*x)",
	"Par" :	("a", "b", "c"),
	"NumPar" : 3,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_rational0
}
fit_sine_dic = {
	"Doc" : "Sine Function",
	"Exp" : "y = a * sin(pi*(x - x0)/w)",
	"Par" :	("a", "x0", "w"),
	"NumPar" : 3,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_sine
}
fit_gaussamp_dic = {
	"Doc" : "Amplitude version of Gaussian Peak Function",
	"Exp" : "y = y0 + a * exp(-(x - x0)**2/(2*w**2))",
	"Par" :	("x0", "y0", "a", "w"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_gaussamp
}
fit_allometric_dic = {
	"Doc" : "Classical Freundlich Model",
	"Exp" : "y = a * x**b",
	"Par" :	("a", "b"),
	"NumPar" : 4,
	"IVar" : "x",
	"DVar" : "y",
	"Function" : fit_allometric
}


functions = {
	"Linear" : fit_linear_dic,
	"Cuadratic" : fit_cuadratic_dic,
	"Gauss" : fit_gauss_dic,
	"Lorentz" : fit_lorentz_dic,
	"Boltzman" : fit_boltzman_dic,
	"Logistic" : fit_logistic_dic,
	"ExpDecay" : fit_expdecay_dic,
	"ExpGrow" : fit_expgrow_dic,
	"ExpAssoc" : fit_expassoc_dic,
	"Hyperbl" : fit_hyperbl_dic,
	"Pulse" : fit_pulse_dic,
	"Rational0" : fit_rational0_dic,
	"Sine" : fit_sine_dic,
	"GaussAmp" : fit_gaussamp_dic,
	"Allometric" : fit_allometric_dic
}

def fit(data):
	main_window = GtkWindow()
	main_window.set_title("Curve Fitting")
	main_window.set_border_width(5)
	main_window.connect("destroy", mainquit)
	main_window.connect("delete_event", mainquit)
	main_box = GtkVBox(FALSE, 5)
	main_window.add(main_box)
	main_frame = GtkFrame("Select Function")
	main_box.pack_start(main_frame)

	table = GtkTable(7, 4, FALSE)
	table.set_col_spacings(10)
	table.set_row_spacings(5)
	table.set_border_width(5)
	main_frame.add(table)

        swindow = GtkScrolledWindow()
        swindow.set_policy(POLICY_AUTOMATIC, POLICY_AUTOMATIC)
	swindow.set_usize(120, 100)
	table.attach(swindow, 0, 1, 0, 6)
        clist = GtkCList(1)
	swindow.add(clist)

	table.attach(GtkVSeparator(), 1, 2, 0, 6)

	text = map(lambda i: str(i), range(20))

	k = functions.keys()
	k.sort()
	for i in k:
		text[0] = i
		clist.append(text)

	label = GtkLabel("Exp:")
	label.set_alignment(1., .5)
	table.attach(label, 2, 3, 0, 1)
        fentry = GtkEntry()
#	fentry.set_editable(FALSE)
	table.attach(fentry, 3, 4, 0, 1)

	label = GtkLabel("Number of Param:")
	label.set_alignment(1., .5)
	table.attach(label, 2, 3, 1, 2)
        nspin = GtkSpinButton(GtkAdjustment(0, 0, 8, 1, 8, 0), 0, 0)
	nspin.set_editable(FALSE)
	nspin.set_state(STATE_INSENSITIVE)
	table.attach(nspin, 3, 4, 1, 2)

	label = GtkLabel("Param:")
	label.set_alignment(1., .5)
	table.attach(label, 2, 3, 2, 3)
        pentry = GtkEntry()
	pentry.set_editable(FALSE)
	pentry.set_state(STATE_INSENSITIVE)
	table.attach(pentry, 3, 4, 2, 3)

	label = GtkLabel("Independent Var:")
	label.set_alignment(1., .5)
	table.attach(label, 2, 3, 3, 4)
        iventry = GtkEntry()
	iventry.set_editable(FALSE)
	iventry.set_state(STATE_INSENSITIVE)
	table.attach(iventry, 3, 4, 3, 4)

	label = GtkLabel("Dependent Var:")
	label.set_alignment(1., .5)
	table.attach(label, 2, 3, 4, 5)
        dventry = GtkEntry()
	dventry.set_editable(FALSE)
	dventry.set_state(STATE_INSENSITIVE)
	table.attach(dventry, 3, 4, 4, 5)

	action_area = GtkHButtonBox()
	action_area.set_layout(BUTTONBOX_END)
	action_area.set_spacing(5)
	main_box.pack_start(action_area)
        fit_button = GtkButton("Fit")
	action_area.pack_start(fit_button)
        close_button = GtkButton("Close")
	action_area.pack_start(close_button)

	lframe = GtkFrame()
	lframe.set_shadow_type(SHADOW_IN)
	main_box.pack_start(lframe)
	explabel = GtkLabel("Choose a Fitting Function")
	lframe.add(explabel)



#	CALLBACK FUNCTIONS
	def select_function(_clist, row, col, event, functions = functions, label = explabel, fentry = fentry, pentry = pentry, nspin = nspin, iventry = iventry, dventry = dventry):
		k = _clist.get_text(row, col)
                f = functions[k]
		label.set_text(f["Doc"])
		fentry.set_text(f["Exp"])
		nspin.set_value(f["NumPar"])
		iventry.set_text(f["IVar"])
		dventry.set_text(f["DVar"])
		s = ""
		for i in f["Par"]:
			s = s + i + ", "
		pentry.set_text(s[:len(s)-2])

	def open_fit_dialog(_button, functions = functions, clist = clist, data = data):
		a = clist.__getattr__("selection")
		k = clist.get_text(a[0], 0)
		f = functions[k]
		param = (1, 1)
		fit_dialog(f, data)


# 	CONNECT OBJECTS

	clist.connect("select_row", select_function)

	fit_button.connect("clicked", open_fit_dialog)
	close_button.connect("clicked", main_window.destroy)

        clist.select_row(0, 0)
	main_window.show_all()
	mainloop()

def fit_dialog(f, data):
	main_window = GtkWindow()
	main_window.set_title("Fit")
	main_window.set_border_width(5)
	main_window.connect("destroy", mainquit)
	main_window.connect("delete_event", mainquit)
	table = GtkTable(len(f["Par"])+3, 2, FALSE)
	table.set_col_spacings(10)
	table.set_row_spacings(5)
	main_window.add(table)

	table.attach(GtkLabel("Variable"), 0, 1, 0, 1)
	table.attach(GtkLabel("Value"), 1, 2, 0, 1)
	table.attach(GtkHSeparator(), 0, 2, 1, 2)
	r = 2
	entries = []
	for i in f["Par"]:
#		check = GtkCheckButton(i+":")
#		table.attach(check, 0, 1, r, r+1)
		table.attach(GtkLabel(i+":"), 0, 1, r, r+1)
		entry = GtkEntry()
		entries = entries + [entry]
		entry.set_text("0.0")
		table.attach(entry, 1, 2, r, r+1)
		r = r + 1

	table.attach(GtkHSeparator(), 0, 2, r, r + 1)
	r = r + 1
	table.attach(GtkLabel("Chi_Sqr:"), 0, 1, r, r + 1)
	err_entry = GtkEntry()
	table.attach(err_entry, 1, 2, r, r + 1)
	r = r + 1
	table.attach(GtkHSeparator(), 0, 2, r, r + 1)

	def run_fit(_button, f = f, data = data, entries = entries, err_entry = err_entry):
		n = 0
		p = ()
		for i in entries:
			s = entries[n].get_text()
			p = p + (atof(s),)
			n = n + 1

		fit, error = leastSquaresFit(f["Function"], p, data)

		n = 0
		for i in entries:
			entries[n].set_text(str(fit[n]))
			n = n + 1

		err_entry.set_text(str(error))
#		print "Fitted parameters: ", fit
#		print "Fit error: ", error

		return


	action_area = GtkHButtonBox()
	action_area.set_layout(BUTTONBOX_SPREAD)
	action_area.set_spacing(5)
	run_button = GtkButton("Run")
	close_button = GtkButton("Close")
	action_area.pack_start(run_button)
	action_area.pack_start(close_button)
	table.attach(action_area, 0, 2, r + 1, r + 2)

#	CONNECT OBJECTS

	run_button.connect("clicked", run_fit)
	close_button.connect("clicked", main_window.destroy)

	main_window.show_all()
	mainloop()


# Test for linear fit:
#
# from gtkLeastSquares import *
# data = [ (0., 0.), (1., 1.1), (2., 1.98), (3., 3.05) ]
# fit(data)






More information about the Numpy-discussion mailing list