數據集背景介紹
- 2009年的《紐約市基准法律》要求對建築的能源和水的使用信息進行說明和評分。 涵蓋的建築包括具有單個建築物的總建築面積超過50,000平方英尺（平方英尺），和群建築面積超過100,000平方英尺。指標是由環境保護署的工具ENERGY STAR Portfolio Manager計算的，並且數據由建築物所有者自行報告。(回歸問題)
字段說明
- 目標數據：
- ENERGY STAR Score：指定建築物類型的1到100百分位排名，在投資組合管理器中根據自我報告計算報告年度的能源使用情況。
數據導入
import pandas as pd import numpy as np import matplotlib.pyplot as plt %matplotlib inline import seaborn as sns from sklearn.model_selection import train_test_split data = pd.read_csv('data/Energy_and_Water_Data_Disclosure_for_Local_Law_84_2017__Data_for_Calendar_Year_2016_.csv') data.head(1)
查看數據形狀
# 查看數據形狀 data.shape # (11746, 60)
查看數據字段類型和是否存在缺失值
# 查看數據字段類型和是否存在缺失值 data.info()
數據集中的字段中存有Not Available值，該值表示無效值可以視為缺失值被處理
# 首先將"Not Available"替換為 np.nan data = data.replace({'Not Available': np.nan})
數據集有的字段顯示為數值型數據，但是實際類型為str，再將部分數值型數據轉換成float
# 數據集有的字段顯示為數值型數據，但是實際類型為str，再將部分數值型數據轉換成float for col in list(data.columns): if ('ft²' in col or 'kBtu' in col or 'Metric Tons CO2e' in col or 'kWh' in col or 'therms' in col or 'gal' in col or 'Score' in col): data[col] = data[col].astype(float)
通過 describe 和 matplotlib 可視化查看數據的相關統計量（柱狀圖）3
data.describe()
# 通過 describe 和 matplotlib 可視化查看數據的相關統計量（柱狀圖） data_desc = data.describe() # 查看數據描述 cols = data_desc.columns # 取得列縮影 index = data_desc.index[1:] # 去除count行 plt.figure(figsize=(30, 30)) # 控制畫布大小 for i in range(len(cols)): ax = plt.subplot(10,6,i+1) # 繪制10x6的表格，當前數據特征維度為60 ax.set_title(cols[i]) # 設置標題 for j in range(len(index)): plt.bar(index[j], data_desc.loc[index[j], cols[i]]) # 對每個特征繪制describe柱狀圖 plt.show() # Order的圖形比較正常，因為最小值，中位數，最大值是錯落分布，正常分布的，且均值和標准差分布也正常 # DOF Gross Floor Area圖形可能有問題，顯示最大值比其他的值都大很多(離均點,異常值)，如果最大值的數據數量較少，則考慮將其刪除 # 發現：經度，維度特征的std極低，且數值分布特別均勻，說明這倆列特征對結果影響幾乎為0，適當考慮過濾該特征
查看每列的缺失值比例
封裝查看每列缺失值比例函數
查看每列的缺失值比例，估計大家有很多種方法可以做，可以給大家提供一個函數，聰明的同學把他當做自己的工具吧，以后你會常用到的
# 查看缺失值 def missing_values_table(df): # 計算每一列缺失值的個數 mis_val = df.isnull().sum(axis=0) # 計算每列缺失值占該列總數據的百分比 mis_val_percent = 100 * df.isnull().sum(axis=0) / data.shape[0] # 將每一列缺失值的數量和缺失值的百分比級聯到一起,形成一個新的表格 mis_val_table = pd.concat([mis_val, mis_val_percent], axis=1) # 重新給上步表格的列命名 mis_val_table_ren_columns = mis_val_table.rename( columns = {0 : 'Missing Values', 1 : '% of Total Values'}) # 將百分比不為0的行數據根據百分比進行降序排序 mis_val_table_ren_columns = mis_val_table_ren_columns[ mis_val_table_ren_columns.iloc[:,1] != 0].sort_values( '% of Total Values', ascending=False).round(1) # 打印概述 print ("Your selected dataframe has " + str(df.shape[1]) + " columns.\n" "There are " + str(mis_val_table_ren_columns.shape[0]) + " columns that have missing values.") # Return the dataframe with missing information return mis_val_table_ren_columns
調用函數查看缺失值及比例
missing_df = missing_values_table(data); missing_df.head(3) Your selected dataframe has 60 columns. There are 46 columns that have missing values. # 下面缺失比例比較大
設置閾值將缺失比例超過百分之50的列刪除
# 設置閾值將缺失比例超過百分之50的列刪除 # 找出超過閾值的列 missing_df = missing_values_table(data); missing_columns = list(missing_df.loc[missing_df['% of Total Values'] > 50].index) print('We will remove %d columns.' % len(missing_columns)) data = data.drop(columns = list(missing_columns)) Your selected dataframe has 60 columns. There are 46 columns that have missing values. We will remove 11 columns.
中位數填充剩下的空值
# 中位數填充剩下的空值 np.median獲取中位數，如果原始數據存在空值就會返回空nan for x in data.columns: # 去除object類型的列(object列不存在中位數) if str(data[x].dtypes) == 'object': continue if data[x].isnull().sum() > 0: # 取出每列非空元素求得中位數進行填充 data[x] = data[x].fillna(value=np.median(data.loc[~data[x].isnull(),x]))
查看目標數據的分布情況
# 查看目標數據的分布情況 data['ENERGY STAR Score'].hist(bins=20)
plt.figure(figsize=(40,20)) plt.scatter(data['ENERGY STAR Score'].index,data['ENERGY STAR Score'].values) # 由圖看到集中到60多出現一條線，說明數據集中在60多，沒有找到離群數據
data['ENERGY STAR Score'].value_counts().sort_values().tail(1) 65.0 2216 Name: ENERGY STAR Score, dtype: int64
發現分布在65的數據量是最多的，沒有發現離群數據
什么是離群點數據？
存在着一些特別大或者特別小的值，這些可能是離群點或記錄錯誤點，對我們結果會有一些影響的。那我們是需要將離群點數據進行過濾的。
- 離群點：離群點是指一個數據序列中，遠離序列的一般水平的極端大值和極端小值，且這些值會對整個數據的分析產生異常的影響。
- 傳統的過濾方式：
- Q1 - 3 * IQ：Q1為序列中25%的中位數，IQ為Q3-Q1
- Q3 + 3 * IQ：Q3為序列中75%的中位數，IQ為Q3-Q1
- 離群點判定：
- 極小的離群點數據：x < (q1 - 3 * iq)
- 極大的離群點數據：x > (q3 + 3 * iq)
假設我們的目標數據為Site EUI (kBtu/ft²)
尋找離群點
data['Site EUI (kBtu/ft²)'].hist(bins=20)
plt.figure(figsize=(15,8)) plt.scatter(data['Site EUI (kBtu/ft²)'].index,data['Site EUI (kBtu/ft²)'].values)
離群點數據過濾
# 離群點數據過濾 q1 = data['Site EUI (kBtu/ft²)'].describe()['25%'] q3 = data['Site EUI (kBtu/ft²)'].describe()['75%'] iq = q3 - q1 # data_copy就是離群的數據 data_copy = data[(data['Site EUI (kBtu/ft²)'] > (q1 - 3 * iq)) & (data['Site EUI (kBtu/ft²)'] < (q3 + 3 * iq))] # 之后我們就可以對離群點做處理，替換replace還是刪除 data_copy['Site EUI (kBtu/ft²)'].hist(bins=30)
plt.scatter(data_copy['Site EUI (kBtu/ft²)'].index,data_copy['Site EUI (kBtu/ft²)'].values)
特征數據的分布
- 可以觀測到特征的取值范圍
- 可以觀測到特征不同數值的分布的密度
- 可以觀測到特征是連續性還是離散型
查看特征的不同取值數量
# 查看特征的不同取值數量 for x in data.columns: print('*'*50) print(x,data[x].nunique()) # 取值少的可能為類別性數據，取值多的為連續性數據
圖形查看所有特征數據分布
for col in data.columns: if 'int' in str(data[col].dtypes) or 'float' in str(data[col].dtypes): # plt.hist(data[col],bins=50) sns.distplot(data.loc[~data[col].isnull(),col]) plt.title(col) plt.show() # 發現有很多特征都是長尾分布的，需要將其轉換為正太或者近正太分布，長尾分布說明特征中少數的數值是離群點數據
長尾分布說明特征中少數的數值是離群點數據，需要將其轉換為正太或者近正太分布了，log操作轉換為近正太分布
log操作轉換為近正太分布
# log操作轉換為近正太分布 data['DOF Gross Floor Area'] sns.distplot(np.log(data.loc[~data['DOF Gross Floor Area'].isnull(),'DOF Gross Floor Area']))
直方圖查看數據分布
# 直方圖 for col in data.columns: if 'int' in str(data[col].dtypes) or 'float' in str(data[col].dtypes): plt.hist(data[col],bins=50) plt.title(col) plt.show()
類別型特征的特征值化
- 結合着項目的目的和對非數值型字段特征的理解等手段，我們只選取出兩個代表性特征Borough和Largest Property Use Type對其進行onehot編碼實現特征值化
features = data.loc[:,data.columns != 'ENERGY STAR Score']# 提取特征數據 fea_name = features.select_dtypes('number').columns # 提取數值型特征名稱 features = features[fea_name] # 提取數值型特征 # 提取指定的兩個類別型特征 categorical_subset = data[['Borough', 'Largest Property Use Type']] categorical_subset = pd.get_dummies(categorical_subset) features = pd.concat([features, categorical_subset], axis = 1)
探索特征之間的相關性
corr查看特征與特征的相關性
features.corr() # corr查看特征與特征的相關性
plt.subplots(figsize=(30,15)) #指定窗口尺寸（單位英尺） features_corr=features.corr().abs() # 返回列與列之間的相關系數 abs求得是絕對值，相關系數與正負無關 # 數據為相關系數，顯示數值，顯示顏色條 這里需要導入模快 import seaborn as sns 也是一個繪圖模塊 sns.heatmap(features_corr, annot=True)
去除相關性強的冗余特征，工具包封裝如下
cols = features.columns # 獲取列的名稱 corr_list = [] size = features.shape[1] high_corr_fea = [] # 存儲相關系數大於0.5的特征名稱 for i in range(0,size): for j in range(i+1, size): if(abs(features_corr.iloc[i,j])>= 0.5): corr_list.append([features_corr.iloc[i,j], i, j]) # features_corr.iloc[i,j]：按位置選取數據 sorted_corr_list = sorted(corr_list, key=lambda xx:-abs(xx[0])) # print(sorted_corr_list) for v,i,j in sorted_corr_list: high_corr_fea.append(cols[i]) print("%s and %s = %.2f" % (cols[i], cols[j],v)) # cols: 列名
刪除相關性強的特征
# 刪除特征 features.drop(labels=high_corr_fea,axis=1,inplace=True) features.shape # (11746, 69)
查看特征和目標之間的相關性
- 如果特征和標簽之間是存在線性關系的才可以采用如下方式
target = data['ENERGY STAR Score'] target = pd.DataFrame(data=target,columns=['ENERGY STAR Score']) # 級聯target&features new_data = pd.concat((features,target),axis=1) # 計算相關性,之后我們就可以選擇特征與目標相關性較大的數據進行特征的選取 fea_target_corr = abs(new_data.corr()['ENERGY STAR Score'][:-1]) fea_target_corr
保存數據
# 改名字 new_data = new_data.rename(columns = {'ENERGY STAR Score': 'score'}) new_data.to_csv('./data/eda_data.csv')
建模 (該案例是回歸問題)
數據集導入及切分
# Pandas and numpy for data manipulation import pandas as pd import numpy as np import matplotlib.pyplot as plt %matplotlib inline import seaborn as sns from sklearn.preprocessing import MinMaxScaler from sklearn.linear_model import LinearRegression from sklearn.ensemble import RandomForestRegressor from sklearn.svm import SVR from sklearn.neighbors import KNeighborsRegressor from lightgbm import LGBMRegressor # Hyperparameter tuning from sklearn.model_selection import GridSearchCV,train_test_split data = pd.read_csv('./data/eda_data.csv').drop(labels='Unnamed: 0',axis=1) fea_name = [x for x in data.columns if x not in ['score']] feature = data[fea_name] target = data['score'] x_train,x_test,y_train,y_test = train_test_split(feature,target,test_size=0.2,random_state=2020)
選擇的機器學習算法（回歸問題）
- Linear Regression
- Support Vector Machine Regression
- Random Forest Regression
- lightGBM (忽略)
- xgboost
- 我們只選擇其默認的參數，這里先不進行調參工作，后續再來調參。
- 由於模型的調用和評價方式是一樣的，則封裝如下工具函數
# Function to calculate mean absolute error def mae(y_true, y_pred): return np.mean(abs(y_true - y_pred)) # Takes in a model, trains the model, and evaluates the model on the test set def fit_and_evaluate(model): # Train the model model.fit(x_train, y_train) # Make predictions and evalute model_pred = model.predict(x_test) model_mae = mae(y_test, model_pred) # Return the performance metric return model_mae
嘗試各種模型
# 線性回歸 lr = LinearRegression() lr_mae = fit_and_evaluate(lr) print('Linear Regression Performance on the test set: MAE = %0.4f' % lr_mae) Linear Regression Performance on the test set: MAE = 20.8585 # svm支向量機 svm = SVR(C = 1000, gamma = 0.1) svm_mae = fit_and_evaluate(svm) print('Support Vector Machine Regression Performance on the test set: MAE = %0.4f' % svm_mae) Support Vector Machine Regression Performance on the test set: MAE = 21.3965 # 隨機森林 random_forest = RandomForestRegressor(random_state=60) random_forest_mae = fit_and_evaluate(random_forest) print('Random Forest Regression Performance on the test set: MAE = %0.4f' % random_forest_mae) Random Forest Regression Performance on the test set: MAE = 11.7919
條形圖顯示mae
plt.style.use('fivethirtyeight') model_comparison = pd.DataFrame({'model': ['Linear Regression', 'Support Vector Machine', 'Random Forest', ], 'mae': [lr_mae, svm_mae, random_forest_mae,]}) model_comparison.sort_values('mae', ascending = False).plot(x = 'model', y = 'mae', kind = 'barh', color = 'red', edgecolor = 'black') plt.ylabel(''); plt.yticks(size = 14); plt.xlabel('Mean Absolute Error'); plt.xticks(size = 14) plt.title('Model Comparison on Test MAE', size = 20);
看起來集成算法更占優勢一些，這里存在一些不公平，因為參數只用了默認，但是對於SVM來說參數可能影響會更大一些。
模型調參 (這里用網格搜索找尋最佳參數)
- 以隨機森林參數調優為例
from sklearn.model_selection import GridSearchCV parameters = { 'n_estimators':[150,200,250] ,'min_samples_split':[2, 3, 5, 10, 15] ,"max_depth":[2, 3, 5, 10, 15] ,'min_samples_leaf':[1, 2, 4, 6, 8] ,'min_samples_split':[2, 4, 6, 10] } model = RandomForestRegressor() GS = GridSearchCV(estimator=model,param_grid=parameters,cv=5,scoring='neg_mean_absolute_error') GS.fit(x_train,y_train) GS.best_params_