截断分布曲线估计：右尾样本的分布识别与参数计算

当然可以，但必须采用截断分布的参数估计方法——你手上的样本是原分布在≥35阈值下的条件样本，直接用普通样本统计量（比如样本均值、标准差）得到的是截断后分布的参数，而非原分布的。下面是具体实现思路和代码示例：

1. 先探索样本的分布形态 #

先通过可视化快速判断样本的大致分布特征，比如直方图、核密度曲线：

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

# 加载样本
sample = np.array([53,53,51,49,49,49,48,48,48,47,47,47,47,47,47,46,46,46,46,46,45,44,44,44,43,43,43,43,43,43,43,43,43,43,43,43,42,42,42,42,42,42,42,42,42,42,41,41,41,41,41,41,41,41,41,41,41,41,40,40,40,40,40,40,40,40,40,40,40,40,39,39,39,39,39,39,39,39,39,39,39,39,38,38,38,38,38,38,38,38,38,38,38,38,38,37,37,37,37,37,37,37,37,37,37,37,37,37,37,37,37,36,36,36,36,36,36,36,36,36,36,36,36,36,36,36,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35,35])

# 画直方图+核密度估计
plt.hist(sample, bins=10, density=True, alpha=0.6, label='样本直方图')
kernel = stats.gaussian_kde(sample)
xvals = np.linspace(35, 53, 1000)
plt.plot(xvals, kernel(xvals), 'r', label='核密度曲线')
plt.xlabel('取值')
plt.ylabel('密度')
plt.legend()
plt.show()

从图中可见样本呈近似对称的钟形，先尝试用截断正态分布拟合原分布。

2. 截断分布的参数估计 #

以正态分布为例，假设原分布为$X \sim N(\mu, \sigma^2)$，样本是$X \geq 35$的条件样本。用Scipy的truncnorm类可自动考虑截断阈值，输出原正态分布的参数：

# 截断正态分布的参数转换：truncnorm的a、b是标准化后的截断阈值
a = (35 - np.mean(sample)) / np.std(sample)  # 初始估计值
b = np.inf  # 仅截断左尾，右边界无限
# 拟合得到原正态分布的参数
mu, sigma = stats.truncnorm.fit(sample, a=a, b=b)[:2]

print(f"原分布均值估计值：{mu:.2f}")
print(f"原分布标准差估计值：{sigma:.2f}")

运行后得到的$\mu$和$\sigma$是原分布的参数，而非截断后样本的均值（样本均值约为41.0，原分布均值会更低，因为左尾被截断）。

若尝试其他分布（比如离散负二项分布、连续伽马分布），可自定义似然函数，用scipy.optimize.minimize最大化似然：

from scipy.optimize import minimize

# 示例：假设原分布是伽马分布，截断在35
def neg_log_likelihood(params):
    shape, scale = params
    # 伽马分布的PDF和生存函数
    pdf = stats.gamma.pdf(sample, a=shape, scale=scale)
    sf = stats.gamma.sf(35, a=shape, scale=scale)  # P(X≥35)
    # 负对数似然
    return -np.sum(np.log(pdf / sf))

# 初始参数
init_params = [20, 2]
result = minimize(neg_log_likelihood, init_params, method='L-BFGS-B')
gamma_shape, gamma_scale = result.x
gamma_mean = gamma_shape * gamma_scale
gamma_std = np.sqrt(gamma_shape) * gamma_scale

print(f"伽马分布均值估计值：{gamma_mean:.2f}")
print(f"伽马分布标准差估计值：{gamma_std:.2f}")

3. 选择最优分布类型 #

通过拟合优度检验比较不同分布的拟合效果，比如Kolmogorov-Smirnov检验：

# 生成截断正态分布的理论分布
trunc_dist = stats.truncnorm(a=(35 - mu)/sigma, b=np.inf, loc=mu, scale=sigma)
# KS检验
stat, p_value = stats.kstest(sample, trunc_dist.cdf)
print(f"KS检验统计量：{stat:.4f}，p值：{p_value:.4f}")

若p值大于0.05，说明不能拒绝“样本来自该截断分布”的假设，即该分布拟合效果较好。

注意事项 #

• 由于仅拥有右尾样本，原分布左半部分信息缺失，参数估计的不确定性比完整样本大，建议用Bootstrap方法生成置信区间。
• 若原分布左尾有特殊形态（比如极端值集中），仅用右尾样本可能导致估计偏差，但这是基于现有数据的最合理方案。

内容的提问来源于stack exchange，提问作者Silvio Duarte