質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

ただいまの
回答率

88.04%

Pythonで、複数の座標情報から、その座標を通る三角関数や指数・対数関数の式あるいは近似式を導出したい

受付中

回答 2

投稿 編集

  • 評価
  • クリップ 2
  • VIEW 7,164

score 15

数学寄りの質問で大変申し訳ないのですが、
Pythonで、複数の座標情報(座標情報は、いくつでも取得できることを想定。もちろん、利用する座標情報はなるべく少ない方が良い)を元に、その座標を通る三角関数や指数・対数関数の式あるいは近似式を導出したいと考えております(その式の係数を求めることになります)。なお、三角関数、指数・対数関数のどの関数の式であるかがわかっている前提とします。

今のところ、うまく導き出す方法を考えられておらず、
・座標を用いて計算を行い、導き出す
・Pythonで利用出来る関数に、導出してくれる関数があり、それを用いることで導き出す
・wolframのようなサイトがあって、そこのAPIに導出してくれる関数があり、それを用いることで導き出す
・例えば正弦波なら、y = asin(bx) + cの式の、a, b, cのとりうる範囲を限定し、二分探索などで、座標が一致しているかチェックすることで導き出す
という4つのぼんやりした方法しか思いつきません。

ネット上で関連した記事を探してみたのですが、Excelのソルバーや物理学(正弦波、余弦波、正接波など)などの記事しか見つかりませんでした。

なにか他に良い、式の導出(式の係数を求める)方法はありませんでしょうか?
プログラム関連のQ&Aサイトにこういった質問を投稿するのはお門違いな気もするのですが、どなたかご教授頂けないでしょうか?
  • 気になる質問をクリップする

    クリップした質問は、後からいつでもマイページで確認できます。

    またクリップした質問に回答があった際、通知やメールを受け取ることができます。

    クリップを取り消します

  • 良い質問の評価を上げる

    以下のような質問は評価を上げましょう

    • 質問内容が明確
    • 自分も答えを知りたい
    • 質問者以外のユーザにも役立つ

    評価が高い質問は、TOPページの「注目」タブのフィードに表示されやすくなります。

    質問の評価を上げたことを取り消します

  • 評価を下げられる数の上限に達しました

    評価を下げることができません

    • 1日5回まで評価を下げられます
    • 1日に1ユーザに対して2回まで評価を下げられます

    質問の評価を下げる

    teratailでは下記のような質問を「具体的に困っていることがない質問」、「サイトポリシーに違反する質問」と定義し、推奨していません。

    • プログラミングに関係のない質問
    • やってほしいことだけを記載した丸投げの質問
    • 問題・課題が含まれていない質問
    • 意図的に内容が抹消された質問
    • 過去に投稿した質問と同じ内容の質問
    • 広告と受け取られるような投稿

    評価が下がると、TOPページの「アクティブ」「注目」タブのフィードに表示されにくくなります。

    質問の評価を下げたことを取り消します

    この機能は開放されていません

    評価を下げる条件を満たしてません

    評価を下げる理由を選択してください

    詳細な説明はこちら

    上記に当てはまらず、質問内容が明確になっていない質問には「情報の追加・修正依頼」機能からコメントをしてください。

    質問の評価を下げる機能の利用条件

    この機能を利用するためには、以下の事項を行う必要があります。

質問への追記・修正、ベストアンサー選択の依頼

  • swordone

    2015/07/26 13:46

    与えられた座標を通る(または近似する)三角関数なり指数・対数関数なりの係数を求めるという意味に見えるのですが,違うのですか?

    キャンセル

  • tukejonnyBot

    2015/07/26 14:01

    はい。そうなります。
    説明不足でした。申し訳ありません・・・

    キャンセル

回答 2

+1

最小二乗法が使いやすいと思います.
いくつか追加のモジュールが必要です.
Macで試したのですが,NumPyとMatplotLibとSciPyをインストールすれば良いと思います.
$ pip3 install numpy
$ pip3 install matplotlib
$ pip3 install scipy

その上で,以下のコードを実行してみてください.
テストに使ったデータは y = x ** 2 + 1 の座標です.
2次関数の一般形 y = a * x ** 2 + b * x + c を定義して,
そのパラメータを推定しています.
a_fit,b_fit,c_fitにそれぞれ推定したパラメータが格納されます.
三角関数,指数・対数関数などでも利用できると思うので,試してみてください.

import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize

x = np.array([0, 1, 2, 3, 4, 5, 6, 7])
y = np.array([1, 2, 5, 10, 17, 26, 37, 50])

def fit_func(parameter,x,y):
    a = parameter[0]
    b = parameter[1]
    c = parameter[2]
    residual = y-(a*x**2+b*x+c)
    return residual

parameter0 = [0.,0.,0.]
result = optimize.leastsq(fit_func,parameter0,args=(x,y))
a_fit=result[0][0]
b_fit=result[0][1]
c_fit=result[0][2]
print(a_fit,b_fit,c_fit)

#PLot
plt.figure(figsize=(8,5))
plt.plot(x,y,'bo', label='Exp.')
plt.plot(x,a_fit*x**2+b_fit*x+c_fit,'k-', label='fitted parabora', linewidth=10, alpha=0.3)
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.legend(loc='best',fancybox=True, shadow=True)
plt.grid(True)
plt.show()

参考
フィッティング3(最小二乗法,optimizeの利用) — Mechanical Design Lab. of TUMSAT 1.0 documentation

投稿

編集

  • 回答の評価を上げる

    以下のような回答は評価を上げましょう

    • 正しい回答
    • わかりやすい回答
    • ためになる回答

    評価が高い回答ほどページの上位に表示されます。

  • 回答の評価を下げる

    下記のような回答は推奨されていません。

    • 間違っている回答
    • 質問の回答になっていない投稿
    • スパムや攻撃的な表現を用いた投稿

    評価を下げる際はその理由を明確に伝え、適切な回答に修正してもらいましょう。

  • 2015/07/26 23:18

    回答ありがとうございます!!
    最小二乗法という方法があるのですね・・・

    勉強不足で大変申し訳ないのですが、一つわからないところがあります。
    optimize.leastsq()の返り値が二次元配列になっているのですが、値はどのように格納されているのでしょうか?
    単純にパラメータ1つ目からresult[0][1], result[0][2], ..., result[1][0], result[1][1],...となっているのでしょうか?

    実際に三角関数(sin)のグラフを描画してみました。描画されているグラフはあっていそうなのですが、パラメータが少しおかしいように思えます

    ```Python
    import numpy as np
    import matplotlib.pyplot as plt
    from scipy import optimize
    from math import*

    """ test data """
    x = np.array([r for r in range(1, 1000)])
    y = np.array([sin(pi * (r) / 180.0) for r in range(1, 1000)])
    """ test data """

    def calcSin(angle):
    rad = pi * angle / 180.0
    for r in range(0, len(angle)):
    angle[r] = sin(rad[r])
    return(angle)

    def fit_func(parameter, x, y):
    a = parameter[0]
    b = parameter[1]
    c = parameter[2]

    residual = y-(a*calcSin(b * x)+c)
    return(residual)

    parameter0 = [0.,0.,0.]
    result = optimize.leastsq(fit_func, parameter0, args=(x, y))
    print "result is " + str(result)
    a_fit = result[0][0]
    b_fit = result[0][1]
    c_fit = result[0][2]

    print "y = (" + str(a_fit) + ")sin((" + str(b_fit) + ")x) + " + str(c_fit)

    #PLot
    plt.figure(figsize=(8,5))
    plt.plot(x,y,'bo', label='Exp.')
    plt.plot(x,a_fit*calcSin(b_fit * x)+c_fit,'k-', label='fitted parabora', linewidth=10, alpha=0.3)
    plt.xlabel('X-axis')
    plt.ylabel('Y-axis')
    plt.legend(loc='best',fancybox=True, shadow=True)
    plt.grid(True)
    plt.show()
    ```

    出力結果(グラフは除く)
    y = (0.0)sin((0.0)x) + 0.0478855562909

    これだとxの値に関わらず、yが0.0478855562909となるように見えてしまいます。
    出力結果として、これは正しいのでしょうか?

    キャンセル

  • 2015/07/27 01:24

    scipy.optimize.leastsqの返り値がresultだったとき,
    これは2次元配列になっていて,
    result[0]が推定パラメータ,result[1]がパラメータの誤差だそうです.
    result[0]ばかり気にしていたので,result[1]が真に何を意味するのかはわかりません.

    キャンセル

0

コメントだとコードが見にくくなるので,別の所へ.
確かに,色々やってみましたが上手く行きませんでした.
入力データを変えた所上手く動いたのですが,なぜだかわからずじまいです.
もしかしたらどこかで,浮動小数の情報落ち,桁落ちなどで
値の更新がされなくなってしまったのかもしれません.
(最小二乗法は,計算を何度も繰り返して最小になる値を求めるので)

pylabのsin関数は配列に対しても行えると聞いたので,使っています.
よくよく考えたら,yの式をmyfitか何かで宣言しておいて,
mydiffか何かで入力データとの差を取れば良いかなと思いました.

また,初期値も[0.,0.,0.]でやるとエラーになることを確認しています.
最小二乗法では最小値を推定する計算があるのですが,
その方法を使う時にいくつか制限があるということです.
この場合は1番目の値が0だとだめなようです.
(ヤコビ行列が陽,と私にもさっぱりの定義でした)

import matplotlib.pyplot as plt
from scipy import optimize
from math import *
import pylab

""" test data """
x = [(pi * (r) / 180.0) for r in range(1, 1000)]
y = pylab.sin(x)
""" test data """

def myfit(p, x):
    return(p[0]*pylab.sin(p[1]*x)+p[2])

def mydiff(p, x, y):
    return(y-myfit(p,x))

p0 = [1.0,1.0,1.0]
result = optimize.leastsq(mydiff, p0, args=(x,y))
pr = result[0]
print pr

#Plot
plt.figure(figsize=(8,5))
plt.plot(x,y,'bo', label='Exp.')
plt.plot(x,myfit(pr,x),'k-', label='fitted curb', linewidth=10, alpha=0.3)
plt.xlabel('X-axis')
plt.ylabel('Y-axis')
plt.legend(loc='best',fancybox=True, shadow=True)
plt.grid(True)
plt.show()

投稿

編集

  • 回答の評価を上げる

    以下のような回答は評価を上げましょう

    • 正しい回答
    • わかりやすい回答
    • ためになる回答

    評価が高い回答ほどページの上位に表示されます。

  • 回答の評価を下げる

    下記のような回答は推奨されていません。

    • 間違っている回答
    • 質問の回答になっていない投稿
    • スパムや攻撃的な表現を用いた投稿

    評価を下げる際はその理由を明確に伝え、適切な回答に修正してもらいましょう。

15分調べてもわからないことは、teratailで質問しよう!

  • ただいまの回答率 88.04%
  • 質問をまとめることで、思考を整理して素早く解決
  • テンプレート機能で、簡単に質問をまとめられる

同じタグがついた質問を見る

  • トップ
  • Pythonに関する質問
  • Pythonで、複数の座標情報から、その座標を通る三角関数や指数・対数関数の式あるいは近似式を導出したい