首页 > 解决方案 > 将零膨胀的未知分布拟合到 python 中的模型

问题描述

我试图找出消费者购买的价值分布是什么。它是零膨胀的,因为大多数消费者不会在给定的时间限制内进行任何购买。我用蟒蛇。由于它是所购买物品的价值,因此我的数据不会像泊松分布那样被亵渎,而是始终非负且连续,这可能意味着对数正态分布、指数分布、伽玛分布、逆伽玛分布等

我的问题归结为如何将分布拟合到零膨胀数据并检查哪个更适合

我找到了很多关于如何进行零膨胀泊松回归的信息,但我的目标是找出过程的基本分布是什么,而不是做出预测,因为我想知道方差。

什么是未知的

  1. 膨胀零的概率 - 并非所有零都膨胀,因为它们也可能是基础分布的结果
  2. 产生购买价值的分布族是什么
  3. 产生购买价值的分布参数是什么

我创建了一个示例代码来生成示例数据并尝试拟合两个分布。不幸的是,真正的 SSE 高于替代方案。

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


N = 1000 * 1000
p_of_inflated_zeros = 0.20

#generation of data
Data = pd.DataFrame({"Prob_bought" : np.random.uniform(0, 1, N) })
Data["If_bought"] = np.where(Data["Prob_bought"] > p_of_inflated_zeros , 1 , 0)
Data["Hipotetical_purchase_value"] = scipy.stats.expon.rvs(scale = 50, loc = 10, size = N) 
#Data["Hipotetical_purchase_value"] = scipy.stats.lognorm.rvs(s = 1, scale = 50, loc = 10, size = N)
Data["Hipotetical_purchase_value"] = np.where(Data["Hipotetical_purchase_value"] < 0 ,0 , Data["Hipotetical_purchase_value"]) 
Data["Purchase_value"] = Data["If_bought"]  * Data["Hipotetical_purchase_value"] 

# fit distribiution
# based on https://stackoverflow.com/questions/6620471/fitting-empirical-distribution-to-theoretical-ones-with-scipy-python
#create 
#x = np.linspace(min(gr_df_trans_tmp), max(gr_df_trans_tmp), 200)
y, x = np.histogram(Data["Purchase_value"], bins = 1000, density = True)
x = (x + np.roll(x, -1))[:-1] / 2.0

#lognormal
FIT_lognorm_sape, FIT_lognorm_loc, FIT_lognorm_scale = scipy.stats.lognorm.fit(Data["Purchase_value"])  
FIT_lognorm_pdf = scipy.stats.lognorm.pdf(x, s = FIT_lognorm_sape, loc = FIT_lognorm_loc, scale = FIT_lognorm_scale)
SSE_lognorm = np.sum(np.power(y - FIT_lognorm_pdf, 2.0))
print(SSE_lognorm)
# 0.036408827144038584

#exponental
FIT_expo_loc, FIT_expo_scale = scipy.stats.expon.fit(Data["Purchase_value"])  
FIT_expo_pdf = scipy.stats.expon.pdf(x, FIT_expo_loc, FIT_expo_scale)
SSE_expo = np.sum(np.power(y - FIT_expo_pdf, 2.0))
print(SSE_expo)
# 0.07564960702319487

# chart
# wykres histogram
axes = plt.gca()
axes.set_xlim([-2, 200])
plt.hist(Data["Purchase_value"], bins = 1000, alpha = 1, density = True)
    
# Plot the PDFs
plt.plot(x, FIT_lognorm_pdf, 'k', linewidth = 1, alpha = 0.5, color = 'red', label = 'lognormal')  
plt.plot(x, FIT_expo_pdf,    'k', linewidth = 1, alpha = 0.5, color = 'blue', label = 'exponental')   
plt.legend(loc='upper right', title = "")

plt.title("Fitting distribiution to ilustrativ data")
plt.xlabel("Hipotetical purchase value")
plt.ylabel('Density')

说明性数据的直方图

标签: pythonstatisticsdistribution

解决方案


推荐阅读