泊松噪声污染图像的直方图泊松拟合问题及报错求助

问题背景 #

我用泊松噪声添加代码给样本图像添加泊松噪声（图1为污染后图像），选取图像中红色框选的同质ROI分析强度值分布，预期该分布符合泊松分布，但遇到以下三个问题：

问题1：拟合曲线与直方图匹配效果差 #

用以下代码拟合泊松曲线时，曲线和直方图匹配效果不好（图2），但直接对np.random.poisson生成的样本数据拟合时，曲线匹配正常（图3及对应代码）。

# 读取图像
clean_img = cv2.imread("my_image.png") 
clean_gray = cv2.cvtColor(clean_img, cv2.COLOR_BGR2GRAY)
# 转换到0-1范围
clean_gray = (clean_gray-clean_gray.min())/(clean_gray.max()-clean_gray.min())

# 用泊松过程污染图像并转换回0-255范围
peak = 10
noisy_gray = np.random.poisson(clean_gray*peak)/peak
noisy_gray = (((noisy_gray-noisy_gray.min())/(noisy_gray.max()-noisy_gray.min()))*255).astype(np.uint8)

# 定义拟合函数
def fit_function(x, lambda_p):
         # 确保x是整数数组
        pois_pmf = (np.exp(-lambda_p)* (lambda_p ** x) )/ factorial(x+0.00001)
        # print(pois_pmf.dtype)
        return pois_pmf

def calculate_fit_curve(noisy_roi):
    # 绘制选中ROI的直方图
    plt.figure(figsize=(6, 4))
    bins = np.arange(150) - 0.05
    count, bin_edges, _ = plt.hist(noisy_roi.ravel(), bins=bins, density=True, facecolor='steelblue', label='Data')    
    x = ((bin_edges[:-1] + bin_edges[1:]) / 2)

    # 用curve_fit将样本数据拟合成泊松分布
    pois_params, _ = curve_fit(fit_function, x, count, maxfev=1000)
    pois_opt = fit_function(x, *pois_params)
    plt.plot(x, pois_opt,label='Poisson', color='r')
    plt.xlabel('强度值')
    plt.ylabel('归一化密度')
    plt.title('带拟合分布的直方图')
    plt.legend()
    plt.show()

问题2：全强度范围拟合时出现数值溢出，曲线严重偏离 #

增大bins数量（比如设为bins = np.arange(171) - 0.05），尝试绘制0-255全强度范围时，代码行pois_params, _ = curve_fit(fit_function, x, count, maxfev=1000)抛出“overflow encountered in power”和“invalid value encountered in true_divide”警告，拟合曲线严重偏离直方图（图4）。

问题3：增大peak值时拟合失败 #

增大peak值（比如设为30）时，出现“Optimal parameters not found: Number of calls to function has reached maxfev = 1000”报错，即使调大maxfev至10000或取消限制，仍会出现数值溢出警告。

针对问题1：拟合匹配差 #

核心问题是噪声后的归一化步骤破坏了泊松分布特性——你对噪声图像再次做了min-max归一化，拉伸/压缩了原始泊松分布的形态，导致拟合失效。

修正步骤：

1. 去掉噪声生成后的归一化，直接将泊松输出转换为0-255范围：

# 替换原噪声生成代码
peak = 10
noisy_gray = np.random.poisson(clean_gray * peak) * (255 / peak)
noisy_gray = noisy_gray.astype(np.uint8)

2. 拟合时用bin的左边界整数（泊松分布变量为整数），替换原浮点数中点：

# 在calculate_fit_curve函数中修改x的计算
x = bin_edges[:-1].astype(int)

针对问题2：数值溢出与曲线偏离 #

自定义拟合函数的幂运算会在数值较大时溢出，改用scipy内置的泊松PMF函数（已做数值稳定处理），同时过滤空bins避免无效拟合：

from scipy.stats import poisson

# 替换自定义拟合函数
def fit_function(x, lambda_p):
    return poisson.pmf(x, lambda_p)

# 在calculate_fit_curve函数中修改拟合逻辑
def calculate_fit_curve(noisy_roi):
    plt.figure(figsize=(6, 4))
    bins = np.arange(256) - 0.05
    count, bin_edges, _ = plt.hist(noisy_roi.ravel(), bins=bins, density=True, facecolor='steelblue', label='Data')    
    x = bin_edges[:-1].astype(int)

    # 过滤空bins，只拟合有数据的部分
    valid_idx = count > 0
    x_valid = x[valid_idx]
    count_valid = count[valid_idx]

    pois_params, _ = curve_fit(fit_function, x_valid, count_valid, maxfev=1000)
    pois_opt = np.zeros_like(x)
    pois_opt[valid_idx] = fit_function(x_valid, *pois_params)
    
    plt.plot(x, pois_opt,label='Poisson', color='r')
    plt.xlabel('强度值')
    plt.ylabel('归一化密度')
    plt.title('带拟合分布的直方图')
    plt.legend()
    plt.show()

针对问题3：peak值增大后拟合失败 #

peak值增大后lambda参数也会变大，自定义函数的数值稳定性极差，改用scipy内置的poisson.fit直接估计参数，效率和稳定性更高：

from scipy.stats import poisson

def calculate_fit_curve(noisy_roi):
    plt.figure(figsize=(6, 4))
    bins = np.arange(256) - 0.05
    count, bin_edges, _ = plt.hist(noisy_roi.ravel(), bins=bins, density=True, facecolor='steelblue', label='Data')    
    x = bin_edges[:-1].astype(int)

    # 直接用poisson.fit估计lambda参数
    lambda_p = poisson.fit(noisy_roi.ravel())[0]
    pois_opt = poisson.pmf(x, lambda_p)
    
    plt.plot(x, pois_opt,label='Poisson', color='r')
    plt.xlabel('强度值')
    plt.ylabel('归一化密度')
    plt.title('带拟合分布的直方图')
    plt.legend()
    plt.show()

内容的提问来源于stack exchange，提问作者Architha Shastry