pythonで任意の範囲の変曲点を求める

Question

前提・実現したいこと

pythonを使用して測定したデータの変曲点を求めるシステムを作成したいと思っています
測定データは加算平均処理を行っているので、測定時間のいずれの中から雑音を気にせず単なる変曲点を求めたいです。

発生している問題・エラーメッセージ

resultの配列が負になるポイントが変曲点とした場合、変曲点にマーカーをつける方法を教えてください。
!](b2b89115a7abe40e5756be2cf2c9401e.png)
11番目の区間が負になっているので12番目の点（グラフ右）が変曲点ということでしょうか？
resultの配列→[0.20171338950000001, 0.21015402949999998, 0.19502000350000001, 0.20433273499999999, 0.16845394799999999, 0.168081231, 0.14231983350000002, 0.12233689899999997, 0.08080066050000001, 0.041403474500000009, -0.0056071255000000042, -0.039149066000000003, -0.07469213949999999, -0.11799488900000002, -0.14328871300000001, -0.18012878500000001, -0.17326162, -0.20117560449999999, -0.20206969749999998, -0.21288156799999999, -0.2059022455, -0.2068159155, -0.207493024, -0.20301948749999998, -0.18277358499999999]

該当のソースコード

import numpy as np
import matplotlib.pyplot as plt

data = np.loadtxt("60－70kasannomi.txt", unpack="True")
h = 1
result=[]

for k in range(125,150):
    point=data[k]-data[k-1]/2*h
    result.append(point)
    k = k+1

plt.subplot(1,2,1)
plt.plot(data,"o-")
plt.ylim(-1,3)

plt.subplot(1,2,2)
plt.plot(data[125:150],"o-")
plt.ylim(-0.5,0.5)

plt.show()
print(result)

試したこと

補足情報（FW/ツールのバージョンなど）

python3.6
開発環境はpycharmを使用しています

Comments

[tiitoi]

1次元配列から傾きの正負が変わる場所を求めたいということでしょうか？

[origami_]

そうです。

[tiitoi]

でしたら、numpy を使うと簡単に計算できます。回答に記載しました。

Answer 1

微分の定義は

f(x)について

(f(t+Δt)-f(t))/Δt　のlim Δt->0のことです。

f(x)が一時微分可能な数式で表現できるようなものなら公式から微分が出てくるわけですが・・・
f(x)が未知であって単にあるx(i)(i=0,1,2,...)の値での測定値しかわからないようなものは数式による微分は求めようがないです。そういう場合は

(f(t+Δt)-f(t))/Δt

を

(f(x(i+1))-f(x(i)))/(x(i+1)-x(i))

とみなして「観測可能なx(i)の地点での大まかな微分相当の値」として考える他に方法はありません。
まずこの点をどう認識しておられるのかご質問からはよくわかりませんでした。

ちなみに測定値は普通雑音を伴うものですので、「雑音を気にせず単なる変曲点を求めればよい」のか「一定レベルの雑音の影響を排除したいのか」で方法は異なります。通常は後者を考えるでしょう。

また「測定しているものが本来は滑らかな曲線であるとみなし、本来の変曲点を推定する」のが目的なのか「測定時間のいずれかの中から変曲点を求めるのか」などでも方法論は異なってきます。

あなたが「何を目的にどういう精度の変曲点を求めようとしているのか」がわからなければ適切な回答は難しいと思います。

Comments

[origami_]

回答ありがとうございます。
微分については数値微分という方法で考えるのですね。不勉強でした。
測定データは加算平均処理を行っているので、測定時間のいずれの中から雑音を気にせず単なる変曲点を求めたいです。

[KSwordOfHaste]

そうでしたか。それならばf(x-1)とf(x),f(x)とf(x+1)の変量を求めその符号が逆転する点を変曲点と考えるのを基本とすればよいと思います。f(x)が0,1,1,0なんて変化をするときどこを変曲点とみなすかなんてことも考えるべきかも知れません。

[origami_]

返信ありがとうございます。勉強になりました。
いろいろ試してまた質問させていただきます。

Answer 2

以下の手順で求めればよいと思います。

1. numpy.gradient() で傾きを求める。
2. numpy.sign() で傾きの正負をとる。
3. numpy.diff() をとると、傾きの正負が変わる場所は2、そうでない場所は0の配列が得られる。
4. indices = numpy.where() で傾きの正負が変わる場所のインデックスを取得する。
5. xs[indices] で傾きの正負が変わる場所を求める。

import matplotlib.pyplot as plt
import numpy as np

def f(x):
    return x**3 - 3 * x + 1

xs = np.linspace(-3, 3, 30)
ys = f(xs)

f_prime = np.gradient(ys)  # 差分近似
indices = np.where(np.diff(np.sign(f_prime)))[0]  # 変曲点を求める。
infections = xs[indices]

print('変曲点', infections)  # 変曲点 [-1.13793103  0.93103448]

fig, ax = plt.subplots(figsize=(7, 7))
ax.plot(xs, ys, 'bo-', ms=2)
ax.plot(infections, f(infections), 'ro', ms=5)
plt.show()

赤の点が変曲点

Comments

[origami_]

回答ありがとうございます。
16行目のaは何の変数でしょうか？
測定点の中から変曲点を求めたいのですが、ファイルを読みこんだ場合xs,ysの変数はどうすればよいでしょうか。

[tiitoi]

a は転記ミスです。回答を修正しました。
xs, ys のデータ点をファイルから読み込む場合は np.loadtxt() を使ってください。
使い方は公式リファレンスか https://deepage.net/features/numpy-loadsavetxt.html のサイトが参考になるかと思います。

[origami_]

解決しました。各種関数の説明は大変助かりました。ありがとうございます。