pythonでの強度分布のプロット方法

Question

前提・実現したいこと

例えば，次のような値を持つデータがあるときに，

5
5
5
5
9
9
16

任意の幅で強度をガウス分布をプロットする方法がわからないです．

具体的には化学でのMSスペクトルやXPSスペクトルのシュミレーションで得られた値をプロットしたいイメージです．

最終的なイメージは以下の画像のようなプロットをしたいです．
青が上で与えられるようなデータで，赤がガウス分布のプロットのイメージです．

Comments

[ppaul]

その近い値ごとに分布を持ったプロット
の、「近い」はどう判断したいのですか？

[hisanori61]

サンプルデータとしては整数値を与えていますが，実際のデータは0.00??ぐらいで近い値が存在し，離散的な状態を取ります．
この場合ですと，txtのデータの整数値が同じものを近いと考えていただけたらと思います．

[ppaul]

四捨五入か切り上げがか切り捨てすれば良いということですか。
それだと、サンプルでないデータの時には想定外のことがおきるかもしれませんが、それでも良いのでしょうか？

[hisanori61]

データとして別物になるので四捨五入する必要性がよくわかりません．桁としてどこで丸めるのかがわからないのですが，想定外のことが起こる場合は丸める桁を変えられれば問題ないかと思います．

[ppaul]

1.1, 1.6, 2.1, 2.6, 3.1, 3.6というような数字があたえられたとします。
「整数値が同じものを近い」だと、[1.1, 1.6], [2.1, 2.6], [3.1, 3.6]と3個に分類されますね。それでよいのですか、という質問です。

[hisanori61]

ややこしくてすみません．近いの判断が機械的というよりは人間が見たときの視覚的に判断できるようにしたいのでその近さだと一概にいえないです．

実データの一部だと，以下のものがあります．この場合だと，整数値部分の[739], [719], [707, 706]の3つに分類したいです．
[739.2065427, 739.2055331, 739.2053889, 739.2041970, 739.2032718, 739.2025915, 739.2019031, 739.2010051, 739.2003792, 739.1997615, 739.1985642, 739.1977914, 739.1666096, 739.1658532, 739.1657307, 739.1649525, 739.1618585, 739.1610313, 739.1596652, 739.1588380, 739.1567835, 739.1560951, 739.1536841, 739.1529821, 719.6762059, 719.2190156, 707.9640013, 707.5060300, 706.9127546, 706.3545658]

[ppaul]

このように分類することをクラスタリングといいます。
クラスタリングは難しい問題です。
https://aizine.ai/wp-content/uploads/2019/05/image5-%E2%80%93-30.jpg
のクラスタリングについて、
https://aizine.ai/clustering0522/
に書かれているのでお読みください。
直前のコメントに書かれている例で提示された「実データの一部」についてはうまく動くコードを書くのは簡単です。しかし、実データの一部を分類できたとしても、実データの全体ではどうすれば期待されるクラスタリングになるのかは、通常のプログラムでは分からないのです。
そのために、機械学習の分野ではさまざまなクラスタリング技法が考案されていますが、どれを選ぶかは、それを作る人がどういうクラスタリングを期待しているかによります。
例を出されるけれども、「近い」はどう判断したいのですか？という質問への回答が返ってこないということは、「期待するクラスタリング」が言語化できない状態だということのように思います。

[hisanori61]

趣旨がずれてきてる気がするので改めていいますが，これは図示をしたいのが主目的です．

いかに簡単にプロットするかのために，判断を機械的にしたいのはわからなくもないです．近いの判断は重要かもしれないのですが，そこは場合によっては人間が分離するでもいいと思っています．（予めデータを分けるなどで対応するのも仕方ないかと思っています）

近いに関してはデータによって変わるので明言できないだけであって，言語化できてないというのはよくわからないです．先に例示した例のような基準で近いを判断したいのですが伝わらないでしょうか？その基準を明確にするというのは，データの間隔や分布を作るときの幅でその近いに含まれる範囲が変わると思うので難しいと思います．

[ppaul]

では、入力データはクラスタリング済みである、つまり、
[[706.3545658, 706.9127546, 707.50603, 707.9640013], [719.2190156, 719.6762059], [739.1529821, 739.1536841, 739.1560951, 739.1567835, 739.158838, 739.1596652, 739.1610313, 739.1618585, 739.1649525, 739.1657307, 739.1658532, 739.1666096, 739.1977914, 739.1985642, 739.1997615, 739.2003792, 739.2010051, 739.2019031, 739.2025915, 739.2032718, 739.204197, 739.2053889, 739.2055331, 739.2065427]]
が与えられていると考えて表示だけすれば良いのですか？

[hisanori61]

そうです！この分けてあるデータを最初の図のような棒＋分布のようなプロットで図示できないかと考えています．

Answer 1

入力データは分類済みと言うことなら以下のようにすれば良いでしょう。

import numpy as np
from scipy.stats import norm
import matplotlib.pyplot as plt

data = [[706.3545658, 706.9127546, 707.50603, 707.9640013], [719.2190156, 719.6762059], [739.1529821, 739.1536841, 739.1560951, 739.1567835, 739.158838, 739.1596652, 739.1610313, 739.1618585, 739.1649525, 739.1657307, 739.1658532, 739.1666096, 739.1977914, 739.1985642, 739.1997615, 739.2003792, 739.2010051, 739.2019031, 739.2025915, 739.2032718, 739.204197, 739.2053889, 739.2055331, 739.2065427]]
colors = ['r', 'b', 'g']

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)

for i, datum in enumerate(data):
    print(i, datum)
    a = np.array(datum)
    average = a.mean()
    std = np.std(a)
    X = np.linspace(average-5, average+5, 1000)
    Y = norm.pdf(X, loc=average, scale=1.6/a.shape[0])
    ax.plot(X[Y>0.0001], Y[Y>0.0001], color=colors[i])

plt.show()

実行結果

Comments

[hisanori61]

ありがとうございます！
イメージしていた感じです．ここに棒（垂線など？）として各データの存在を示すことは可能でしょうか？

Answer 2

参考データ（txt）の意味がよくわからなかったのですが、このようなイメージでしょうか？

import numpy as np
import matplotlib.pyplot as plt

from scipy import stats


mu = [-0.0, 1.0, 2.0, 3.0]  # 平均
std = [0.1, 0.5, 0.2, 0.5]  # 標準偏差
intensity = [0.5, 1, 1.5, 2.0]  # 強度

x = np.linspace(-5,5,100)

for m, s, i in zip(mu, std, intensity):
  norm = stats.norm.pdf(x=x, loc=m, scale=s)
  plt.plot(x, norm*i)


plt.show()

Comments

[hisanori61]

回答していただいてありがとうございます．

自分の説明が分かりづらかったのですが，txtの値がいわゆるX軸の値として存在していて，Y軸方向の値を持っていない状態です（txtの値がすべて同じ強度あるいは数）．

いわゆる頻度分布のようなヒストグラムに近いのかもしれないのです．任意の区間で仕切ってプロットするわけではなく，各値（質問の図での棒）が存在するようにしておいて，その近い値ごとに分布を持ったプロットを作成するイメージです．

[taront]

あー、そういう意味でしたか。
データをグループ分けした後にフィッティングすれば正規分布のパラメータが得られると思いますが、サンプルデータのようにグループ毎の分布が無いデータでは適用できませんね（標準偏差が決まらないと思います）。

from scipy.stats import norm
param = norm.fit(sample)

ちょっと良い方法が思いつきません。失礼いたしました。

[taront]

補足ですが、正規分布を描くには平均値と標準偏差が必要ですので、「その近い値ごとに分布を持ったプロット」に関して、自分なりの定義を作らないとコード化は出来ないと思いますよ。

[hisanori61]

tarontさん，ありがとうございます．
近しい値を選んでそれごとにプロットを作らないといけないかと思いました．また，分布に関しては，ガウス分布を取るのを考えていたのですが，任意に分布の幅をつけようと考えていました．
自分の理解がついていってなかったかもしれないです．