下載數據
import os import tarfile # 用於壓縮和解壓文件 import urllib DOWNLOAD_ROOT = "https://raw.githubusercontent.com/ageron/handson-ml/master/" HOUSING_PATH = "datasets/housing" HOUSING_URL = DOWNLOAD_ROOT + HOUSING_PATH + "/housing.tgz" # 下載數據 def fetch_housing_data(housing_url=HOUSING_URL, housing_path=HOUSING_PATH): if not os.path.isdir(housing_path): os.makedirs(housing_path) tgz_path = os.path.join(housing_path, "housing.tgz") # urlretrieve()方法直接將遠程數據下載到本地 urllib.request.urlretrieve(housing_url, tgz_path) housing_tgz = tarfile.open(tgz_path) housing_tgz.extractall(path=housing_path) # 解壓文件到指定路徑,不指定就是解壓到當前路徑 housing_tgz.close() fetch_housing_data()
加載數據
import pandas as pd def load_housing_data(housing_path=HOUSING_PATH + "/"): csv_path = os.path.join(housing_path, "housing.csv") return pd.read_csv(csv_path)
housing = load_housing_data()
housing.head()
查看數據結構
info()
info()方法可以快速查看數據的描述,特別是總行數、每個屬性的類型和非空值的數量
housing.describe()
housing.info() # 分析:數據集中共有 20640 個實例,按照機器學習的標准這個數據量很小,但是非常適合入門。 # 我們注意到總卧室數只有 20433 個非空值,這意味着有 207 個街區缺少這個值。我們將在后面對它進行處理。
value_counts()
所有的屬性都是數值的,除了離大海距離這項。它的類型是對象,因此可以包含任意 Python 對象,但是因為該項是從 CSV 文件加載的,所以必然是文本類型。在剛才查看數據前五項時,你可能注意到那一列的值是重復的,意味着它可能是一項表示類別的屬性。可以使用value_counts()
方法查看該項中都有哪些類別,每個類別中都包含有多少個街區:
housing["ocean_proximity"].value_counts()
describe()
describe()方法展示了數值屬性的概括
housing.describe()
圖形描述
使用matplotlib的hist()將屬性值畫成柱狀圖,更直觀
import matplotlib.pyplot as plt housing.hist(bins=50, figsize=(10,10)) plt.show() # 不是必要的
- 房屋年齡中位數和房屋價值中位數也被設了上限,因此圖中末尾為一條直線。這種情況解決辦法有兩種
- 1是對於被設置了上線的數據重新收集
- 2是將這些數據從訓練集中移除
- 有些柱狀圖尾巴很長,離中位數過遠。這會使得檢測規律變難,因此會后面后嘗試變換屬性使其變為正太分布。
創建測試集
在這個階段就要分割數據。如果你查看了測試集,就會不經意地按照測試集中的規律來選擇某個特定的機器學習模型。再當你使用測試集來評估誤差率時,就會導致評估過於樂觀,而實際部署的系統表現就會差。這稱為數據透視偏差。
下面3種切分方法:
1.下面的方法,再次運行程序,就會產生一個不同的測試集。解決的辦法之一是保存第一次運行得到的測試集,並在隨后的過程加載。另一種方法是在調用np.random.permutation()
之前,設置隨機數生成器的種子(比如np.random.seed(42)
),以產生總是相同的洗牌指數(shuffled indices)
但是仍舊不完美
import numpy as np def split_train_test(data, test_ratio): shuffled_indices = np.random.permutation(len(data)) # permutation中文排列,輸入數字x,將x以內的數字隨機打散 test_set_size = int(len(data)*test_ratio) test_indices = shuffled_indices[:test_set_size] train_indices = shuffled_indices[test_set_size:] return data.iloc[train_indices], data.iloc[test_indices] train_set, test_set = split_train_test(housing, 0.2) print(len(train_set), "train +", len(test_set), "test")
2.通過實例的哈希值切分
import hashlib def test_set_check(identifier, test_ratio, hash): return hash(np.int64(identifier)).digest()[-1] < 256 * test_ratio def split_train_test_by_id(data, test_ratio, id_column, hash=hashlib.md5): ids = data[id_column] in_test_set = ids.apply(lambda id_: test_set_check(id_, test_ratio, hash)) return data.loc[~in_test_set], data.loc[in_test_set] housing_with_id = housing.reset_index() # adds an `index` column train_set, test_set = split_train_test_by_id(housing_with_id, 0.2, "index") housing_with_id["id"] = housing["longitude"] * 1000 + housing["latitude"] train_set, test_set = split_train_test_by_id(housing_with_id, 0.2, "id") print(len(train_set), "train +", len(test_set), "test")
3.sklearn切分函數
Scikit-Learn 提供了一些函數,可以用多種方式將數據集分割成多個子集。最簡單的函數是`train_test_split`,它的作用和之前的函數`split_train_test`很像,並帶有其它一些功能。比如它有一個`random_state`參數,可以設定前面講過的隨機生成器種子。
from sklearn.model_selection import train_test_split train_set, test_set = train_test_split(housing, test_size=0.2, random_state=42) print(len(train_set), "train +", len(test_set), "test")
sklearn切分函數2
train_test_split屬於純隨機采樣,樣本數量大時很適合。但是如果數據集不大,就會出現采樣偏差的風險。進行分層采樣。可以使用 Scikit-Learn 的StratifiedShuffleSplit類
housing["income_cat"] = np.ceil(housing["median_income"] / 1.5) # ceil對值舍入(以產生離散的分類)除以1.5是為了限制收入分類的數量 housing["income_cat"].where(housing["income_cat"] < 5, 5.0, inplace=True) # 將所有大於 5的分類歸入到分類 5 from sklearn.model_selection import StratifiedShuffleSplit split = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=42) for train_index, test_index in split.split(housing, housing["income_cat"]): strat_train_set = housing.loc[train_index] strat_test_set = housing.loc[test_index] # 記得剔除`income_cat`屬性 for set in (strat_train_set, strat_test_set): set.drop(["income_cat"], axis=1, inplace=True)
數據探索和可視化、發現規律
查看過數據了,現在需要對數據進行了解。
只研究訓練集,如果訓練集非常大的時候還需要再弄一個探索集用來加快運行速度。(本例子不需要)
創建一個副本,以免損傷訓練集
housing = strat_train_set.copy()
地理數據可視化
housing.plot(kind="scatter", x="longitude", y="latitude", alpha=0.4, s=housing["population"]/100, label="population", c="median_house_value", cmap=plt.get_cmap("jet"), colorbar=True,) # 每個圈的半徑表示街區的人口(選項`s`),顏色代表價格(選項`c`) plt.legend()
查找關聯
相關系數1
數據集不是很大時,可以很容易地使用corr()
方法計算出每對屬性間的標准相關系數(standard correlation coefficient,也稱作皮爾遜相關系數)
corr_matrix = housing.corr() corr_matrix["median_house_value"].sort_values(ascending=False)
相關系數2
scatter_matrix
函數,它能畫出每個數值屬性對每個其它數值屬性的圖。
attributes = ["median_house_value", "median_income", "total_rooms", "housing_median_age"] pd.plotting.scatter_matrix(housing[attributes], figsize=(12, 8))
屬性組合實驗
某些屬性本身沒有用,與其它屬性結合起來就有用了。比如下面房間數與房主數本身沒有用,相除得出的每戶的房間數更有用。所以實際操作中需要各種屬性組合,然后比較相關系數。
housing["rooms_per_household"] = housing["total_rooms"]/housing["households"] housing["bedrooms_per_room"] = housing["total_bedrooms"]/housing["total_rooms"] housing["population_per_household"]=housing["population"]/housing["households"] corr_matrix = housing.corr() corr_matrix["median_house_value"].sort_values(ascending=False)
不要手工來做,需要一些函數。原因:
- 函數可以在任何數據集上方便重復地數據轉換
- 慢慢建立一個函數庫,在未來的項目中重復使用
- 可以方便嘗試多種數據轉換
第一步,將特征和標簽分開
housing = strat_train_set.drop("median_house_value", axis=1) housing_labels = strat_train_set["median_house_value"]
數據清洗
缺失值(本例使用total_bedrooms
):
- 去掉對應的街區
dropna()
- 去掉整個屬性
drop()
- 進行賦值(0,平均值,中位數等)
fillna()
使用該方法時記得保存平均值或者中位值等,后面測試集也要填充
housing.dropna(subset=["total_bedrooms"]) # 選項1 housing.drop("total_bedrooms", axis=1) # 選項2 median = housing["total_bedrooms"].median() housing["total_bedrooms"].fillna(median) # 選項3
scikit-learn提供了一個類來處理缺失值:Imputer
from sklearn.preprocessing import Imputer imputer = Imputer(strategy="median") # 因為只有數值屬性才能算出中位數,我們需要創建一份不包括文本屬性`ocean_proximity`的數據副本 housing_num = housing.drop("ocean_proximity", axis=1) # 用`fit()`方法將`imputer`實例擬合到訓練數據 imputer.fit(housing_num)
X = imputer.transform(housing_num)
housing_tr = pd.DataFrame(X, columns=housing_num.columns)
Scikit-Learn 設計
Scikit-Learn 設計的 API 設計的非常好。它的主要設計原則是:
-
一致性:所有對象的接口一致且簡單:
- 估計器(estimator)。任何可以基於數據集對一些參數進行估計的對象都被稱為估計器(比如,
imputer
就是個估計器)。估計本身是通過fit()
方法,只需要一個數據集作為參數(對於監督學習算法,需要兩個數據集;第二個數據集包含標簽)。任何其它用來指導估計過程的參數都被當做超參數(比如imputer
的strategy
),並且超參數要被設置成實例變量(通常通過構造器參數設置)。 - 轉換器(transformer)。一些估計器(比如
imputer
)也可以轉換數據集,這些估計器被稱為轉換器。API也是相當簡單:轉換是通過transform()
方法,被轉換的數據集作為參數。返回的是經過轉換的數據集。轉換過程依賴學習到的參數,比如imputer
的例子。所有的轉換都有一個便捷的方法fit_transform()
,等同於調用fit()
再transform()
(但有時fit_transform()
經過優化,運行的更快)。 - 預測器(predictor)。最后,一些估計器可以根據給出的數據集做預測,這些估計器稱為預測器。例如,上一章的
LinearRegression
模型就是一個預測器:它根據一個國家的人均 GDP 預測生活滿意度。預測器有一個predict()
方法,可以用新實例的數據集做出相應的預測。預測器還有一個score()
方法,可用於評估測試集(如果是監督學習算法的話,還要給出相應的標簽)的預測質量。
- 估計器(estimator)。任何可以基於數據集對一些參數進行估計的對象都被稱為估計器(比如,
-
可檢驗。所有估計器的超參數都可以通過實例的public變量直接訪問(比如,
imputer.strategy
),並且所有估計器學習到的參數也可以通過在實例變量名后加下划線來訪問(比如,imputer.statistics_
)。 -
類不可擴散。數據集被表示成 NumPy 數組或 SciPy 稀疏矩陣,而不是自制的類。超參數只是普通的 Python 字符串或數字。
-
可組合。盡可能使用現存的模塊。例如,用任意的轉換器序列加上一個估計器,就可以做成一個流水線,后面會看到例子。
-
合理的默認值。Scikit-Learn 給大多數參數提供了合理的默認值,很容易就能創建一個系統。
處理文本和類別屬性
from sklearn.preprocessing import LabelEncoder encoder = LabelEncoder() housing_cat = housing["ocean_proximity"] housing_cat_encoded = encoder.fit_transform(housing_cat) housing_cat_encoded
OneHotEncoder
原理是創建一個二元屬性,當分類是<1H OCEAN
,該屬性為 1(否則為 0),當分類是INLAND
,另一個屬性等於 1(否則為 0),以此類推。這稱作獨熱編碼(One-Hot Encoding),因為只有一個屬性會等於 1(熱),其余會是 0(冷)
注意:fit_transform()
用於 2D 數組,而
housing_cat_encoded`是一個 1D 數組,所以需要將其變形
from sklearn.preprocessing import OneHotEncoder encoder = OneHotEncoder() housing_cat_1hot = encoder.fit_transform(housing_cat_encoded.reshape(-1,1)) # housing_cat_1hot結果是一個稀疏矩陣,即只存儲非零項。這是為了當分類很多時節省內存。 將其轉換成numpy數組需要用到toarray函數 housing_cat_1hot.toarray()
LabelBinarizer
使用類LabelBinarizer
,可以用一步執行這兩個轉換(從文本分類到整數分類,再從整數分類到獨熱向量)
from sklearn.preprocessing import LabelBinarizer encoder = LabelBinarizer() housing_cat_1hot = encoder.fit_transform(housing_cat) housing_cat_1hot
注意默認返回的結果是一個密集 NumPy 數組。向構造器LabelBinarizer
傳遞sparse_output=True
,就可以得到一個稀疏矩陣。
自定義轉換器
盡管 Scikit-Learn 提供了許多有用的轉換器,你還是需要自己動手寫轉換器執行任務,比如自定義的清理操作,或屬性組合。你需要讓自制的轉換器與 Scikit-Learn 組件(比如流水線)無縫銜接工作,因為 Scikit-Learn 是依賴鴨子類型的(而不是繼承),你所需要做的是創建一個類並執行三個方法:fit()
(返回self
),transform()
,和fit_transform()
。通過添加TransformerMixin
作為基類,可以很容易地得到最后一個。另外,如果你添加BaseEstimator
作為基類(且構造器中避免使用*args
和**kargs
),你就能得到兩個額外的方法(get_params()
和set_params()
),二者可以方便地進行超參數自動微調。例如,一個小轉換器類添加了上面討論的屬性:
from sklearn.base import BaseEstimator, TransformerMixin rooms_ix, bedrooms_ix, population_ix, household_ix = 3, 4, 5, 6 class CombinedAttributesAdder(BaseEstimator, TransformerMixin): def __init__(self, add_bedrooms_per_room = True): # no *args or **kargs self.add_bedrooms_per_room = add_bedrooms_per_room def fit(self, X, y=None): return self # nothing else to do def transform(self, X, y=None): rooms_per_household = X[:, rooms_ix] / X[:, household_ix] population_per_household = X[:, population_ix] / X[:, household_ix] if self.add_bedrooms_per_room: bedrooms_per_room = X[:, bedrooms_ix] / X[:, rooms_ix] return np.c_[X, rooms_per_household, population_per_household, bedrooms_per_room] else: return np.c_[X, rooms_per_household, population_per_household] attr_adder = CombinedAttributesAdder(add_bedrooms_per_room=False) housing_extra_attribs = attr_adder.transform(housing.values)
特征縮放
屬性的量度相差很大時,模型性能會很差。比如某個屬性范圍是0-1,另外一個屬性的值是10000-50000.這種情況就需要進行特征縮放。有兩種方法:
- 線性函數歸一化:cikit-Learn 提供了一個轉換器
MinMaxScaler
來實現這個功能。 - 標准化:Scikit-Learn 提供了一個轉換器
StandardScaler
來進行標准化。警告:縮放器只用於訓練集擬合。
轉換流水線
由上可以看出,很多轉換步驟需要按照一定的數據。Scikit-Learn 提供了類Pipeline
,來進行這一系列的轉換。
只對數值的流水線
from sklearn.pipeline import Pipeline from sklearn.preprocessing import StandardScaler num_pipeline = Pipeline([ ('imputer', Imputer(strategy="median")), ('attribs_adder', CombinedAttributesAdder()), ('std_scaler', StandardScaler()), ]) housing_num_tr = num_pipeline.fit_transform(housing_num)
scikit-learn中沒有工具來處理Pandas的DataFrame,所以我們需要來寫一個簡單的自定義轉換器來做這項工作
from sklearn.base import BaseEstimator, TransformerMixin class DataFrameSelector(BaseEstimator, TransformerMixin): def __init__(self, attribute_names): self.attribute_names = attribute_names def fit(self, X, y=None): return self def transform(self, X): return X[self.attribute_names].values
一個完整的處理數值以及類別屬性的流水線
# 由於sklearn更新使得LabelBinarizer的fit_transform只能接受兩個參數,直接運行的話會報錯,所以重寫一個轉化器,增加一個參數。 from sklearn.base import TransformerMixin #gives fit_transform method for free class MyLabelBinarizer(TransformerMixin): def __init__(self, *args, **kwargs): self.encoder = LabelBinarizer(*args, **kwargs) def fit(self, x, y=0): self.encoder.fit(x) return self def transform(self, x, y=0): return self.encoder.transform(x)
from sklearn.pipeline import FeatureUnion num_attribs = list(housing_num) cat_attribs = ["ocean_proximity"] num_pipeline = Pipeline([ ('selector', DataFrameSelector(num_attribs)), ('imputer', Imputer(strategy="median")), ('attribs_adder', CombinedAttributesAdder()), ('std_scaler', StandardScaler()), ]) cat_pipeline = Pipeline([ ('selector', DataFrameSelector(cat_attribs)), ('Mylabel_binarizer', MyLabelBinarizer()), ]) full_pipeline = FeatureUnion(transformer_list=[ ("num_pipeline", num_pipeline), ("cat_pipeline", cat_pipeline), ])
housing_prepared = full_pipeline.fit_transform(housing)
housing_prepared[0]
選擇並訓練模型
線性模型使用與評估
from sklearn.linear_model import LinearRegression lin_reg = LinearRegression() lin_reg.fit(housing_prepared, housing_labels)
#結果0.6558010255907188
評估方式1:RMSE
用mean_squared_error`函數計算一下RMSE(歐幾里得范數的平方根的和的根)
結果68628萬美元的誤差顯然不能讓人滿意,欠擬合
from sklearn.metrics import mean_squared_error housing_predictions = lin_reg.predict(housing_prepared) lin_mse = mean_squared_error(housing_labels, housing_predictions) lin_rmse = np.sqrt(lin_mse) lin_rmse
評估方式2:交叉驗證法
將數據集分成k個大小相似的互斥子集,每個子集保持數據分布的一致。然后,每次用k-1個子集的並集作為訓練集,剩下的一個作為測試集,最終返回k個測試結果的均值。k最長用的取值是10,另外5和20也比較常用。
from sklearn.model_selection import cross_val_score lin_scores = cross_val_score(lin_reg, housing_prepared, housing_labels,scoring="neg_mean_squared_error", cv=10) lin_rmse_scores = np.sqrt(-lin_scores) # lin_scores越大越好,取負值開方后的lin_rmse_scores越小越好
def display_scores(scores): print("Scores:", scores) print("Mean:", scores.mean()) print("Standard deviation:", scores.std())
display_scores(lin_rmse_scores) # Mean就是RMSE
決策樹
也可以換其它的模型試試,下面換DecisionTreeRegressor模型
from sklearn.tree import DecisionTreeRegressor tree_reg = DecisionTreeRegressor() tree_reg.fit(housing_prepared, housing_labels)
評估方式1:RMSE
housing_predictions = tree_reg.predict(housing_prepared) tree_mse = mean_squared_error(housing_labels, housing_predictions) tree_rmse = np.sqrt(tree_mse) tree_rmse
#結果 0.0
評估方式2:交叉驗證
scores = cross_val_score(tree_reg, housing_prepared, housing_labels, scoring="neg_mean_squared_error", cv=10) tree_rmse_scores = np.sqrt(-scores) display_scores(tree_rmse_scores)
隨機森林
RandomForestRegressor
from sklearn.ensemble import RandomForestRegressor forest_reg = RandomForestRegressor() forest_reg.fit(housing_prepared, housing_labels)
評估方式1:RMSE
housing_predictions = forest_reg.predict(housing_prepared) forest_mse = mean_squared_error(housing_labels, housing_predictions) forest_rmse = np.sqrt(forest_mse) forest_rmse
#22088.901578494966
評估方式2:交叉驗證
scores = cross_val_score(forest_reg, housing_prepared, housing_labels, scoring="neg_mean_squared_error", cv=10) forest_rmse_scores = np.sqrt(-scores) display_scores(forest_rmse_scores)
比較下來,隨機森林的評分52779.8955803107最小,性能最好,但是仍舊不夠好。還需要給模型加一些限制,或者更多地訓練數據來提高准確率。
模型微調
找到三五個還不錯的模型,列成一個列表,然后對它們進行微調
網格搜索
微調就是手動調整超參數,但是做起來會非常繁雜。應該用Scikit-Learn 的GridSearchCV
來做這項搜索工作
你所需要做的是告訴GridSearchCV
要試驗有哪些超參數,要試驗什么值,GridSearchCV
就能用交叉驗證試驗所有可能超參數值的組合。
from sklearn.model_selection import GridSearchCV param_grid = [ {'n_estimators': [3, 10, 30], 'max_features': [2, 4, 6, 8]}, {'bootstrap': [False], 'n_estimators': [3, 10], 'max_features': [2, 3, 4]}, ] forest_reg = RandomForestRegressor() grid_search = GridSearchCV(forest_reg, param_grid, cv=5, scoring='neg_mean_squared_error') grid_search.fit(housing_prepared, housing_labels)
grid_search.best_params_
提示:因為 30 是n_estimators
的最大值,你也應該估計更高的值,因為評估的分數可能會隨n_estimators
的增大而持續提升。
你還能直接得到最佳的估計器:
grid_search.best_estimator_
# 得到評估得分 cvres = grid_search.cv_results_ for mean_score, params in zip(cvres["mean_test_score"], cvres["params"]): print(np.sqrt(-mean_score), params)
通過網格搜索,找到'max_features': 8, 'n_estimators': 30下的RMSE為49897,低於全部默認值時的52779。微調成功
隨機搜索
搜索組合較多時,網格搜索就不太適用了。這時最好使用RandomizedSearchCV
它不是嘗試所有可能的組合,而是通過選擇每個超參數的一個隨機值的特定數量的隨機組合。
集成方法
另一種微調系統的方法是將表現最好的模型組合起來。組合(集成)之后的性能通常要比單獨的模型要好
分析最佳模型和它們的誤差
可以指出每個屬性對於做出准確預測的相對重要性
feature_importances = grid_search.best_estimator_.feature_importances_
feature_importances
將重要性分數和屬性名放到一起
extra_attribs = ["rooms_per_hhold", "pop_per_hhold", "bedrooms_per_room"] cat_one_hot_attribs = list(encoder.classes_) attributes = num_attribs + extra_attribs + cat_one_hot_attribs sorted(zip(feature_importances,attributes), reverse=True)
由上可以看出,字段中ISLAND最重要,其它的幾個可以刪掉
用測試集評估系統
final_model = grid_search.best_estimator_ X_test = strat_test_set.drop("median_house_value", axis=1) y_test = strat_test_set["median_house_value"].copy() X_test_prepared = full_pipeline.transform(X_test) final_predictions = final_model.predict(X_test_prepared) final_mse = mean_squared_error(y_test, final_predictions) final_rmse = np.sqrt(final_mse) display_scores(final_rmse)
啟動、監控、維護系統
你還需要編寫監控代碼,以固定間隔檢測系統的實時表現,當發生下降時觸發報警。這對於捕獲突然的系統崩潰和性能下降十分重要。做監控很常見,是因為模型會隨着數據的演化而性能下降,除非模型用新數據定期訓練。
評估系統的表現需要對預測值采樣並進行評估。這通常需要人來分析。分析者可能是領域專家,或者是眾包平台(比如 Amazon Mechanical Turk 或 CrowdFlower)的工人。不管采用哪種方法,你都需要將人工評估的流水線植入系統。
你還要評估系統輸入數據的質量。有時因為低質量的信號(比如失靈的傳感器發送隨機值,或另一個團隊的輸出停滯),系統的表現會逐漸變差,但可能需要一段時間,系統的表現才能下降到一定程度,觸發警報。如果監測了系統的輸入,你就可能盡量早的發現問題。對於線上學習系統,監測輸入數據是非常重要的。
最后,你可能想定期用新數據訓練模型。你應該盡可能自動化這個過程。如果不這么做,非常有可能你需要每隔至少六個月更新模型,系統的表現就會產生嚴重波動。如果你的系統是一個線上學習系統,你需要定期保存系統狀態快照,好能方便地回滾到之前的工作狀態。