用python将零膨胀未知分布拟合到模型中

Question

我试图找出消费者购买的价值的分布情况。它是零膨胀的，因为大多数消费者在给定的时间内没有进行任何购买。我使用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')