图7.1,出自《统计学习导论》
我目前正在学习一本名为《统计学习导论:R语言应用》的书,并将书中的解决方案转换为Python语言。
我无法理解如何获取置信区间并像上图中显示的那样绘制它们(虚线)。我已经绘制了线条。以下是我的代码 – (我使用的是多项式回归,预测变量是’age’,响应变量是’wage’,次数为4)
poly = PolynomialFeatures(4)X = poly.fit_transform(data['age'].to_frame())y = data['wage']# X.shapemodel = sm.OLS(y,X).fit()print(model.summary())# 所以,我们不仅想要最终的线条,还想要与线条相关的标准误差# 为了找到它,我们需要计算一些年龄值的预测test_ages = np.linspace(data['age'].min(),data['age'].max(),100)X_test = poly.transform(test_ages.reshape(-1,1))pred = model.predict(X_test)plt.figure(figsize = (12,8))plt.scatter(data['age'],data['wage'],facecolors='none', edgecolors='darkgray')plt.plot(test_ages,pred)
这里的数据是R语言中可用的WAGE数据。这是我的结果图 –
回答:
我使用自助法(bootstraping)来计算置信区间,为此我使用了一个自定义模块 –
import numpy as npimport pandas as pdfrom tqdm import tqdmclass Bootstrap_ci: def boot(self,X_data,y_data,R,test_data,model): predictions = [] for i in tqdm(range(R)): predictions.append(self.alpha(X_data,y_data,self.get_indices(X_data,200),test_data,model)) return np.percentile(predictions,2.5,axis = 0),np.percentile(predictions,97.5,axis = 0) def alpha(self,X_data,y_data,index,test_data,model): X = X_data.loc[index] y = y_data.loc[index] lr = model lr.fit(pd.DataFrame(X),y) return lr.predict(pd.DataFrame(test_data)) def get_indices(self,data,num_samples): return np.random.choice(data.index, num_samples, replace=True)
上述模块的使用方法如下 –
poly = PolynomialFeatures(4)X = poly.fit_transform(data['age'].to_frame())y = data['wage']X_test = np.linspace(min(data['age']),max(data['age']),100)X_test_poly = poly.transform(X_test.reshape(-1,1))from bootstrap import Bootstrap_cibootstrap = Bootstrap_ci()li,ui = bootstrap.boot(pd.DataFrame(X),y,1000,X_test_poly,LinearRegression())
这将给我们提供置信区间的下限和上限。要绘制图形 –
plt.scatter(data['age'],data['wage'],facecolors='none', edgecolors='darkgray')plt.plot(X_test,pred,label = '拟合线')plt.plot(X_test,ui,linestyle = 'dashed',color = 'r',label = '置信区间')plt.plot(X_test,li,linestyle = 'dashed',color = 'r')
结果图如下