Python Scipy 中的麦克斯韦分布答案

【问题标题】：Maxwellian Distribution in Python ScipyPython Scipy 中的麦克斯韦分布
【发布时间】：2020-11-27 17:37:36
【问题描述】：

在我感兴趣的文章中，它指出数据很好地表示为麦克斯韦分布，它还提供了平均速度 (307 km/s) 和 1 sigma 不确定性 (47 km/s) 的分布.

使用提供的值，我尝试重新生成数据，然后使用 python scipy.stats 将其与麦克斯韦分布拟合。

正如here 中所述，scipy 中的 maxwell 函数需要两个输入，1) "loc" 移动 x 变量和 2) "a" 参数，对应于 maxwell-Boltzmann 方程中的参数 "a" .

在我的例子中，我没有这两个参数，因此使用wiki page 中的均值和方差 (sigma^2) 描述，我尝试计算“a”和“loc”参数。 mean 和 sigma 参数都只依赖于“a”参数。

我遇到的第一个问题是我从 Mean (a = 192.4) 和 sigma (a = 69.8) 得到的“a”参数彼此不同。第二个问题是我不知道如何从 Mean 和 sigma 中获得准确的 loc (shift) 值。

根据分布的形状（图中平均速度值落在图中，请查看图 2），我尝试猜测“loc”值以及从 sigma 获得的“a”值（a = 69.8），我已经生成并拟合了数据。大约看起来是正确的，但我不知道我上面提到的问题的答案，我需要一些专家的指导。感谢您的帮助。

import matplotlib.pyplot as plt
import math
from scipy.stats import norm
import random
import numpy as np
import scipy.optimize
from scipy.stats import maxwell

samplesize = 100000

mean = 307
sigma = 47
loc = 175 #my guess
a_value = np.sqrt((sigma**2 * math.pi)/(3*math.pi - 8)) #calculated based on wiki description

fig, axs = plt.subplots(1)
v_2d = maxwell.rvs(loc, a_value, size=samplesize) #array corresponding to 2D proper motion obtained from Hubbs
mean, var, skew, kurt = maxwell.stats(moments='mvsk')

N, bins, patches = plt.hist(v_2d, bins=100, density=True, alpha=0.5, histtype='bar', ec='black')
maxx = np.linspace(min(v_2d), max(v_2d), samplesize)

axs.plot(maxx, maxwell.pdf(maxx, loc, a_value), color=colorset[6], lw=2, label= r'$\mathdefault{\mu}$ = '+'{:0.1f}'.format(mean)+r' , '+r'$\mathdefault{\sigma}$ = '+'{:0.1f}'.format(sigma))

axs.set(xlabel=r'2-D Maxwellian speed (km s$^{-1}$)')
axs.set(ylabel='Frequency')
plt.legend(loc='upper right')

【问题讨论】：

我相信，由于您的“loc”不为零，您需要重新评估期望值和 sigma 的公式。在 MB 分布中，平均值为 2*asqrt(2/pi)。不同的位置不一样。在这种情况下，您需要将 xf(x) 集成到您的域中，其中 f(x) 是带有“loc”的 MB 分布。
你能提供文章的链接吗？您说“它还为分布提供了平均速度（307 km/s）和 1 sigma 不确定性（47 km/s）”； “1 sigma 不确定性”到底指的是什么？
@WarrenWeckesser 他们表示：该值表示最低有效数字的误差（在 68% 的置信水平上）。所以这是一个 1 sigma 的拟合误差。如果你愿意，我仍然可以提供论文。
我仍然不明白值 47 是什么“不确定性”。如果那是样本均值 307 的不确定性，那么这与分布的标准差不一样。
您能说出实际的问题吗？似乎 sigma 是拟合均值的置信区间，它是拟合参数的统计参数，而不是分布参数。你文章的参考文献是什么？

标签： python numpy scipy statistics distribution

【解决方案1】：

好吧，平均值受位置影响，而 sigma 不会。所以从 sigma 计算a，计算平均值，如果 loc=0，找到差异并将其分配给位置，采样 100K RV 以检查是否采样均值/标准差足够接近。

代码、Python 3.8、Windows 10 x64

import numpy as np

from scipy.stats import maxwell

σ = 47
μ = 307

a = σ * np.sqrt(np.pi/(3.0*np.pi - 8.0))
print(a)

m = 2.0*a*np.sqrt(2.0/np.pi)
print(m) # as if loc=0

loc = μ - m
print(loc)

print("----------Now test--------------------")

# sampling
q = maxwell.rvs(loc=loc, scale=a, size=100000)

print(np.mean(q))
print(np.std(q))

作为我得到的输出

306.9022249667151
47.05319429681308

够好吗？

【讨论】：