非参数化的自助重采样方法在Logistic回归应用及模型差异分析|附数据代码

 本文探讨了计算逻辑回归参数抽样分布的不同方法，包括非参数化的自助重采样方法、参数化的自助方法以及一种混合模式。通过对这些方法的实施和分析，揭示了它们基于不同的建模假设，产生不同的结果。同时指出，在实际应用中，应根据差异的重要性选择最适合的方法，而不是认为某种方法是唯一正确的。

---

---

关键词：重采样；逻辑回归；参数化自助法；非参数化自助法

一、引言

计算抽样分布的不同方法会产生不同结果，但在实践中差异通常较小，我们可以选择方便的方法。本文由 Reddit 上的一个问题引发，探讨了逻辑回归的参数化引导问题中随机误差项的确定，并介绍了两种计算逻辑回归参数抽样分布的方法。本文将实施这些方法并解释其假设，提出混合方法，并给出选择方法的标准。

1. 为了对二次关系进行建模，添加包含AGE和EDUC平方值的列，并删除有缺失值的行。

gss['AGE2'] = gss['AGE']**2gss['EDUC2'] = gss['EDUC']**2data = gss.dropna(subset=['AGE', 'EDUC', 'SEX', 'GUNLAW', 'GRASS'])

三、逻辑回归结果

使用 StatsModels 进行逻辑回归，得到回归结果。绘制男性和女性受访者支持大麻合法化的预测概率随年龄变化的曲线。import statsmodels.formula.api as smfformula = 'GRASS ~ AGE + AGE2 + EDUC + EDUC2 + C(SEX)'result_hat = smf.logit(formula, data=data).fit()result_hat.summary()df = pd.DataFrame()df['AGE'] = np.linspace(18, 89)df['EDUC'] = 16df['AGE2'] = df['AGE']**2df['EDUC2'] = df['EDUC']**2df['SEX'] = 1pred1 = result_hat.predict(df)pred1.index = df['AGE']df['SEX'] = 2pred2 = result_hat.predict(df)pred2.index = df['AGE']pred1.plot(label='Male', alpha=0.6)pred2.plot(label='Female', alpha=0.6)plt.xlabel('Age')plt.ylabel('Fraction')plt.title('Support for legal marijuana')plt.legend();

结果显示男性比女性更有可能支持合法化，年轻人比老年人更有可能支持合法化。 

四、非参数化自助重采样方法

（一）方法描述
基于自助重采样，对data的行进行有放回抽样，并对重新采样的数据运行回归模型。

options = dict(disp=False, start_params=result_hat.params) def bootstrap(i):    bootstrapped = data.sample(n=len(data), replace=True)    results = smf.logit(formula, data=bootstrapped).fit(**options)    return results.params

（二）结果分析
通过多次运行该过程，从抽样分布中生成样本。分析其中一个参数的抽样分布，发现抽样分布的平均值接近使用原始数据集估计的参数，标准差接近 StatsModels 计算的标准误差。

estimates = pqdm(range(101), bootstrap, n_jobs=4)sampling_dist = pd.DataFrame(estimates)ci90 = sampling_dist['C(SEX)[T.2]'].quantile([0.05, 0.95])sns.kdeplot(sampling_dist['C(SEX)[T.2]'])[plt.axvline(x, ls=':') for x in ci90]plt.title('Sampling distribution of a parameter');pd.DataFrame({"Sampling mean": sampling_dist.mean(), "Estimates": result_hat.params})def standard_errors(sampling_dist, result_hat): df = pd.DataFrame({"Sampling std": sampling_dist.std(), "Standard error": result_hat.bse}) num, den = df.values.T df['Percent diff'] = (num / den - 1) * 100 return dfstandard_errors(sampling_dist, result_hat)

五、参数化自助方法

（一）方法描述
假设从原始数据中估计的参数是正确的，使用回归模型计算每个受访者的预测概率，然后用这些概率为每个受访者生成有偏差的抛硬币，将模拟值作为因变量运行回归模型。

pi_hat = result_hat.predict(data)from scipy.stats import bernoullisimulated = bernoulli.rvs(pi_hat.values)def bootstrap2(i):    flipped = data.assign(GRASS=bernoulli.rvs(pi_hat.values))    results = smf.logit(formula, data=flipped).fit(**options)    return results.params

（二）结果分析
多次运行该方法得到抽样分布样本，分析参数的抽样分布，其标准差也接近标准误差。

estimates = pqdm(range(101), bootstrap2, n_jobs=4)sampling_dist2 = pd.DataFrame(estimates)ci90 = sampling_dist2['C(SEX)[T.2]'].quantile([0.05, 0.95])sns.kdeplot(sampling_dist2['C(SEX)[T.2]'])[plt.axvline(x, ls=':') for x in ci90]plt.title('Sampling distribution of a parameter');standard_errors(sampling_dist2, result_hat)

六、混合模式方法

（一）方法描述
结合非参数化和参数化方法，先使用自助抽样模拟调查不同的人，再用参数化方法模拟他们的响应。

def bootstrap3(i):    bootstrapped = data.sample(n=len(data), replace=True)    pi_hat = result_hat.predict(bootstrapped)    flipped = bootstrapped.assign(GRASS=bernoulli.rvs(pi_hat.values))    results = smf.logit(formula, data=flipped).fit(**options)    return results.params

（二）结果分析
得到抽样分布样本，分析参数的抽样分布，标准误差与 StatsModels 计算的误差相当。

estimates = pqdm(range(101), bootstrap3, n_jobs=4)sampling_dist3 = pd.DataFrame(estimates)ci90 = sampling_dist3['C(SEX)[T.2]'].quantile([0.05, 0.95])sns.kdeplot(sampling_dist3['C(SEX)[T.2]'])[plt.axvline(x, ls=':') for x in ci90]plt.title('Sampling distribution of a parameter');standard_errors(sampling_dist3, result_hat)

七、讨论与结论

（一）不同方法的比较
本文介绍了三种计算抽样分布的方法，它们基于不同的抽样过程模型，产生不同的结果。在某些情况下，不同方法可能渐近收敛于相同结果，但对于有限数据集通常不同。
（二）选择方法的标准
如果差异较小，在实践中可能无关紧要，可以选择最容易实现、计算最快或方便的方法。不能认为分析方法的结果是唯一正确的，它们也基于建模假设和近似值。

综上所述，在实际应用中应根据具体情况选择合适的方法来计算逻辑回归参数的抽样分布。

---

关于分析师

在此对 Anting Li 对本文所作的贡献表示诚挚感谢，她在中央财经大学完成了应用统计学专业的硕士学位，专注统计学领域。擅长 R 语言、Python、Matlab。