Ajuste de la función por partes en Python

Estoy tratando de ajustar una función definida por partes a un conjunto de datos en Python. He buscado durante bastante tiempo ahora, pero no he encontrado una respuesta si es posible o no.

Para obtener una impresión de lo que estoy tratando de hacer, mire el siguiente ejemplo (que no funciona para mí). Aquí estoy tratando de ajustar una función de valor absoluto desplazado (f (x) = | xp |) a un conjunto de datos con p como el parámetro de ajuste.

import scipy.optimize as so import numpy as np def fitfunc(x,p): if x>p: return xp else: return -(xp) fitfunc = np.vectorize(fitfunc) #vectorize so you can use func with array x=np.arange(1,10) y=fitfunc(x,6)+0.1*np.random.randn(len(x)) popt, pcov = so.curve_fit(fitfunc, x, y) #fitting routine that gives error 

¿Hay alguna forma de lograr esto en Python?

Una forma de hacer esto en R es:

 # Fit of a absolute value function f(x)=|xp| f.lr p, xp,-(xp)) } x <- seq(0,10) # y <- f.lr(x,6) + rnorm (length(x),0,2) plot(y ~ x) fit.lr <- nls(y ~ f.lr(x,p), start = list(p = 0), trace = T, control = list(warnOnly = T,minFactor = 1/2048)) summary(fit.lr) coefficients(fit.lr) p.fit <- coefficients(fit.lr)["p"] x_fine <- seq(0,10,length.out=1000) lines(x_fine,f.lr(x_fine,p.fit),type='l',col='red') lines(x,f.lr(x,6),type='l',col='blue') 

Después de más investigación encontré una manera de hacerlo. En esta solución, no me gusta el hecho de que yo mismo deba definir la función de error. Además, no estoy realmente seguro de por qué tiene que ser en este estilo lambda. Por lo tanto, cualquier tipo de sugerencias o soluciones más sofisticadas son muy bienvenidas.

     import scipy.optimize as so import numpy as np import matplotlib.pyplot as plt def fitfunc(p,x): return x - p if x > p else p - x def array_fitfunc(p,x): y = np.zeros(x.shape) for i in range(len(y)): y[i]=fitfunc(x[i],p) return y errfunc = lambda p, x, y: array_fitfunc(p, x) - y # Distance to the target function x=np.arange(1,10) x_fine=np.arange(1,10,0.1) y=array_fitfunc(6,x)+1*np.random.randn(len(x)) #data with noise p1, success = so.leastsq(errfunc, -100, args=(x, y), epsfcn=1.) # -100 is the initial value for p; epsfcn sets the step width plt.plot(x,y,'o') # fit data plt.plot(x_fine,array_fitfunc(6,x_fine),'r-') #original function plt.plot(x_fine,array_fitfunc(p1[0],x_fine),'b-') #fitted version plt.show() 

    Para terminar esto aquí, compartiré mi propia solución final al problema. Para permanecer cerca de mi pregunta original, solo tiene que definir la función vectorizada y no usar np.vectorize .

     import scipy.optimize as so import numpy as np def fitfunc(x,p): if x>p: return xp else: return -(xp) fitfunc_vec = np.vectorize(fitfunc) #vectorize so you can use func with array def fitfunc_vec_self(x,p): y = np.zeros(x.shape) for i in range(len(y)): y[i]=fitfunc(x[i],p) return y x=np.arange(1,10) y=fitfunc_vec_self(x,6)+0.1*np.random.randn(len(x)) popt, pcov = so.curve_fit(fitfunc_vec_self, x, y) #fitting routine that gives error print popt print pcov 

    Salida:

     [ 6.03608994] [[ 0.00124934]] 

    ¿No podría simplemente reemplazar fitfunc con

     def fitfunc2(x, p): return np.abs(xp) 

    que luego produce algo como

     >>> x = np.arange(1,10) >>> y = fitfunc2(x,6) + 0.1*np.random.randn(len(x)) >>> >>> so.curve_fit(fitfunc2, x, y) (array([ 5.98273313]), array([[ 0.00101859]])) 

    Usando una función de conmutación y / o bloques de construcción como where reemplazar twigs, esto debería escalar hasta expresiones más complicadas sin necesidad de llamar vectorize .

    [PS: el errfunc en tu ejemplo de mínimos cuadrados no necesita ser un lambda. Usted podría escribir

     def errfunc(p, x, y): return array_fitfunc(p, x) - y 

    en cambio, si te ha gustado.]