5

区分的に定義された関数をPythonのデータセットに適合させようとしています。かなり前から探していましたが、可能かどうかはわかりません。

私がやろうとしていることの印象をつかむために、次の例を見てください（これは私にとってはうまくいきません）。ここでは、シフトされた絶対値関数（f（x）= | xp |）を、適合パラメーターとしてpを使用してデータセットに適合させようとしています。

import scipy.optimize as so
import numpy as np

def fitfunc(x,p):
   if x>p:
      return x-p
   else:
      return -(x-p)

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

Pythonでこれを実現する方法はありますか？

Rでこれを行う方法は次のとおりです。

# Fit of a absolute value function f(x)=|x-p|

f.lr <- function(x,p) {
    ifelse(x>p, x-p,-(x-p))
}
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')

さらに調査した後、私はそれを行う方法を見つけました。このソリューションでは、エラー関数を自分で定義する必要があるという事実は好きではありません。さらに、なぜそれがこのラムダスタイルでなければならないのかよくわかりません。したがって、あらゆる種類の提案やより洗練されたソリューションを歓迎します。

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()
4

2 に答える 2

5

ここでこれを完了するために、私は問題に対する私自身の最終的な解決策を共有します。私の元の質問に近づくには、ベクトル化された関数を自分で定義するだけで、を使用しないでくださいnp.vectorize。

import scipy.optimize as so
import numpy as np

def fitfunc(x,p):
   if x>p:
      return x-p
   else:
      return -(x-p)

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

出力： 

[ 6.03608994]
[[ 0.00124934]]
于 2013-08-09T14:41:59.727 に答える
2

単にfitfuncを次のように置き換えることはできませんか

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

次に、次のようなものを生成します

>>> 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]]))

スイッチ関数やwhereブランチの置き換えなどのビルディングブロックを使用すると、を呼び出すことなく、より複雑な式にスケールアップする必要がありますvectorize。

[追記：errfunc最小二乗法の例では、ラムダである必要はありません。あなたは書くことができます

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

代わりに、気に入ったら。]

于 2012-06-21T02:50:30.003 に答える