# 如何规避线性回归的陷阱（下）

## 使用变量变换或广义线性模型

# Create gamma dataset
np.random.seed(1)
x1 = np.random.uniform(-1, 1, 200)
x2 = np.random.uniform(-1, 1, 200)
mu = np.exp(1 + x1 + 2*x2 + np.random.randn())
y = np.random.gamma(shape = 2, scale = mu/2, size = 200)
gamma_data = pd.DataFrame({'X1':x1, 'X2':x2, 'Y':y})
# Plot data
plt.hist(gamma_data['Y'], bins=30)
plt.ylabel('Count')
plt.xlabel('Y')
plt.show()

# Fit linear regression
non_norm_model = smf.ols(formula='Y ~ X1 + X2',
data=gamma_data).fit()
# Calculate residuals
resid = gamma_data['Y'] - non_norm_model.predict(gamma_data[['X1', 'X2']])
# Plot residuals
plt.scatter(non_norm_model.predict(gamma_data[['X1', 'X2']]), resid, alpha=0.5)
plt.xlabel('fitted')
plt.ylabel('residual')
plt.show()

# Transform Y by taking the log of it
gamma_data.loc[:, 'log_y'] = gamma_data['Y'].apply(lambda x:
np.log(x))
# Refit linear regression to transformed data
non_norm_model2 = smf.ols(formula='log_y ~ X1 + X2',
data=gamma_data).fit()
# Calculate residuals
resid = gamma_data['log_y'] -
non_norm_model2.predict(gamma_data[['X1', 'X2']])
# Plot residuals
plt.scatter(non_norm_model2.predict(gamma_data[['X1', 'X2']]), resid, alpha=0.5)
plt.xlabel('fitted')
plt.ylabel('residual')

# Fit GLM to data
gamma_model = sm.GLM(gamma_data['Y'],
# Calculate residuals
resid = gamma_model.resid_deviance
# Plot residuals
plt.xlabel('fitted')
plt.ylabel('residual')
plt.show()

## 使用时间序列模型处理自相关

# Only keep data for AAL
sandp_data = sandp_data[sandp_data['Name'] == 'AAL']
# Transform date to datetime format
sandp_data['date'] = pd.to_datetime(sandp_data['date'])
# Sort data by date
sandp_data.sort_values(by = ['date'], inplace = True)
# Plot data
plt.plot(sandp_data['date'], sandp_data['close'])
plt.show()

# Create year and month variables
sandp_data['Year'] = sandp_data['date'].map(lambda x: x.year)
sandp_data['Month'] = sandp_data['date'].map(lambda x: x.month)
# Fit linear regression to data
sandp_model = smf.ols(formula='close ~ Year + C(Month)',
data=sandp_data).fit()# Plot data with regression line overlay
plt.plot(sandp_data['date'], sandp_data['close'])
plt.plot(sandp_data['date'], sandp_model.predict(sandp_data[['Year', 'Month']]), color = 'red')
plt.show()

##### 例如，如果d=1，y（t）=y（t）-y（t-1），如果d=2，y（t）=z（t）-z（t-1），其中z（t）=y（t）-y（t-1），依此类推。

# Index dataset by date
sandp_data2 = sandp_data.set_index('date')
# Fit ARIMA model to the data
ts_model = sm.tsa.ARIMA(sandp_data2['close'], (5, 1, 0)).fit()
# Plot actual vs fitted values
ts_model.plot_predict(dynamic=False)

Occam剃刀原理指出，如果“对某一事件存在两种解释，通常需要最少猜测的解释是正确的”。

0 条评论

• ### 数据科学家常犯的十大编程错误

数据科学家是“比任何软件工程师都更擅长统计，比任何软件工程师都更擅长软件工程的的统计学家”。许多数据科学家都有统计学背景却缺乏在软件工程方面的经验。我是资深的数...

• ### 使用LSTM预测比特币价格

本文以“时间序列预测的LSTM神经网络”这篇文章为基础。如果没有阅读，我强烈建议你读一读。 考虑到近期对比特币货币的泡沫的讨论，我写了这篇文章，主要是为了预测比...

• ### Python/PyMC3/ArviZ贝叶斯统计实战（上）

如果你认为贝叶斯定理是反直觉的，那么建立在贝叶斯定理基础上的贝叶斯统计就很难理解。在这一点上我和你的感受完全一致。

• ### 超全的pandas数据分析常用函数总结：下篇

基础知识在数据分析中就像是九阳神功，熟练的掌握，加以运用，就可以练就深厚的内力，成为绝顶高手自然不在话下！

• ### 从0到1掌握R语言网络爬虫

引言 网上的数据和信息无穷无尽，如今人人都用百度谷歌来作为获取知识，了解新鲜事物的首要信息源。所有的这些网上的信息都是直接可得的，而为了满足日益增长的数据需求，...

• ### 使用脚手架应用做单元测试

版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。 ...

• ### 上映4天，票房7.4亿的《海王》，用Python分析数据看大片！

《海王》一部电影带你重温《驯龙高手》《变形金刚》《星球大战》《星河战队》《铁血战士》《安德的游戏》《异形》可能还借鉴了对手的《钢铁侠》与《黑豹》剧情，再稍稍带一...

• ### BootstrapTable的列排序怎么搞

先搞一个table，使用ajax将数据查询出来，然后可以在所有列都加上排序。满足自己的需求。

• ### 使用脚手架应用做单元测试

因为后台service比较复杂，需要三个不同的实例协同工作，所以之前Oliver开发了Scaffolding App这个Angular前端，目的是方便我们随时测...