Подгонка кривой с использованием реального распределения Гаусса в Python ⇐ Python

[Anonymous]

У меня есть набор точек, их разбросанное изображение похоже на нормальное распределение Гаусса. При поиске в Интернете можно найти множество примеров Python, как подогнать кривую к точкам. Они основаны на:

Код: Выделить всё

def Gauss1(X, C, mu, sigma):
return C * np.exp(-(X-mu) ** 2 / (2 * sigma ** 2))

и

Код: Выделить всё

from scipy.optimize import curve_fit

Это дает очень хорошую кривую. Проблема в том, что Gauss1 не является нормальным распределением Гаусса, а должно быть:

Код: Выделить всё

def Gauss2(X, mu, sigma):
return (1 / (sigma * np.sqrt(2 * np.pi))) * np.exp(-(X-mu) ** 2 / (2 * sigma ** 2))

При использовании Gauss2 выводится следующее сообщение без всплывающего окна:

Код: Выделить всё

OptimizeWarning: Covariance of the parameters could not be estimated

Можно ли построить аппроксимирующую кривую, используя реальное распределение Гаусса?

Код: Выделить всё

popt2, pcov2 = curve_fit(Gauss2, x, y)

OptimizeWarning: Covariance of the parameters could not be estimated

Думаю, баллы не нужны, но могу опубликовать 25 строк.

Код: Выделить всё

import matplotlib.pyplot as plt # Follow the convention
import numpy as np              # Library for working with arrays
from scipy.optimize import curve_fit

######################################## G L O B A L S
# X the ratings after grouping by 23 the chess players from FIDE 2024 Oct Standard
# Y is the average of players in the averaged ratings
XY = (
(1411, 231), (1434, 271), (1457, 281), (1480, 287),
(1503, 292), (1526, 298), (1549, 293), (1572, 299),
(1595, 300), (1618, 303), (1641, 304), (1664, 308),
(1687, 301), (1710, 311), (1733, 304), (1756, 291),
(1779, 286), (1802, 283), (1825, 279), (1848, 268),
(1871, 260), (1894, 243), (1917, 229), (1940, 221),
(1963, 203), (1986, 169), (2009, 134), (2032, 112),
(2055, 102), (2078,  94), (2101,  81), (2124,  76),
(2147,  70), (2170,  61), (2193,  54), (2216,  45),
(2239,  42), (2262,  35), (2285,  31), (2308,  27),
(2331,  26), (2354,  22), (2377,  20), (2400,  18),
(2423,  15), (2446,  10), (2469,  10), (2492,   7),
(2515,   6), (2538,   5), (2561,   4), (2584,   3),
(2607,   3), (2630,   2), (2653,   2), (2676,   2),
(2699,   1), (2722,   2), (2745,   1), (2768,   1),
(2791,   1), (2837,   1))

class c:  # Constants
X_LABEL_STEP = 50
Y_LABEL_STEP = 10
A4 = (11.69, 8.27)  # W x H in inches

###################################### F U N C T I O N S

def Graph_Gauss():
global c, XY

X = [] ; Y = []
for e in XY:
X += [e[0]]
Y += [e[1]]

# Set up the limits for X, ratings
minX = min(X) ; maxX = max(X)
# Set up the limits for Y, number of players with that rating
minY = min(Y) ; maxY = max(Y)

fig, ax = plt.subplots()
fig.set_size_inches(c.A4)
fig.suptitle("Rating Distribution", fontsize=18, y=0.935)

x1  = []  # Make the x axle
xlb = []  # Make the labels for x
stp = c.X_LABEL_STEP
for k in range(stp*(minX//stp), stp*(maxX//stp + 2), stp):
x1 += [k]
xlb.append(f"{k}")
ax.set_xticks(ticks=x1, labels=xlb, rotation=270)
ax.set_xlabel("Rating", fontsize=15, labelpad=7)

stp = c.Y_LABEL_STEP
yticks = np.arange(stp*(minY//stp), stp*(maxY//stp + 2), stp)
ax.set_ylabel("Number of Players", fontsize=15, rotation=270, labelpad=18)
ax.set_yticks(ticks=yticks)
ax.grid(which="major", axis="both")

ax.scatter(X, Y, color="#000FFF", marker='o', s=14)

# plt.show()

# Fit a normal distribution, aka Gaussian fitting curve
# https://www.wasyresearch.com/tutorial-python-for-fitting-gaussian-distribution-on-data/
# mu = mean = sum(x) / len(x) ;  In our case sum(x * y) / sum(y)
# sigma = standard deviation = sqrt((sum((x - mean)**2) / len(x))
#
# https://www.scribbr.com/statistics/normal-distribution/
# 1/(sigma*sqrt(2*pi)) * e**(-(x - mean)**2 / (2 * sigma**2))
# Calculating the Gaussian PDF values given Gaussian parameters and random variable X
def Gauss1(X, C, mu, sigma):
return C * np.exp(-(X-mu)**2 / (2 * sigma**2))
def Gauss2(X, mu, sigma):
return (1/(sigma*np.sqrt(2*np.pi))) * np.exp(-(X-mu)**2 / (2 * sigma**2))

x = np.array(X)
y = np.array(Y)
mu    = sum(x * y) / sum(y)
sigma = np.sqrt(sum(y*(x - mu)**2)/sum(y))

print(f"{1/(sigma*np.sqrt(2*np.pi))=:.2f} {mu=:.2f} {sigma=:.2f}")

popt1, pcov1 = curve_fit(Gauss1, x, y, p0=[max(y), mu, sigma], maxfev=5000)
popt2, pcov2 = curve_fit(Gauss2, x, y)  # This generates:
# OptimizeWarning: Covariance of the parameters could not be estimated

yg1 = Gauss1(x, *popt1)
yg2 = Gauss2(x, *popt2)

ax.plot(x, yg1, color="#FF0F00", linewidth=3,
label=f"Normal: mu={popt1[1]:.2f}, sigma={popt1[2]:.2f}")
plt.legend(fontsize=12)

breakpoint()  # ???? DEBUG

plt.show()
pass  # To set a breakpoint

#######################################################################
if __name__ == '__main__':
# breakpoint()  # ???? DEBUG, to set other breakpoints
Graph_Gauss()

Добавлено 28 октября 2024 г.
Привет, Альдегид и Lastchance!
Большое спасибо за ваши ответы, все становится яснее. Теперь у меня два вопроса:
• Что такое функция распределения Гаусса, где находится сигма; Я не смог понять это из множества результатов поиска в Интернете.
• Для использования Gauss2() вы нормализовали значения y на 230 000, и это привело к получению хорошо подходящей кривой. Я попробовал несколько других чисел, 240 000 даже лучше подходят. Откуда взялось это или подобное число. Я попробовал np.sum(y)=8,241 и область под точками np.sum((x[1:]-x[:-1]) * ((y[1:]+y[:-1]) /2)) = 186 898, они плохо подошли.
Так какой метод правильный, как из заданных 62 баллов вычесть 240 000 или что-то подобное. Если вы можете порекомендовать текст или другую книгу, где я могу найти ответ, это тоже хорошо.
Еще раз спасибо,
Вилмос


Подробнее здесь: https://stackoverflow.com/questions/791 ... -in-python