基於模型的特征選擇詳解 (Embedded & Wrapper)
單變量特征選擇方法獨立的衡量每個特征與響應變量之間的關系，另一種主流的特征選擇方法是基於機器學習模型的方法。有些機器學習方法本身就具有對特征進行打分的機制，或者很容易將其運用到特征選擇任務中，例如回歸模型，SVM，決策樹，隨機森林等等。
1. 線性模型和正則化(Embedded方式)
下面將介紹如何用回歸模型的系數來選擇特征。越是重要的特征在模型中對應的系數就會越大，而跟輸出變量越是無關的特征對應的系數就會越接近於0。在噪音不多的數據上，或者是數據量遠遠大於特征數的數據上，如果特征之間相對來說是比較獨立的，那么即便是運用最簡單的線性回歸模型也一樣能取得非常好的效果。
from sklearn.linear_model import LinearRegression
import numpy as np
# A helper method for pretty-printing linear models
def pretty_print_linear(coefs, names=None, sort=False):
if names == None:
names = ["X%s" % x for x in range(len(coefs))]
lst = zip(coefs, names)
if sort:
lst = sorted(lst, key=lambda x: -np.abs(x[0]))
return " + ".join("%s * %s" % (round(coef, 3), name) for coef, name in lst)
np.random.seed(0) # 有了這段代碼，下次再生成隨機數的時候，與上次一樣的結果
size = 5000 # 表示抽取多少個樣本
# 隨機生成3個特征的樣本，每個維度的特征都是服從期望為0，標准差為1的正態分布
X = np.random.normal(0, 1, (size, 3)) # 抽取5000個樣本，每個樣本都是3維的
# Y = X0 + 2*X1 + noise
Y = X[:, 0] + 2 * X[:, 1] + np.random.normal(0, 2, size)
lr = LinearRegression()
lr.fit(X, Y)
print("Linear model:", pretty_print_linear(lr.coef_))
>>> Linear model: 0.984 * X0 + 1.995 * X1 + -0.041 * X2
在這個例子當中，盡管數據中存在一些噪音，但這種特征選擇模型仍然能夠很好的體現出數據的底層結構。當然這也是因為例子中的這個問題非常適合用線性模型來解：特征和響應變量之間全都是線性關系，並且特征之間均是獨立的。
然而，在很多實際的數據當中，往往存在多個互相關聯的特征，這時候模型就會變得不穩定，對噪聲很敏感，數據中細微的變化就可能導致模型的巨大變化（模型的變化本質上是系數，或者叫參數），這會讓模型的預測變得困難，這種現象也稱為多重共線性。例如，假設我們有個數據集，它的真實模型應該是\(Y=X_1+X_2\)，當我們觀察的時候，發現\(Y’=X_1+X_2+e\)，\(e\)是噪音。如果\(X_1\)和\(X_2\)之間存在線性關系，例如\(X_1\)約等於\(X_2\)，這個時候由於噪音\(e\)的存在，我們學到的模型可能就不是\(Y=X_1+X_2\)了，有可能是\(Y=2X_1\)，或者\(Y=-X_1+3X_2\)。【在\(X_1\)約等於\(X_2\)的線性關系下，學到的這三種模型，理論上來說都是正確的，注意他們有個共同的特點是，系數之和為2】
size = 100
# 另外一個種子為5的隨機采樣，若要執行相同的結果，以種子號來區分隨機采樣的結果
np.random.seed(seed=5)
X_seed = np.random.normal(0, 1, size)
X1 = X_seed + np.random.normal(0, .1, size)
X2 = X_seed + np.random.normal(0, .1, size)
X3 = X_seed + np.random.normal(0, .1, size)
X = np.array([X1, X2, X3]).T
Y = X1 + X2 + X3 + np.random.normal(0, 1, size)
lr = LinearRegression()
lr.fit(X, Y)
print("Linear model:", pretty_print_linear(lr.coef_))
>>> Linear model: -1.291 * X0 + 1.591 * X1 + 2.747 * X2
這個例子和上個例子中，三個特征都來源於標准正態分布的隨機采樣，這里系數之和接近3，基本上和上上個例子的結果一致，應該說學到的模型對於預測來說還是不錯的。但是，如果從系數的字面意思上去解釋特征的重要性的話，X2對於輸出變量來說具有很強的正面影響，而X0具有負面影響，而實際上所有特征與輸出變量之間的影響是均等的。
同樣的方法和套路可以用到類似的線性模型上，比如邏輯回歸。
正則化的概念：
正則化就是把額外的約束或者懲罰項加到已有模型（損失函數）上，以防止過擬合並提高泛化能力。損失函數由原來的\(L(X,Y)\)變為$L(X,Y)+\alpha \left| w \right| \(，\)w$是模型系數組成的向量（有些地方也叫參數parameter，coefficients）， \(\left\| \cdot \right\|\)一般是L1或者L2范數，\(\alpha\)是一個可調的參數，控制着正則化的強度。當用在線性模型上時，L1正則化和L2正則化也稱為Lasso(least absolute shrinkage and selection operator)和Ridge。
1.1 L1正則化(Lasso)
L1正則化將系數\(w\)的l1范數作為懲罰項加到損失函數上，由於正則項非零，這就迫使那些弱的特征所對應的系數變成0。因此L1正則化往往會使學到的模型很稀疏（系數w經常為0），這個特性使得L1正則化成為一種很好的特征選擇方法。
Scikit-learn為線性回歸模型提供了Lasso，為分類模型提供了L1邏輯回歸。下面的例子在波士頓房價數據上運行了Lasso，其中參數alpha是通過grid search進行優化的。
from sklearn.linear_model import Lasso
from sklearn.preprocessing import StandardScaler
from sklearn.datasets import load_boston
import numpy as np
boston = load_boston()
scaler = StandardScaler()
X = scaler.fit_transform(boston["data"])
Y = boston["target"]
names = boston["feature_names"]
lasso = Lasso(alpha=.3)
lasso.fit(X, Y)
print("Lasso model: ", pretty_print_linear(lasso.coef_, names, sort=True))
Lasso model: -3.707 * LSTAT + 2.992 * RM + -1.757 * PTRATIO + -1.081 * DIS + -0.7 * NOX + 0.631 * B + 0.54 * CHAS + -0.236 * CRIM + 0.081 * ZN + -0.0 * INDUS + -0.0 * AGE + 0.0 * RAD + -0.0 * TAX
可以看到，很多特征的系數都是0。如果繼續增加\(\alpha\)的值，得到的模型就會越來越稀疏，即越來越多的特征系數會變成0。
然而，L1正則化像非正則化線性模型一樣也是不穩定的，如果特征集合中具有相關聯的特征，當數據發生細微變化時也有可能導致很大的模型差異。
1.2 L2正則化(Ridge Regression)
L2正則化同樣將系數向量的L2范數添加到了損失函數中。由於L2懲罰項中系數是二次方的，這使得L2和L1有着諸多差異，最明顯的一點就是，L2正則化會讓系數的取值變得平均。對於關聯特征，這意味着他們能夠獲得更相近的對應系數。還是以\(Y=X_1+X_2\)為例，假設\(X_1\)和\(X_2\)具有很強的關聯，如果用L1正則化，不論學到的模型是\(Y=X_1+X_2\)還是\(Y=2X_1\)，懲罰都是一樣的，都是2\(\alpha\)。但是對於L2來說，第一個模型的懲罰項是\(2\alpha\)，但第二個模型的是\(4\alpha\)。可以看出，系數(待求參數)之和為常數時，各系數相等時懲罰是最小的，所以才有了L2會讓各個系數趨於相同的特點。
可以看出，L2正則化對於特征選擇來說一種穩定的模型，不像L1正則化那樣，系數會因為細微的數據變化而波動。所以L2正則化和L1正則化提供的價值是不同的，L2正則化對於特征理解來說更加有用：表示能力強的特征對應的系數是非零。
回過頭來看看3個互相關聯的特征的例子，分別以10個不同的種子隨機初始化運行10次，來觀察L1和L2正則化的穩定性。
import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Ridge
size = 100
# We run the method 10 times with different random seeds
for i in range(10):
print("Random seed %s" % i)
np.random.seed(seed=i)
X_seed = np.random.normal(0, 1, size)
X1 = X_seed + np.random.normal(0, .1, size)
X2 = X_seed + np.random.normal(0, .1, size)
X3 = X_seed + np.random.normal(0, .1, size)
Y = X1 + X2 + X3 + np.random.normal(0, 1, size)
X = np.array([X1, X2, X3]).T
lr = LinearRegression()
lr.fit(X, Y)
print("Linear model:", pretty_print_linear(lr.coef_))
ridge = Ridge(alpha=10)
ridge.fit(X, Y)
print("Ridge model:", pretty_print_linear(ridge.coef_))
print()
>>>
Random seed 0
Linear model: 0.728 * X0 + 2.309 * X1 + -0.082 * X2
Ridge model: 0.938 * X0 + 1.059 * X1 + 0.877 * X2
Random seed 1
Linear model: 1.152 * X0 + 2.366 * X1 + -0.599 * X2
Ridge model: 0.984 * X0 + 1.068 * X1 + 0.759 * X2
Random seed 2
Linear model: 0.697 * X0 + 0.322 * X1 + 2.086 * X2
Ridge model: 0.972 * X0 + 0.943 * X1 + 1.085 * X2
Random seed 3
Linear model: 0.287 * X0 + 1.254 * X1 + 1.491 * X2
Ridge model: 0.919 * X0 + 1.005 * X1 + 1.033 * X2
Random seed 4
Linear model: 0.187 * X0 + 0.772 * X1 + 2.189 * X2
Ridge model: 0.964 * X0 + 0.982 * X1 + 1.098 * X2
Random seed 5
Linear model: -1.291 * X0 + 1.591 * X1 + 2.747 * X2
Ridge model: 0.758 * X0 + 1.011 * X1 + 1.139 * X2
Random seed 6
Linear model: 1.199 * X0 + -0.031 * X1 + 1.915 * X2
Ridge model: 1.016 * X0 + 0.89 * X1 + 1.091 * X2
Random seed 7
Linear model: 1.474 * X0 + 1.762 * X1 + -0.151 * X2
Ridge model: 1.018 * X0 + 1.039 * X1 + 0.901 * X2
Random seed 8
Linear model: 0.084 * X0 + 1.88 * X1 + 1.107 * X2
Ridge model: 0.907 * X0 + 1.071 * X1 + 1.008 * X2
Random seed 9
Linear model: 0.714 * X0 + 0.776 * X1 + 1.364 * X2
Ridge model: 0.896 * X0 + 0.903 * X1 + 0.98 * X2
可以看出，不同的數據上線性回歸得到的模型（系數）相差甚遠，但對於L2正則化模型來說，結果中的系數非常的穩定，差別較小，都比較接近於1，能夠反映出數據的內在結構。
2. 基於樹模型的特征選擇(Embedded方式)
隨機森林具有准確率高、魯棒性好、易於使用等優點，這使得它成為了目前最流行的機器學習算法之一。隨機森林提供了兩種特征選擇的方法：平均不純度減少(mean decrease impurity) 和 mean decrease accuracy。
2.1 平均不純度減少 (Mean Decrease Impurity)
隨機森林由多個決策樹構成。決策樹中的每一個節點都是關於某個特征的條件，為的是將數據集按照不同的響應變量一分為二。利用不純度可以確定節點（最優條件），對於分類問題，通常采用 基尼不純度 或者 信息增益 ，對於回歸問題，通常采用的是 方差 或者 最小二乘擬合。當訓練決策樹的時候，可以計算出每個特征減少了多少樹的不純度。對於一個決策樹森林來說，可以算出每個特征平均減少了多少不純度，並把它平均減少的不純度作為特征選擇的值。
下面的例子是sklearn中基於隨機森林的特征重要度度量方法：
from sklearn.datasets import load_boston
from sklearn.ensemble import RandomForestRegressor
import numpy as np
#Load boston housing dataset as an example
boston = load_boston()
X = boston["data"]
Y = boston["target"]
names = boston["feature_names"]
rf = RandomForestRegressor()
rf.fit(X, Y)
print("Features sorted by their score:")
print(sorted(zip(map(lambda x: "%.4f"%x, rf.feature_importances_), names), reverse=True))
Features sorted by their score:
[('0.5128', 'LSTAT'), ('0.2896', 'RM'), ('0.0792', 'DIS'), ('0.0311', 'CRIM'), ('0.0188', 'NOX'), ('0.0179', 'AGE'), ('0.0174', 'TAX'), ('0.0123', 'PTRATIO'), ('0.0086', 'B'), ('0.0074', 'RAD'), ('0.0037', 'INDUS'), ('0.0007', 'ZN'), ('0.0007', 'CHAS')]
這里特征得分實際上采用的是 Gini Importance 。使用基於不純度的方法的時候，要記住：
-
- 這種方法存在偏向 ，對具有更多類別的變量會更有利；
-
- 對於存在_關聯的多個特征_，其中任意一個都可以作為指示器（優秀的特征），並且_一旦某個特征被選擇之后，其他特征的重要度就會急劇下降_(因為不純度已經被選中的那個特征降下來了，其他的特征就很難再降低那么多不純度了，這樣一來，只有先被選中的那個特征重要度很高，其他的關聯特征重要度往往較低)。在理解數據時，這就會造成誤解，導致錯誤的認為先被選中的特征是很重要的，而其余的特征是不重要的，但實際上這些特征對響應變量的作用確實非常接近的（這跟Lasso是很像的）。
特征隨機選擇 方法稍微緩解了這個問題，但總的來說並沒有完全解決。下面的例子中，X0、X1、X2是三個互相關聯的變量，在沒有噪音的情況下，輸出變量是三者之和。
from sklearn.ensemble import RandomForestRegressor
import numpy as np
size = 10000
np.random.seed(seed=10)
X_seed = np.random.normal(0, 1, size)
X0 = X_seed + np.random.normal(0, .1, size)
X1 = X_seed + np.random.normal(0, .1, size)
X2 = X_seed + np.random.normal(0, .1, size)
X = np.array([X0, X1, X2]).T
Y = X0 + X1 + X2
rf = RandomForestRegressor(n_estimators=20, max_features=2)
rf.fit(X, Y)
print('Scores for X0, X1, X2:', ['%.3f'%x for x in rf.feature_importances_])
>>>
Scores for X0, X1, X2 ['0.272', '0.548', '0.179']
當計算特征重要性時，可以看到X1的重要度比X2的重要度要高出3倍，但實際上他們真正的重要度是一樣的。盡管數據量已經很大且沒有噪音，且用了20棵樹來做隨機選擇，但這個問題還是會存在。
需要注意的一點是，關聯特征的打分存在不穩定的現象，這不僅僅是隨機森林特有的，大多數基於模型的特征選擇方法都存在這個問題。
2.2 平均精確率減少 (Mean Decrease Accuracy)
另一種常用的特征選擇方法就是直接度量每個特征對模型精確率的影響。主要思路是打亂每個特征的特征值順序，並且度量順序變動對模型的精確率的影響。很明顯，對於不重要的變量來說，打亂順序對模型的精確率影響不會太大，但是對於重要的變量來說，打亂順序就會降低模型的精確率。
這個方法sklearn中沒有直接提供，但是很容易實現，下面繼續在波士頓房價數據集上進行實現。
from sklearn.cross_validation import ShuffleSplit
from sklearn.metrics import r2_score
from sklearn.datasets import load_boston
from collections import defaultdict
from sklearn.ensemble import RandomForestRegressor
import numpy as np
boston = load_boston()
X = boston["data"]
Y = boston["target"]
names = boston["feature_names"]
rf = RandomForestRegressor()
scores = defaultdict(list)
# crossvalidate the scores on a number of different random splits of the data
for train_idx, test_idx in ShuffleSplit(len(X), 100, .3):
X_train, X_test = X[train_idx], X[test_idx]
Y_train, Y_test = Y[train_idx], Y[test_idx]
r = rf.fit(X_train, Y_train)
acc = r2_score(Y_test, rf.predict(X_test))
for i in range(X.shape[1]):
X_t = X_test.copy()
np.random.shuffle(X_t[:, i])
shuff_acc = r2_score(Y_test, rf.predict(X_t))
scores[names[i]].append((acc - shuff_acc) / acc)
print("Features sorted by their score:")
print(sorted( [(float('%.4f'%np.mean(score)), feat) for
feat, score in scores.items()], reverse=True) )
Features sorted by their score:
[(0.7508, 'LSTAT'), (0.5691, 'RM'), (0.0947, 'DIS'), (0.0396, 'CRIM'), (0.0371, 'NOX'), (0.0223, 'PTRATIO'), (0.0173, 'TAX'), (0.0132, 'AGE'), (0.0071, 'B'), (0.0053, 'INDUS'), (0.0036, 'RAD'), (0.0005, 'CHAS'), (0.0003, 'ZN')]
在這個例子當中，LSTAT和RM這兩個特征對模型的性能有着很大的影響，打亂這兩個特征的特征值使得模型的性能下降了75%和57%。注意，盡管這些我們是在所有特征上進行了訓練得到了模型，然后才得到了每個特征的重要性測試，這並不意味着我們扔掉某個或者某些重要特征后模型的性能就一定會下降很多，因為即便某個特征刪掉之后，其關聯特征一樣可以發揮作用，讓模型性能基本上不變。
3. 頂層特征選擇算法(Wrapper方式)
之所以叫做頂層，是因為他們都是建立在基於模型的特征選擇方法基礎之上的，例如回歸和SVM，在不同的子集上建立模型，然后匯總最終確定特征得分。
3.1 穩定性選擇 (Stability Selection)
在sklearn的官方文檔中，該方法叫做隨機稀疏模型 (Randomized sparse models)：基於L1的稀疏模型的局限在於，當面對一組互相關的特征時，它們只會選擇其中一項特征。為了減輕該問題的影響可以使用隨機化技術，通過_多次重新估計稀疏模型來擾亂設計矩陣_，或通過_多次下采樣數據來統計一個給定的回歸量被選中的次數_。
穩定性選擇是一種基於 二次抽樣 和 選擇算法 相結合較新的方法，選擇算法可以是回歸、SVM或其他類似的方法。它的主要思想是在不同的數據子集和特征子集上運行特征選擇算法，不斷的重復，最終匯總特征選擇結果，比如可以統計某個特征被認為是重要特征的頻率（被選為重要特征的次數除以它所在的子集被測試的次數）。理想情況下，重要特征的得分會接近100%。稍微弱一點的特征得分會是非0的數，而最無用的特征得分將會接近於0。
sklearn在 隨機lasso(RandomizedLasso) 和 隨機邏輯回歸(RandomizedLogisticRegression) 中有對穩定性選擇的實現。
from sklearn.linear_model import RandomizedLasso
from sklearn.datasets import load_boston
boston = load_boston()
#using the Boston housing data.
#Data gets scaled automatically by sklearn's implementation
X = boston["data"]
Y = boston["target"]
names = boston["feature_names"]
rlasso = RandomizedLasso(alpha=0.025)
rlasso.fit(X, Y)
print("Features sorted by their score:")
print(sorted(zip(map(lambda x: format(x, '.4f'), rlasso.scores_), names), reverse=True))
Features sorted by their score:
[('1.0000', 'RM'), ('1.0000', 'PTRATIO'), ('1.0000', 'LSTAT'), ('0.6450', 'CHAS'), ('0.6100', 'B'), ('0.3950', 'CRIM'), ('0.3800', 'TAX'), ('0.2250', 'DIS'), ('0.2000', 'NOX'), ('0.1150', 'INDUS'), ('0.0750', 'ZN'), ('0.0200', 'RAD'), ('0.0100', 'AGE')]
在上邊這個例子當中，最高的3個特征得分是1.0，這表示他們總會被選作有用的特征（當然，得分會收到正則化參數alpha的影響，但是sklearn的隨機lasso能夠自動選擇最優的alpha）。接下來的幾個特征得分就開始下降，但是下降的不是特別急劇，這跟純lasso的方法和隨機森林的結果不一樣。能夠看出穩定性選擇對於克服過擬合和對數據理解來說都是有幫助的：總的來說，好的特征不會因為有相似的特征、關聯特征而得分為0，這跟Lasso是不同的。對於特征選擇任務，在許多數據集和環境下，穩定性選擇往往是性能最好的方法之一。
3.2 遞歸特征消除 (Recursive Feature Elimination (RFE))
遞歸特征消除的主要思想是反復的構建模型（如SVM或者回歸模型）然后選出最好的（或者最差的）的特征（可以根據系數來選），把選出來的特征放到一遍，然后在剩余的特征上重復這個過程，直到所有特征都遍歷了。這個過程中特征被消除的次序就是特征的排序。因此，這是一種尋找最優特征子集的貪心算法。
RFE的穩定性很大程度上取決於在迭代的時候底層用哪種模型。例如，假如RFE采用的普通的回歸，沒有經過正則化的回歸是不穩定的，那么RFE就是不穩定的；假如采用的是Ridge，而用Ridge正則化的回歸是穩定的，那么RFE就是穩定的。
Sklearn提供了 RFE 包，可以用於特征消除，還提供了 RFECV ，可以通過交叉驗證來對的特征進行排序。
from sklearn.feature_selection import RFE
from sklearn.linear_model import LinearRegression
from sklearn.datasets import load_boston
boston = load_boston()
X = boston["data"]
Y = boston["target"]
names = boston["feature_names"]
#use linear regression as the model
lr = LinearRegression()
#rank all features, i.e continue the elimination until the last one
rfe = RFE(lr, n_features_to_select=1)
rfe.fit(X,Y)
print("Features sorted by their rank:")
print(sorted(zip(rfe.ranking_, names)))
Features sorted by their rank:
[(1, 'NOX'), (2, 'RM'), (3, 'CHAS'), (4, 'PTRATIO'), (5, 'DIS'), (6, 'LSTAT'), (7, 'RAD'), (8, 'CRIM'), (9, 'INDUS'), (10, 'ZN'), (11, 'TAX'), (12, 'B'), (13, 'AGE')]
4. 一個完整的例子
下面將本文所有提到的方法進行實驗對比，數據集采用Friedman #1 回歸數據（ 這篇論文 中的數據）。數據是用這個公式產生的：
\(x_1\)到\(x_5\)是由 單變量分布 生成的，\(e\)是 標准正態離差 \(N(0,1)\)。另外，原始的數據集中含有5個噪音變量 \(x_6,…,x_{10}\)，跟響應變量\(y\)是獨立的，並增加了4個額外的變量\(x_{11},…x_{14}\)，分別是\(x_1,…,x_4\)的關聯變量，通過\(f(x) = x + N(0,0.01)\)生成，這將產生大於0.999的關聯系數。這樣生成的數據能夠體現出不同的特征排序方法應對關聯特征時的表現。
接下來在上述數據上運行所有的特征選擇方法，並且將每種方法給出的得分進行歸一化，讓取值都落在0-1之間。對於RFE來說，由於它給出的是順序而不是得分，我們將最好的5個的得分定為1，其他的特征的得分均勻的分布在0-1之間。
__author__ = 'steven'
from sklearn.linear_model import (LinearRegression, Ridge, Lasso, RandomizedLasso)
from sklearn.feature_selection import RFE, f_regression
from sklearn.preprocessing import MinMaxScaler
from sklearn.ensemble import RandomForestRegressor
import numpy as np
from minepy import MINE
np.random.seed(0)
size = 750
X = np.random.uniform(0, 1, (size, 14)) # 產生14個特征變量,序號為0-13;一共750個樣本
# 跟y有直接關系的是變量x1-x5(序號是0-4), x6-x10為噪聲變量,並且獨立(與其他變量不相關)
Y = (10 * np.sin(np.pi * X[:, 0] * X[:, 1]) + 20 * (X[:, 2] - .5) ** 2 +
10 * X[:, 3] + 5 * X[:, 4] + np.random.normal(0, 1))
# 變量x11,x12,x13,x14(序號10-13),這些數據由x1,x2,x3,x4加上隨機噪聲產生
X[:, 10:] = X[:, :4] + np.random.normal(0, .025, (size, 4))
feature_names = ["x%s" % i for i in range(1, 15)]
ranks = {}
# 將每個特征的得分縮放后,以字典形式存儲
def rank_to_dict(score, feature_names, order=1):
nd_array = order * np.array([score]).T # [score]為二維行向量,但是fit_transform是以列來進行規范化的,所以需要轉置
ranks = MinMaxScaler().fit_transform(nd_array).T[0] # T[0]返回ndarray的第一列,為1維數組;.T[0]等價於[:, 0]
ranks = list(map(lambda x: round(x, 2), ranks)) # 注意在python3中，map函數得到的是對象，需要用list()轉化才能得到list中的數據
return dict(zip(feature_names, ranks))
#### 單變量特征選擇
# 線性相關程度: 計算每個特征xi和應變量Y的相關程度;這里的f_regression通過F檢驗用於評估兩個隨機變量的線性相關性
f, pval = f_regression(X, Y, center=True) # 注意X一定為二維ndarray(n_samples, n_features)
ranks["Correlation"] = rank_to_dict(f, feature_names)
# # 皮爾森相關系數(Pearson Correlation Coefficient): Scipy的pearsonr方法能夠同時計算相關系數和p-value
# from scipy.stats import pearsonr
# corrs = []
# for i in range(X.shape[1]):
# corr, pval = pearsonr(X[:, i], Y) # pearsonr的系數是兩個一維ndarray數組,也可以是list
# corrs.append(corr)
# ranks["Correlation"] = rank_to_dict(corrs, feature_names)
# 最大信息系數(Maximal Information Coefficient): 計算每個特征xi和應變量Y的最大信息系數
mine = MINE()
mic_scores = []
for i in range(X.shape[1]): # shape[0]為樣本數,shape[1]為特征數
mine.compute_score(X[:, i], Y)
m = mine.mic()
mic_scores.append(m)
ranks["MIC"] = rank_to_dict(mic_scores, feature_names)
#### 線性回歸和正則化
# 回歸系數: 根據線性回歸的系數判斷特征的重要性
lr = LinearRegression(normalize=True)
lr.fit(X, Y)
ranks["LinearRegression"] = rank_to_dict(np.abs(lr.coef_), feature_names)
# l1正則: Lasso的參數
lasso = Lasso(alpha=.05)
lasso.fit(X, Y)
ranks["Lasso"] = rank_to_dict(np.abs(lasso.coef_), feature_names)
# l2正則: 嶺回歸的參數
ridge = Ridge(alpha=7)
ridge.fit(X, Y)
ranks["Ridge"] = rank_to_dict(np.abs(ridge.coef_), feature_names)
#### 隨機森林特征選擇
# 平均不純度減少(Mean Decrease Impurity): 隨機森林建樹的過程中 根據不純度選擇特征的過程
rf = RandomForestRegressor()
rf.fit(X, Y)
ranks["RandomForestRegressor"] = rank_to_dict(rf.feature_importances_, feature_names)
#### 頂層特征選擇基於基礎模型的特征選擇
# 穩定特征選擇(Stability Selection): 隨機lasso算法中實現穩定特征選擇
rlasso = RandomizedLasso(alpha=0.04)
rlasso.fit(X, Y)
ranks["Stability"] = rank_to_dict(np.abs(rlasso.scores_), feature_names)
# 遞歸特征消除(Recursive Feature Elimination): 普通線性回歸(lr)實現遞歸特征消除
# stop the search when 5 features are left (they will get equal scores)
rfe = RFE(lr, n_features_to_select=5)
rfe.fit(X, Y)
ranks["RFE_lr"] = rank_to_dict(list(map(float, rfe.ranking_)), feature_names, order=-1)
#### 根據前面每個特征選擇的方式的得到每個特征xi的平均得分
xi_mean = {}
for x_i in feature_names:
xi_mean[x_i] = round(np.mean([ranks[method][x_i] for method in ranks.keys()]), 2)
ranks["mean_score"] = xi_mean
if __name__ == '__main__':
method_order = ['Correlation', 'MIC', 'LinearRegression', 'Lasso', 'Ridge',
'RandomForestRegressor', 'Stability', 'RFE_lr', 'mean_score']
print("\t%s" % "\t".join(method_order))
for x_i in feature_names:
score_info = '\t'.join([str(ranks[method][x_i]) for method in method_order])
print("%s\t%s" % (x_i, score_info))
| Correlation | MIC | LinearRegression | Lasso | Ridge | RandomForestRegressor | Stability | RFE_lr | mean_score |
|---|---|---|---|---|---|---|---|---|
| x1 | 0.3 | 0.39 | 1 | 0.79 | 0.77 | 0.37 | 0.65 | 1 |
| x2 | 0.44 | 0.61 | 0.56 | 0.83 | 0.75 | 0.54 | 0.67 | 1 |
| x3 | 0 | 0.34 | 0.5 | 0 | 0.05 | 0.1 | 0 | 1 |
| x4 | 1 | 1 | 0.57 | 1 | 1 | 0.64 | 1 | 1 |
| x5 | 0.1 | 0.2 | 0.27 | 0.51 | 0.88 | 0.25 | 0.63 | 0.78 |
| x6 | 0 | 0 | 0.02 | 0 | 0.05 | 0 | 0 | 0.44 |
| x7 | 0.01 | 0.07 | 0 | 0 | 0.01 | 0.01 | 0 | 0 |
| x8 | 0.02 | 0.05 | 0.03 | 0 | 0.09 | 0.01 | 0 | 0.56 |
| x9 | 0.01 | 0.09 | 0 | 0 | 0 | 0 | 0 | 0.11 |
| x10 | 0 | 0.04 | 0.01 | 0 | 0.01 | 0 | 0 | 0.33 |
| x11 | 0.29 | 0.43 | 0.6 | 0 | 0.59 | 0.61 | 0.32 | 1 |
| x12 | 0.44 | 0.71 | 0.14 | 0 | 0.68 | 0.53 | 0.44 | 0.67 |
| x13 | 0 | 0.23 | 0.48 | 0 | 0.02 | 0.09 | 0 | 0.89 |
| x14 | 0.99 | 1 | 0 | 0.16 | 0.95 | 1 | 0.51 | 0.22 |
從以上結果中可以找到一些有趣的發現：
MIC對特征一視同仁，這一點上和關聯系數有點像，另外，它能夠找出\(x_3\)和響應變量之間的非線性關系。
從LinearRegression的結果中，可發現：特征之間存在線性相關關系，每個特征都是獨立評價的，因此\(x_1,…,x_4\)的得分和\(x_{11},…x_{14}\)的得分非常接近，而噪音特征\(x_6,…,x_{10}\)正如預期的那樣和響應變量之間幾乎沒有關系。由於變量\(x_3\)是二次的，因此\(x_3\)和響應變量之間看不出有關系（除了MIC之外，其他方法都找不到直接關系）。這種方法能夠衡量出特征和響應變量之間的線性關系，但若想選出優質特征來提升模型的泛化能力，這種方法就不是特別給力了，因為所有的優質特征的相關特征也會被挑出來。
Lasso能夠挑出一些優質特征，同時讓其他特征的系數趨於0。當如需要減少特征數的時候它很有用，但是對於數據理解來說不是很好用。例如在結果表中，\(x_{11}, x_{12}, x_{13}\)的得分都是0，好像他們跟輸出變量之間沒有很強的聯系，但實際上不是這樣的。
Ridge將回歸系數均勻的分攤到各個關聯變量上，從表中可以看出，\(x_{11},…x_{14}\)和\(x_1,…,x_4\)的得分非常接近。
隨機森林基於不純度的排序結果非常鮮明，在得分最高的幾個特征之后的特征，得分急劇的下降。而其他的特征選擇算法就沒有下降的這么劇烈。
穩定性選擇常常是一種既能夠有助於理解數據又能夠挑出優質特征的這種選擇，在結果表中就能很好的看出。像Lasso一樣，它能找到那些性能比較好的特征\(x_1,x_2,x_4,x_5\)，同時，與這些特征關聯度很強的變量也得到了較高的得分。
5. 總結
-
- 對於理解數據、數據的結構、特點來說，單變量特征選擇是個非常好的選擇。盡管可以用它對特征進行排序來優化模型，但由於它不能發現冗余（例如假如一個特征子集，其中的特征之間具有很強的關聯，那么從中選擇最優的特征時就很難考慮到冗余的問題）。
-
- 正則化的線性模型對於特征理解和特征選擇來說是非常強大的工具。L1正則化能夠生成稀疏的模型，對於選擇特征子集來說非常有用；相比起L1正則化，L2正則化的表現更加穩定，由於有用的特征往往對應系數非零，因此L2正則化對於數據的理解來說很合適。由於響應變量和特征之間往往是非線性關系，可以采用basis expansion(基展開)的方式將特征轉換到一個更加合適的空間當中，在此基礎上再考慮運用簡單的線性模型。
-
- 隨機森林是一種非常流行的特征選擇方法，它易於使用，一般不需要feature engineering、調參等繁瑣的步驟，並且很多工具包都提供了平均不純度下降方法。它的兩個主要問題，1是重要的特征有可能得分很低（關聯特征問題），2是這種方法對特征變量類別多的特征越有利（偏向問題）。盡管如此，這種方法仍然非常值得在你的應用中試一試。
-
- 特征選擇在很多機器學習和數據挖掘場景中都是非常有用的。在使用的時候要弄清楚自己的目標是什么，然后找到哪種方法適用於自己的任務；當選擇最優特征以提升模型性能的時候，可以采用交叉驗證的方法來驗證某種方法是否比其他方法要好；當用特征選擇的方法來理解數據的時候要留心，特征選擇模型的穩定性非常重要，穩定性差的模型很容易就會導致錯誤的結論；數據進行二次采樣然后在子集上運行特征選擇算法能夠有所幫助，如果在各個子集上的結果是一致的，那就可以說在這個數據集上得出來的結論是可信的，可以用這種特征選擇模型的結果來理解數據。
6. Tips
什么是 卡方檢驗 ？
用方差來衡量某個觀測頻率和理論頻率之間差異性的方法。
什么是 皮爾森卡方檢驗 ？
這是一種最常用的卡方檢驗方法，它有兩個用途：1是計算某個變量對某種分布的擬合程度，2是根據兩個觀測變量的 Contingency table 來計算這兩個變量是否是獨立的。主要有三個步驟：第一步用方差和的方式來計算觀測頻率和理論頻率之間卡方值；第二步算出卡方檢驗的自由度（行數-1乘以列數-1）；第三步比較卡方值和對應自由度的卡方分布，判斷顯著性。
什么是 p-value ？
簡單地說，p-value就是為了驗證假設和實際之間一致性的統計學意義的值，即假設檢驗。有些地方叫右尾概率，根據卡方值和自由度可以算出一個固定的p-value，
什么是 統計能力(statistical power) ?
什么是 零假設(null hypothesis) ?
在相關性檢驗中，一般會取“兩者之間無關聯”作為零假設，而在獨立性檢驗中，一般會取“兩者之間是獨立”作為零假設。與零假設相對的是備擇假設（對立假設），即希望證明是正確的另一種可能。
什么是 多重共線性 ？
7. References
http://blog.datadive.net/selecting-good-features-part-i-univariate-selection/
http://blog.datadive.net/selecting-good-features-part-ii-linear-models-and-regularization/
http://scikit-learn.org/stable/modules/feature_selection.html#univariate-feature-selection
http://www.quora.com/What-are-some-feature-selection-methods
http://www.quora.com/What-are-some-feature-selection-algorithms
http://www.quora.com/What-are-some-feature-selection-methods-for-SVMs
http://www.quora.com/What-is-the-difference-between-principal-component-analysis-PCA-and-feature-selection-in-machine-learning-Is-PCA-a-means-of-feature-selection
文章出處： http://chaoslog.com/te-zheng-xuan-ze.html