pandas 使用 StatsModels 绘制二阶多项式的分位数回归
声明：本页面是StackOverFlow热门问题的中英对照翻译，遵循CC BY-SA 4.0协议，如果您需要使用它，必须同样遵循CC BY-SA许可，注明原文地址和作者信息，同时你必须将它归于原作者（不是我）：StackOverFlow
原文地址: http://stackoverflow.com/questions/35185046/
Warning: these are provided under cc-by-sa 4.0 license. You are free to use/share it, But you must attribute it to the original authors (not me):
StackOverFlow
Using StatsModels to plot quantile regression for 2nd order polynomial
提问by PJW
I am following the StatsModels example hereto plot quantile regression lines. With only slight modification for my data, the example works great, producing this plot (note that I have modified the code to only plot the 0.05, 0.25, 0.5, 0.75, and 0.95 quantiles) :
我正在按照此处的 StatsModels 示例绘制分位数回归线。只对我的数据稍作修改，这个例子就很好用，生成了这个图（请注意，我修改了代码，只绘制了 0.05、0.25、0.5、0.75 和 0.95 分位数）：
However, I would like to plot the OLS fit and corresponding quantiles for a 2nd order polynomial fit (instead of linear). For example, here is the 2nd-order OLS line for the same data:
但是，我想为二阶多项式拟合（而不是线性）绘制 OLS 拟合和相应的分位数。例如，这里是相同数据的二阶 OLS 线：
How can I modify the code in the linked example to produce non-linear quantiles?
如何修改链接示例中的代码以生成非线性分位数？
Here is my relevant code modified from the linked example to produce the 1st plot:
这是我从链接示例修改的相关代码以生成第一个图：
d = {'temp': x, 'dens': y}
df = pd.DataFrame(data=d)
# Least Absolute Deviation
#
# The LAD model is a special case of quantile regression where q=0.5
mod = smf.quantreg('dens ~ temp', df)
res = mod.fit(q=.5)
print(res.summary())
# Prepare data for plotting
#
# For convenience, we place the quantile regression results in a Pandas DataFrame, and the OLS results in a dictionary.
quantiles = [.05, .25, .50, .75, .95]
def fit_model(q):
res = mod.fit(q=q)
return [q, res.params['Intercept'], res.params['temp']] + res.conf_int().ix['temp'].tolist()
models = [fit_model(x) for x in quantiles]
models = pd.DataFrame(models, columns=['q', 'a', 'b','lb','ub'])
ols = smf.ols('dens ~ temp', df).fit()
ols_ci = ols.conf_int().ix['temp'].tolist()
ols = dict(a = ols.params['Intercept'],
b = ols.params['temp'],
lb = ols_ci[0],
ub = ols_ci[1])
print(models)
print(ols)
x = np.arange(df.temp.min(), df.temp.max(), 50)
get_y = lambda a, b: a + b * x
for i in range(models.shape[0]):
y = get_y(models.a[i], models.b[i])
plt.plot(x, y, linestyle='dotted', color='grey')
y = get_y(ols['a'], ols['b'])
plt.plot(x, y, color='red', label='OLS')
plt.scatter(df.temp, df.dens, alpha=.2)
plt.xlim((-10, 40))
plt.ylim((0, 0.4))
plt.legend()
plt.xlabel('temp')
plt.ylabel('dens')
plt.show()
回答by PJW
After a day of looking into this, came up with a solution, so posting my own answer. Much credit to Josef Perktold at StatsModels for assistance.
经过一天的调查，想出了一个解决方案，所以发布我自己的答案。非常感谢 StatsModels 的 Josef Perktold 提供的帮助。
Here is the relevant code and plot:
这是相关的代码和情节：
d = {'temp': x, 'dens': y}
df = pd.DataFrame(data=d)
x1 = pd.DataFrame({'temp': np.linspace(df.temp.min(), df.temp.max(), 200)})
poly_2 = smf.ols(formula='dens ~ 1 + temp + I(temp ** 2.0)', data=df).fit()
plt.plot(x, y, 'o', alpha=0.2)
plt.plot(x1.temp, poly_2.predict(x1), 'r-',
label='2nd order poly fit, $R^2$=%.2f' % poly_2.rsquared,
alpha=0.9)
plt.xlim((-10, 50))
plt.ylim((0, 0.25))
plt.xlabel('mean air temp')
plt.ylabel('density')
plt.legend(loc="upper left")
# with quantile regression
# Least Absolute Deviation
# The LAD model is a special case of quantile regression where q=0.5
mod = smf.quantreg('dens ~ temp + I(temp ** 2.0)', df)
res = mod.fit(q=.5)
print(res.summary())
# Quantile regression for 5 quantiles
quantiles = [.05, .25, .50, .75, .95]
# get all result instances in a list
res_all = [mod.fit(q=q) for q in quantiles]
res_ols = smf.ols('dens ~ temp + I(temp ** 2.0)', df).fit()
plt.figure()
# create x for prediction
x_p = np.linspace(df.temp.min(), df.temp.max(), 50)
df_p = pd.DataFrame({'temp': x_p})
for qm, res in zip(quantiles, res_all):
# get prediction for the model and plot
# here we use a dict which works the same way as the df in ols
plt.plot(x_p, res.predict({'temp': x_p}), linestyle='--', lw=1,
color='k', label='q=%.2F' % qm, zorder=2)
y_ols_predicted = res_ols.predict(df_p)
plt.plot(x_p, y_ols_predicted, color='red', zorder=1)
#plt.scatter(df.temp, df.dens, alpha=.2)
plt.plot(df.temp, df.dens, 'o', alpha=.2, zorder=0)
plt.xlim((-10, 50))
plt.ylim((0, 0.25))
#plt.legend(loc="upper center")
plt.xlabel('mean air temp')
plt.ylabel('density')
plt.title('')
plt.show()
red line: 2nd order polynomial fit
红线：二阶多项式拟合
black dashed lines: 5th, 25th, 50th, 75th, 95th percentiles
黑色虚线：第 5、第 25、第 50、第 75、第 95 个百分位数
