[SciPy-User] How to fit a curve/function?

Wed Jun 8 10:56:25 EDT 2011

-------- Original-Nachricht --------
> Datum: Wed, 8 Jun 2011 10:33:45 -0400
> Von: josef.pktd at gmail.com
> An: SciPy Users List <scipy-user at scipy.org>
> Betreff: Re: [SciPy-User] How to fit a curve/function?

> On Wed, Jun 8, 2011 at 10:27 AM, Johannes Radinger <JRadinger at gmx.at>
> wrote:
> >
> > -------- Original-Nachricht --------
> >> Datum: Wed, 8 Jun 2011 10:12:58 -0400
> >> Von: josef.pktd at gmail.com
> >> An: SciPy Users List <scipy-user at scipy.org>
> >> Betreff: Re: [SciPy-User] How to fit a curve/function?
> >
> >> On Wed, Jun 8, 2011 at 9:41 AM, Johannes Radinger <JRadinger at gmx.at>
> >> wrote:
> >> >
> >> > -------- Original-Nachricht --------
> >> >> Datum: Wed, 8 Jun 2011 07:10:38 -0400
> >> >> Von: josef.pktd at gmail.com
> >> >> An: SciPy Users List <scipy-user at scipy.org>
> >> >> Betreff: Re: [SciPy-User] How to fit a curve/function?
> >> >
> >> >> On Wed, Jun 8, 2011 at 6:52 AM, Johannes Radinger <JRadinger at gmx.at>
> >> >> wrote:
> >> >> > Hello,
> >> >> >
> >> >> > I've got following function describing any kind of animal
> dispersal
> >> >> kernel:
> >> >> >
> >> >> > def pdf(x,s1,s2):
> >> >> >    return
> >> >>
> >>
> (p/(math.sqrt(2*math.pi*s1**2))*numpy.exp(-((x-0)**(2)/(2*s1**(2)))))+((1-p)/(s2*math.sqrt(2*math.pi))*numpy.exp(-((x-0)**(2)/(2*s2**(2)))))
> >> >> >
> >> >> > On the other hand I've got data from literature with which I want
> to
> >> fit
> >> >> the function so that I get s1, s2 and x.
> >> >> > Ususally the data in the literature are as follows:
> >> >> >
> >> >> > Example 1: 50% of the animals are between -270m and +270m and 90%
> >>  are
> >> >> between -500m and + 500m
> >> >> >
> >> >> > Example 2: 84% is between - 5000 m and +5000m, and 73% are between
> >> >> -1000m and +1000m
> >> >> >
> >> >> > So far as I understand an integration of the function is needed to
> >> solve
> >> >> for s1 and s2 as all the literature data give percentage (area under
> >> the
> >> >> curve) Can that be used to fit the curve or can that create ranges
> for
> >> s1
> >> >> and s2.
> >> >>
> >> >> I don't see a way around integration.
> >> >>
> >> >> If you have exactly 2 probabilities, then you can you a solver like
> >> >> scipy.optimize.fsolve to match the probabilites
> >> >> eg.
> >> >> 0.5 = integral pdf from -270 to 270
> >> >> 0.9 = integral pdf from -500 to 500
> >> >>
> >> >> If you have more than 2 probabilities, then using optimization of a
> >> >> weighted function of the moment conditions would be better.
> >> >>
> >> >> Josef
> >> >
> >> >
> >> >
> >> > Hello again
> >> >
> >> > I tried following, but without success so far. What do I have to do
> >> excactly...
> >> >
> >> > import numpy
> >> > from scipy import stats
> >> > from scipy import integrate
> >> > from scipy.optimize import fsolve
> >> > import math
> >> >
> >> > p=0.3
> >> >
> >> > def pdf(x,s1,s2):
> >> >    return
> >>
> (p/(math.sqrt(2*math.pi*s1**2))*numpy.exp(-((x-0)**(2)/(2*s1**(2)))))+((1-p)/(s2*math.sqrt(2*math.pi))*numpy.exp(-((x-0)**(2)/(2*s2**(2)))))
> >> >
> >> > def equ(s1,s2):
> >> >    0.5==integrate.quad(pdf,-270,270,args=(s1,s2))
> >> >    0.9==integrate.quad(pdf,-500,500,args=(s1,s2))
> >> >
> >> > result=fsolve(equ, 1,500)
> >> >
> >> > print result
> >>
> >> equ needs to return the deviation of the equations (I changed some
> >> details for s1 just to try it)
> >>
> >> import numpy
> >> from scipy import stats
> >> from scipy import integrate
> >> from scipy.optimize import fsolve
> >> import math
> >>
> >> p=0.3
> >>
> >> def pdf(x,s1,s2):
> >>     return
> >>
> (p/(math.sqrt(2*math.pi*s1**2))*numpy.exp(-((x-0)**(2)/(2*s1**(2)))))+((1-p)/(math.sqrt(2*math.pi*s2**2))*numpy.exp(-((x-0)**(2)/(2*s2**(2)))))
> >>
> >> def equ(arg):
> >>     s1,s2 = numpy.abs(arg)
> >>     cond1 = 0.5 - integrate.quad(pdf,-270,270,args=(s1,s2))[0]
> >>     cond2 = 0.9 - integrate.quad(pdf,-500,500,args=(s1,s2))[0]
> >>     return [cond1, cond2]
> >>
> >> result=fsolve(equ, [200., 1200])

thank you for your last reply...seems that the parameters of the two normals are nearly identical... anyway just two small addtional questions:

1)in fsolve(equ, [200., 1200]) the 200 and 1200 are kind of start values so far as I understand...how should these be choosen? what is recommended?

2) How can that be solve if I have I third condition (overfitted) can that be used as well or how does the alternative look like?

/johannes

> >>
> >> print result
> >>
> >> but in the results I get the parameters are very close to each other
> >> [-356.5283675   353.82544075]
> >>
> >> the pdf looks just like a mixture of 2 normals both with loc=0, then
> >> maybe the cdf of norm can be used directly
> >
> >
> > Thank you for that hint... First yes these are 2 superimposed normals
> but for other reasons I want to use the original formula instead of the
> stats.functions...
> >
> > anyway there is still a thing...the locator s1 and s2 are like the scale
> parameter of stats.norm so the are both + and -. For fsolve above it seems
> that I get only one parameter (s1 or s2) but for the positive and negative
> side of the distribution. So in actually there are four parameters -s1,
> +s1, -s2, +s2. How can I solve that? Maybe I can restrict the fsolve to look
> for the two values only in the positive range...
> 
> It doesn't really matter, if the scale only shows up in quadratic
> terms, or as in my initial change I added a absolute value, so whether
> it's positive or negative, it's still only one value, and we
> interprete it as postive scale
> 
> s1 = sqrt(s1**2)
> 
> Josef
> 
> >
> > any guesses?
> >
> > /J
> >
> >>
> >> >>> from scipy import stats
> >> >>> stats.norm.cdf(270, scale=350) - stats.norm.cdf(-270, scale=350)
> >> 0.55954705470577526
> >> >>>
> >> >>> stats.norm.cdf(270, scale=354) - stats.norm.cdf(-270, scale=354)
> >> 0.55436474670960978
> >> >>> stats.norm.cdf(500, scale=354) - stats.norm.cdf(-500, scale=354)
> >> 0.84217642881921018
> >>
> >> Josef
> >> >
> >> >
> >> > /Johannes
> >> >>
> >> >> >
> >> >> > /Johannes
> >> >> >
> >> >> > --
> >> >> > NEU: FreePhone - kostenlos mobil telefonieren!
> >> >> > Jetzt informieren: http://www.gmx.net/de/go/freephone
> >> >> > _______________________________________________
> >> >> > SciPy-User mailing list
> >> >> > SciPy-User at scipy.org
> >> >> > http://mail.scipy.org/mailman/listinfo/scipy-user
> >> >> >
> >> >> _______________________________________________
> >> >> SciPy-User mailing list
> >> >> SciPy-User at scipy.org
> >> >> http://mail.scipy.org/mailman/listinfo/scipy-user
> >> >
> >> > --
> >> > NEU: FreePhone - kostenlos mobil telefonieren!
> >> > Jetzt informieren: http://www.gmx.net/de/go/freephone
> >> > _______________________________________________
> >> > SciPy-User mailing list
> >> > SciPy-User at scipy.org
> >> > http://mail.scipy.org/mailman/listinfo/scipy-user
> >> >
> >> _______________________________________________
> >> SciPy-User mailing list
> >> SciPy-User at scipy.org
> >> http://mail.scipy.org/mailman/listinfo/scipy-user
> >
> > --
> > NEU: FreePhone - kostenlos mobil telefonieren!
> > Jetzt informieren: http://www.gmx.net/de/go/freephone
> > _______________________________________________
> > SciPy-User mailing list
> > SciPy-User at scipy.org
> > http://mail.scipy.org/mailman/listinfo/scipy-user
> >
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user

-- 
NEU: FreePhone - kostenlos mobil telefonieren!			
Jetzt informieren: http://www.gmx.net/de/go/freephone