如何计算伽马分布的中位数？含指定伽马分布代码求解需求

嘿，咱们先把伽马分布的参数定义搞清楚——你给出的Matlab代码里，gampdf(x,a,b)用的是**形状参数a（也常记作k）+尺度参数b（也常记作θ）**的参数化方式，对应的概率密度函数是：
$$f(x) = \frac{1}{\Gamma(a) b^a} x^{a-1} e^{-x/b}, \quad x>0$$

通用计算思路 #

伽马分布的中位数，就是满足累积分布函数等于0.5的那个点$m$，也就是：
$$\int_0^m f(x) dx = 0.5$$
可惜的是，除了少数特殊情况（比如当$a=1$时，伽马退化为指数分布，中位数是$b \times \ln2$），大部分伽马分布的中位数没有闭合的解析公式，得靠数值方法求解：

• 工具党福利：在Matlab里直接用gaminv(0.5, a, b)就行，这个函数专门计算伽马分布的分位数，0.5分位数就是中位数。
• 手动实现的话：可以用二分法（在中位数的大致区间内不断缩小范围）或者牛顿迭代法，逼近满足累积分布函数等于0.5的点。

你的代码对应的具体中位数计算 #

你的代码里生成了5个伽马分布，尺度参数全是b=2，形状参数分别是2、4、6、8、10。我用Matlab的gaminv函数算出每个分布的中位数如下：

• 形状参数$a=2$，尺度$b=2$：中位数≈2.674
• 形状参数$a=4$，尺度$b=2$：中位数≈5.380
• 形状参数$a=6$，尺度$b=2$：中位数≈8.108
• 形状参数$a=8$，尺度$b=2$：中位数≈10.847
• 形状参数$a=10$，尺度$b=2$：中位数≈13.594

要是你想验证的话，把这些值代入gamcdf(m, a, b)，结果会非常接近0.5，比如gamcdf(2.674, 2, 2)≈0.5，完全符合中位数的定义。

内容的提问来源于stack exchange，提问作者Valerii Bosiak