3816.sagews

import numpy as np

data=np.loadtxt(DATA+'data10.txt')

data

array([[ 0.00000000e+00, 1.06303624e-07], [ 1.00000000e+00, 2.78865755e-07], [ 2.00000000e+00, 5.43302390e-07], [ 3.00000000e+00, 1.12037774e-06], [ 4.00000000e+00, 2.23183048e-06], [ 5.00000000e+00, 5.39501842e-06], [ 6.00000000e+00, 1.05425434e-05], [ 7.00000000e+00, 1.83262304e-05], [ 8.00000000e+00, 3.35136176e-05], [ 9.00000000e+00, 6.00607833e-05], [ 1.00000000e+01, 1.21112469e-04], [ 1.10000000e+01, 2.05524920e-04], [ 1.20000000e+01, 3.93775246e-04], [ 1.30000000e+01, 6.98180160e-04], [ 1.40000000e+01, 1.08503928e-03], [ 1.50000000e+01, 2.06738256e-03], [ 1.60000000e+01, 3.26789359e-03], [ 1.70000000e+01, 5.35184608e-03], [ 1.80000000e+01, 9.24101365e-03], [ 1.90000000e+01, 1.25659508e-02], [ 2.00000000e+01, 1.93323342e-02], [ 2.10000000e+01, 2.64702201e-02], [ 2.20000000e+01, 4.08125355e-02], [ 2.30000000e+01, 5.81097319e-02], [ 2.40000000e+01, 9.44573011e-02], [ 2.50000000e+01, 9.62788915e-02], [ 2.60000000e+01, 1.40389823e-01], [ 2.70000000e+01, 2.00509760e-01], [ 2.80000000e+01, 2.22924049e-01], [ 2.90000000e+01, 3.04465393e-01], [ 3.00000000e+01, 3.55468980e-01], [ 3.10000000e+01, 4.48517768e-01], [ 3.20000000e+01, 5.95769069e-01], [ 3.30000000e+01, 6.19332979e-01], [ 3.40000000e+01, 7.95076984e-01], [ 3.50000000e+01, 6.45010715e-01], [ 3.60000000e+01, 8.81946795e-01], [ 3.70000000e+01, 9.56774135e-01], [ 3.80000000e+01, 1.07420930e+00], [ 3.90000000e+01, 1.12315905e+00], [ 4.00000000e+01, 8.49608389e-01], [ 4.10000000e+01, 8.74969831e-01], [ 4.20000000e+01, 9.36188001e-01], [ 4.30000000e+01, 9.95437847e-01], [ 4.40000000e+01, 1.02838136e+00], [ 4.50000000e+01, 7.45478522e-01], [ 4.60000000e+01, 7.34565760e-01], [ 4.70000000e+01, 7.01208185e-01], [ 4.80000000e+01, 5.17418172e-01], [ 4.90000000e+01, 5.18380873e-01], [ 5.00000000e+01, 3.28983872e-01], [ 5.10000000e+01, 2.75376214e-01], [ 5.20000000e+01, 2.70140279e-01], [ 5.30000000e+01, 1.98010599e-01], [ 5.40000000e+01, 1.24825434e-01], [ 5.50000000e+01, 9.99630383e-02], [ 5.60000000e+01, 7.48163653e-02], [ 5.70000000e+01, 4.19234766e-02], [ 5.80000000e+01, 3.43473013e-02], [ 5.90000000e+01, 2.46199534e-02], [ 6.00000000e+01, 1.83980407e-02], [ 6.10000000e+01, 1.08684577e-02], [ 6.20000000e+01, 9.03862746e-03], [ 6.30000000e+01, 5.36881679e-03], [ 6.40000000e+01, 3.30836161e-03], [ 6.50000000e+01, 1.75524197e-03], [ 6.60000000e+01, 1.13850872e-03], [ 6.70000000e+01, 7.64783751e-04], [ 6.80000000e+01, 4.00069374e-04], [ 6.90000000e+01, 1.76444304e-04], [ 7.00000000e+01, 1.39110285e-04], [ 7.10000000e+01, 5.76224429e-05], [ 7.20000000e+01, 3.55358016e-05], [ 7.30000000e+01, 1.50822422e-05], [ 7.40000000e+01, 9.85036345e-06], [ 7.50000000e+01, 4.83466645e-06], [ 7.60000000e+01, 2.32112152e-06], [ 7.70000000e+01, 1.02468328e-06], [ 7.80000000e+01, 5.73216367e-07], [ 7.90000000e+01, 2.35057778e-07], [ 8.00000000e+01, 1.23921353e-07], [ 8.10000000e+01, 4.62800773e-08], [ 8.20000000e+01, 2.20658347e-08], [ 8.30000000e+01, 8.77742622e-09], [ 8.40000000e+01, 3.05402531e-09], [ 8.50000000e+01, 2.05262014e-09], [ 8.60000000e+01, 5.89336258e-10], [ 8.70000000e+01, 2.06633675e-10], [ 8.80000000e+01, 1.14723040e-10], [ 8.90000000e+01, 4.75526804e-11], [ 9.00000000e+01, 1.59046236e-11], [ 9.10000000e+01, 6.16295746e-12], [ 9.20000000e+01, 1.72555772e-12], [ 9.30000000e+01, 6.08521100e-13], [ 9.40000000e+01, 2.28635426e-13], [ 9.50000000e+01, 6.53443947e-14], [ 9.60000000e+01, 2.10213937e-14], [ 9.70000000e+01, 8.40599566e-15], [ 9.80000000e+01, 2.70328686e-15], [ 9.90000000e+01, 6.93012174e-16]])

list_plot([(data[i,0],data[i,1]) for i in range(100)])

import numpy as np

data=np.loadtxt(DATA+'data10.txt')

dataplot=list_plot([(data[i,0],data[i,1]) for i in range(100)])

var("a b c d")
model1(x)=a*exp(-(x-b)**2/(c)**2)
model2(x)=a*exp(-(x-b)**2/(c)**2)+d

param1=find_fit(data,model1,solution_dict=True)
param2=find_fit(data,model2,solution_dict=True,parameters=[a,b,c,d],variables=[x],initial_quess=[1,40,10,0])

Traceback (most recent call last): File "<stdin>", line 1, in <module> File "_sage_input_175.py", line 10, in <module> exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("cGFyYW0xPWZpbmRfZml0KGRhdGEsbW9kZWwxLHNvbHV0aW9uX2RpY3Q9VHJ1ZSkKcGFyYW0yPWZpbmRfZml0KGRhdGEsbW9kZWwyLHNvbHV0aW9uX2RpY3Q9VHJ1ZSxwYXJhbWV0ZXJzPVthLGIsYyxkXSx2YXJpYWJsZXM9W3hdLGluaXRpYWxfcXVlc3M9WzEsNDAsMTAsMF0p"),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))' + '\n', '', 'single') File "", line 1, in <module> File "/tmp/tmpBDZQtI/___code___.py", line 4, in <module> exec compile(u'param2=find_fit(data,model2,solution_dict=True,parameters=[a,b,c,d],variables=[x],initial_quess=[_sage_const_1 ,_sage_const_40 ,_sage_const_10 ,_sage_const_0 ])' + '\n', '', 'single') File "", line 1, in <module> TypeError: find_fit() got an unexpected keyword argument 'initial_quess'

plotmodel1=plot(model1(x).subs_expr((param1)),(x,0,100),color='red')
plotmodel2=plot(model2(x).subs_expr((param2)),(x,0,100),color='purple')
show(plotmodel1+plotmodel2+dataplot)

import scipy.optimize as so

so.leastsq

<function leastsq at 0xd1aff7c>

###TD 2.estimation approximation...............

### calcule de a,b,c ? et" 1 sinefie qu'il a reusie"

import scipy.optimize as so

def erreur(p,y,x):
    sigma=(0.05/1.96)*y
    a,b,c=p
    #a,b,c,d=p
    model=a*exp(-(x-b)**2/(c)**2)
    #model=a*exp(-(x-b**2/(c)**2)+d
    err=(y-model)/sigma
    return(err)

def FitChi2(p0,y,x):
    p=so.leastsq(erreur,p0,(y,x))
    return (p)

p0=(1,40,10)
y=data[:,1]
x=data[:,0]
FitChi2(p0,y,x)

(array([ 0.98361202, 39.96540068, 10.00427897]), 1)

### calcule de a,b,c et x**2?

import scipy.optimize as so

def model1(x,a,b,c):
    return(a*exp(-(x-b)**2/(c)**2))
def model2(x,a,b,c,d):
    return(a*exp(-(x-b)**2/(c)**2)+d)

def erreur(p,y,x):
    sigma=(0.05/1.96)*y
    a,b,c=p
    #a,b,c,d=p
    model=model1(x,a,b,c)
    #model=model2(x,a,b,c,d)
    err=(y-model)/sigma
    return(err)

def FitChi2(p0,y,x):
    p,v=so.leastsq(erreur,p0,(y,x))
    return (p)

def Chi2(p,y,x):
    a,b,c=p
    SIGMA=(0.05/1.96)*y
    #a,b,c,d=p
    Chi2=((y-model1(x,a,b,c))/SIGMA)**2
    #Chi2=((y-model2(x,a,b,c,d))/SIGMA)**2
    return(Chi2.sum())

p0=(1,40,10)
y=data[:,1]
x=data[:,0]
p=FitChi2(p0,y,x)
print(p)
Chi2(p,y,x)

[ 0.98361202 39.96540068 10.00427897] 1786.6314531010312

##Chi2=1786.6314531010312

erf?

File: /usr/local/sage/local/lib/python2.6/site-packages/sage/functions/other.py

Type: <class ‘sage.functions.other.Function_erf’>

Definition: erf(*args, coerce=True, hold=False)

Docstring:

The error function, defined as erf(x) = \frac{2}{\sqrt{\pi}} \int_0^x e^{-t^2} dt.

Sage currently only implements the error function (via a call to PARI) when the input is real.

EXAMPLES:

sage: erf(2)
erf(2)
sage: erf(2).n()
0.995322265018953
sage: loads(dumps(erf))
erf

The following fails because we haven’t implemented erf yet for complex values:

sage: complex(erf(3*I))
Traceback (most recent call last):
...
TypeError: unable to simplify to complex approximation

TESTS:

Check if conversion from maxima elements work:

sage: merf = maxima(erf(x)).sage().operator()
sage: merf == erf
True