基于模型的特征選擇詳解 (Embedded & Wrapper)

Z2ty6osc12zs6c 2018-04-17

展開全文

[toc]

　　單變量特征選擇方法獨立的衡量每個特征與響應變量之間的關(guān)系，另一種主流的特征選擇方法是基于機器學習模型的方法,。有些機器學習方法本身就具有對特征進行打分的機制,，或者很容易將其運用到特征選擇任務(wù)中，例如回歸模型,，SVM,，決策樹，隨機森林等等,。

1. 線性模型和正則化(Embedded方式)

　　下面將介紹如何用==回歸==模型的系數(shù)來選擇特征,。越是重要的特征在模型中對應的系數(shù)就會越大，而跟輸出變量越是無關(guān)的特征對應的系數(shù)就會越接近于0,。在噪音不多的數(shù)據(jù)上,，或者是數(shù)據(jù)量遠遠大于特征數(shù)的數(shù)據(jù)上，如果特征之間相對來說是比較獨立的,，那么即便是運用最簡單的線性回歸模型也一樣能取得非常好的效果,。

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)  # 有了這段代碼，下次再生成隨機數(shù)的時候,，與上次一樣的結(jié)果
size = 5000  # 表示抽取多少個樣本
# 隨機生成3個特征的樣本,，每個維度的特征都是服從期望為0，標準差為1的正態(tài)分布
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

在這個例子當中,，盡管數(shù)據(jù)中存在一些噪音，但這種特征選擇模型仍然能夠很好的體現(xiàn)出數(shù)據(jù)的底層結(jié)構(gòu),。當然這也是因為例子中的這個問題非常適合用線性模型來解：特征和響應變量之間全都是線性關(guān)系,，并且特征之間均是獨立的。

　　然而,，在很多實際的數(shù)據(jù)當中,，往往存在多個互相關(guān)聯(lián)的特征，這時候模型就會變得不穩(wěn)定,，對噪聲很敏感,，數(shù)據(jù)中細微的變化就可能導致模型的巨大變化（模型的變化本質(zhì)上是系數(shù),，或者叫參數(shù)），這會讓模型的預測變得困難,，這種現(xiàn)象也稱為多重共線性,。例如，假設(shè)我們有個數(shù)據(jù)集,，它的真實模型應該是 $Y = X_{1} + X_{2}$ ,，當我們觀察的時候，發(fā)現(xiàn) $Y^{'} = X_{1} + X_{2} + e$ ,， $e$ 是噪音,。如果 $X_{1}$ 和 $X_{2}$ 之間存在線性關(guān)系，例如 $X_{1}$ 約等于 $X_{2}$ ,，這個時候由于噪音 $e$ 的存在,，我們學到的模型可能就不是 $Y = X_{1} + X_{2}$ 了，有可能是 $Y = 2 X_{1}$ ,，或者 $Y = - X_{1} + 3 X_{2}$ ,。【在 $X_{1}$ 約等于 $X_{2}$ 的線性關(guān)系下,，學到的這三種模型,，理論上來說都是正確的，注意他們有個共同的特點是,，系數(shù)之和為2】

size = 100
# 另外一個種子為5的隨機采樣，若要執(zhí)行相同的結(jié)果,，以種子號來區(qū)分隨機采樣的結(jié)果
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

　　這個例子和上個例子中,，三個特征都來源于標準正態(tài)分布的隨機采樣，這里系數(shù)之和接近3,，基本上和上上個例子的結(jié)果一致,，應該說學到的模型對于預測來說還是不錯的。但是,，如果從系數(shù)的字面意思上去解釋特征的重要性的話,，X2對于輸出變量來說具有很強的正面影響，而X0具有負面影響,，而實際上所有特征與輸出變量之間的影響是均等的,。
　　同樣的方法和套路可以用到類似的線性模型上，比如邏輯回歸,。

正則化的概念：

　　正則化就是把額外的約束或者懲罰項加到已有模型（損失函數(shù)）上,，以防止過擬合并提高泛化能力。損失函數(shù)由原來的 $L (X, Y)$ 變?yōu)?L(X,Y)+\alpha \left| w \right| $，$ w$是模型系數(shù)組成的向量（有些地方也叫參數(shù)parameter,，coefficients）， $‖ \cdot ‖$ 一般是L1或者L2范數(shù),， $α$ 是一個可調(diào)的參數(shù),，控制著正則化的強度。當用在線性模型上時,，L1正則化和L2正則化也稱為Lasso(least absolute shrinkage and selection operator)和Ridge,。

1.1 L1正則化(Lasso)

　　L1正則化將系數(shù) $w$ 的l1范數(shù)作為懲罰項加到損失函數(shù)上，由于正則項非零,，這就迫使那些弱的特征所對應的系數(shù)變成0,。因此L1正則化往往會使學到的模型很稀疏（系數(shù)w經(jīng)常為0），這個特性使得L1正則化成為一種很好的特征選擇方法,。

　　Scikit-learn為線性回歸模型提供了Lasso,，為分類模型提供了L1邏輯回歸。下面的例子在波士頓房價數(shù)據(jù)上運行了Lasso,，其中參數(shù)alpha是通過grid search進行優(yōu)化的,。

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

可以看到，很多特征的系數(shù)都是0,。如果繼續(xù)增加 $α$ 的值,，得到的模型就會越來越稀疏，即越來越多的特征系數(shù)會變成0,。

　　然而,，L1正則化像非正則化線性模型一樣也是不穩(wěn)定的，如果特征集合中具有相關(guān)聯(lián)的特征,，當數(shù)據(jù)發(fā)生細微變化時也有可能導致很大的模型差異,。

1.2 L2正則化(Ridge Regression)

　　L2正則化同樣將系數(shù)向量的L2范數(shù)添加到了損失函數(shù)中。由于L2懲罰項中系數(shù)是二次方的,，這使得L2和L1有著諸多差異,，最明顯的一點就是，L2正則化會讓系數(shù)的取值變得平均,。對于關(guān)聯(lián)特征,，這意味著他們能夠獲得更相近的對應系數(shù)。還是以 $Y = X_{1} + X_{2}$ 為例,，假設(shè) $X_{1}$ 和 $X_{2}$ 具有很強的關(guān)聯(lián),，如果用L1正則化，不論學到的模型是 $Y = X_{1} + X_{2}$ 還是 $Y = 2 X_{1}$ ,，懲罰都是一樣的,，都是2 $α$ 。但是對于L2來說,，第一個模型的懲罰項是 $2 α$ ,，但第二個模型的是 $4 α$ ,。可以看出，系數(shù)(待求參數(shù))之和為常數(shù)時,，各系數(shù)相等時懲罰是最小的,，所以才有了L2會讓各個系數(shù)趨于相同的特點。

　　可以看出,，L2正則化對于特征選擇來說一種穩(wěn)定的模型,，不像L1正則化那樣，系數(shù)會因為細微的數(shù)據(jù)變化而波動,。所以L2正則化和L1正則化提供的價值是不同的,，L2正則化對于特征理解來說更加有用：表示能力強的特征對應的系數(shù)是非零。

　　回過頭來看看3個互相關(guān)聯(lián)的特征的例子,，分別以10個不同的種子隨機初始化運行10次,，來觀察L1和L2正則化的穩(wěn)定性。

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

可以看出,，不同的數(shù)據(jù)上線性回歸得到的模型（系數(shù)）相差甚遠,，但對于L2正則化模型來說，結(jié)果中的系數(shù)非常的穩(wěn)定,，差別較小,，都比較接近于1，能夠反映出數(shù)據(jù)的內(nèi)在結(jié)構(gòu),。

2. 基于樹模型的特征選擇(Embedded方式)

　　隨機森林具有準確率高,、魯棒性好、易于使用等優(yōu)點,，這使得它成為了目前最流行的機器學習算法之一,。隨機森林提供了兩種特征選擇的方法：平均不純度減少(mean decrease impurity) 和 mean decrease accuracy。

2.1 平均不純度減少 (Mean Decrease Impurity)

　　隨機森林由多個決策樹構(gòu)成,。決策樹中的每一個節(jié)點都是關(guān)于某個特征的條件，為的是將數(shù)據(jù)集按照不同的響應變量一分為二,。利用不純度可以確定節(jié)點（最優(yōu)條件）,，對于分類問題，通常采用 ==基尼不純度== 或者 ==信息增益== ,，對于回歸問題,，通常采用的是 ==方差== 或者 ==最小二乘擬合==。當訓練決策樹的時候,，可以計算出每個特征減少了多少樹的不純度,。對于一個決策樹森林來說，可以算出每個特征平均減少了多少不純度,，并把它平均減少的不純度作為特征選擇的值,。

下面的例子是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 ,。使用基于不純度的方法的時候，要記?。?/p>

1. 這種方法存在偏向 ,，對具有更多類別的變量會更有利；
1. 對于存在_關(guān)聯(lián)的多個特征_,，其中任意一個都可以作為指示器（優(yōu)秀的特征）,，并且_一旦某個特征被選擇之后，其他特征的重要度就會急劇下降_(因為不純度已經(jīng)被選中的那個特征降下來了,，其他的特征就很難再降低那么多不純度了,，這樣一來，只有先被選中的那個特征重要度很高,，其他的關(guān)聯(lián)特征重要度往往較低),。在理解數(shù)據(jù)時，這就會造成誤解,，導致錯誤的認為先被選中的特征是很重要的,，而其余的特征是不重要的，但實際上這些特征對響應變量的作用確實非常接近的（這跟Lasso是很像的）,。

　　特征隨機選擇方法稍微緩解了這個問題,，但總的來說并沒有完全解決。下面的例子中,，X0,、X1、X2是三個互相關(guān)聯(lián)的變量,，在沒有噪音的情況下,，輸出變量是三者之和。

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倍,，但實際上他們真正的重要度是一樣的。盡管數(shù)據(jù)量已經(jīng)很大且沒有噪音,，且用了20棵樹來做隨機選擇,，但這個問題還是會存在。

　　需要注意的一點是,，關(guān)聯(lián)特征的打分存在不穩(wěn)定的現(xiàn)象,，這不僅僅是隨機森林特有的，大多數(shù)基于模型的特征選擇方法都存在這個問題,。

2.2 平均精確率減少 (Mean Decrease Accuracy)

　　另一種常用的特征選擇方法就是直接度量每個特征對模型精確率的影響,。主要思路是打亂每個特征的特征值順序，并且度量順序變動對模型的精確率的影響,。很明顯,，對于不重要的變量來說,，打亂順序?qū)δＰ偷木_率影響不會太大，但是對于重要的變量來說,，打亂順序就會降低模型的精確率,。

　　這個方法sklearn中沒有直接提供，但是很容易實現(xiàn),，下面繼續(xù)在波士頓房價數(shù)據(jù)集上進行實現(xiàn),。

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%,。注意，盡管這些我們是在所有特征上進行了訓練得到了模型,，然后才得到了每個特征的重要性測試,，這并不意味著我們?nèi)拥裟硞€或者某些重要特征后模型的性能就一定會下降很多,，因為即便某個特征刪掉之后,，其關(guān)聯(lián)特征一樣可以發(fā)揮作用,，讓模型性能基本上不變,。

3. 頂層特征選擇算法(Wrapper方式)

　　之所以叫做頂層，是因為他們都是建立在基于模型的特征選擇方法基礎(chǔ)之上的,，例如回歸和SVM,，在不同的子集上建立模型,，然后匯總最終確定特征得分,。

3.1 穩(wěn)定性選擇 (Stability Selection)

　　在sklearn的官方文檔中,，該方法叫做隨機稀疏模型 (Randomized sparse models)：基于L1的稀疏模型的局限在于，當面對一組互相關(guān)的特征時,，它們只會選擇其中一項特征,。為了減輕該問題的影響可以使用隨機化技術(shù)，通過_多次重新估計稀疏模型來擾亂設(shè)計矩陣_,，或通過_多次下采樣數(shù)據(jù)來統(tǒng)計一個給定的回歸量被選中的次數(shù)_,。
　　穩(wěn)定性選擇是一種基于二次抽樣和選擇算法相結(jié)合較新的方法，選擇算法可以是回歸,、SVM或其他類似的方法,。它的主要思想是在不同的數(shù)據(jù)子集和特征子集上運行特征選擇算法，不斷的重復,，最終匯總特征選擇結(jié)果，比如可以統(tǒng)計某個特征被認為是重要特征的頻率（被選為重要特征的次數(shù)除以它所在的子集被測試的次數(shù)）,。理想情況下,，重要特征的得分會接近100%。稍微弱一點的特征得分會是非0的數(shù),，而最無用的特征得分將會接近于0,。

　　sklearn在隨機lasso(RandomizedLasso) 和隨機邏輯回歸(RandomizedLogisticRegression) 中有對穩(wěn)定性選擇的實現(xiàn),。

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,，這表示他們總會被選作有用的特征（當然,，得分會收到正則化參數(shù)alpha的影響，但是sklearn的隨機lasso能夠自動選擇最優(yōu)的alpha）,。接下來的幾個特征得分就開始下降,，但是下降的不是特別急劇，這跟純lasso的方法和隨機森林的結(jié)果不一樣,。能夠看出穩(wěn)定性選擇對于克服過擬合和對數(shù)據(jù)理解來說都是有幫助的：總的來說,，好的特征不會因為有相似的特征、關(guān)聯(lián)特征而得分為0,，這跟Lasso是不同的,。對于特征選擇任務(wù)，在許多數(shù)據(jù)集和環(huán)境下,，穩(wěn)定性選擇往往是性能最好的方法之一,。

3.2 遞歸特征消除 (Recursive Feature Elimination (RFE))

　　遞歸特征消除的主要思想是反復的構(gòu)建模型（如SVM或者回歸模型）然后選出最好的（或者最差的）的特征（可以根據(jù)系數(shù)來選），把選出來的特征放到一遍,，然后在剩余的特征上重復這個過程,，直到所有特征都遍歷了。這個過程中特征被消除的次序就是特征的排序,。因此,，這是一種尋找最優(yōu)特征子集的貪心算法。

　　RFE的穩(wěn)定性很大程度上取決于在迭代的時候底層用哪種模型,。例如,，假如RFE采用的普通的回歸，沒有經(jīng)過正則化的回歸是不穩(wěn)定的,，那么RFE就是不穩(wěn)定的,；假如采用的是Ridge，而用Ridge正則化的回歸是穩(wěn)定的,，那么RFE就是穩(wěn)定的,。

　　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. 一個完整的例子

　　下面將本文所有提到的方法進行實驗對比,，數(shù)據(jù)集采用Friedman #1 回歸數(shù)據(jù)（這篇論文中的數(shù)據(jù)）,。數(shù)據(jù)是用這個公式產(chǎn)生的：

y = 10 \sin (π x_{1} x_{2}) + 20 (x_{3} - 0.5)^{2} + 10 x_{4} + 5 x_{5} + e

$x_{1}$ 到 $x_{5}$ 是由單變量分布生成的， $e$ 是標準正態(tài)離差 $N (0, 1)$ 。另外,，原始的數(shù)據(jù)集中含有5個噪音變量 $x_{6}, \dots, x_{10}$ ,，跟響應變量 $y$ 是獨立的，并增加了4個額外的變量 $x_{11}, \dots x_{14}$ ,，分別是 $x_{1}, \dots, x_{4}$ 的關(guān)聯(lián)變量,，通過 $f (x) = x + N (0, 0.01)$ 生成，這將產(chǎn)生大于0.999的關(guān)聯(lián)系數(shù),。這樣生成的數(shù)據(jù)能夠體現(xiàn)出不同的特征排序方法應對關(guān)聯(lián)特征時的表現(xiàn),。

　　接下來在上述數(shù)據(jù)上運行所有的特征選擇方法，并且將每種方法給出的得分進行歸一化,，讓取值都落在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))   # 產(chǎn)生14個特征變量,序號為0-13;一共750個樣本
# 跟y有直接關(guān)系的是變量x1-x5(序號是0-4), x6-x10為噪聲變量,并且獨立(與其他變量不相關(guān))
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),這些數(shù)據(jù)由x1,x2,x3,x4加上隨機噪聲產(chǎn)生
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是以列來進行規(guī)范化的,所以需要轉(zhuǎn)置
    ranks = MinMaxScaler().fit_transform(nd_array).T[0]   # T[0]返回ndarray的第一列,為1維數(shù)組;.T[0]等價于[:, 0]
    ranks = list(map(lambda x: round(x, 2), ranks))   # 注意在python3中,，map函數(shù)得到的是對象,，需要用list()轉(zhuǎn)化才能得到list中的數(shù)據(jù)
    return dict(zip(feature_names, ranks))

#### 單變量特征選擇
# 線性相關(guān)程度: 計算每個特征xi和應變量Y的相關(guān)程度;這里的f_regression通過F檢驗用于評估兩個隨機變量的線性相關(guān)性
f, pval = f_regression(X, Y, center=True)  # 注意X一定為二維ndarray(n_samples, n_features)
ranks["Correlation"] = rank_to_dict(f, feature_names)

# # 皮爾森相關(guān)系數(shù)(Pearson Correlation Coefficient): Scipy的pearsonr方法能夠同時計算相關(guān)系數(shù)和p-value
# from scipy.stats import pearsonr
# corrs = []
# for i in range(X.shape[1]):
#     corr, pval = pearsonr(X[:, i], Y)   # pearsonr的系數(shù)是兩個一維ndarray數(shù)組,也可以是list
#     corrs.append(corr)
# ranks["Correlation"] = rank_to_dict(corrs, feature_names)

# 最大信息系數(shù)(Maximal Information Coefficient): 計算每個特征xi和應變量Y的最大信息系數(shù)
mine = MINE()
mic_scores = []
for i in range(X.shape[1]):  # shape[0]為樣本數(shù),shape[1]為特征數(shù)
    mine.compute_score(X[:, i], Y)
    m = mine.mic()
    mic_scores.append(m)
ranks["MIC"] = rank_to_dict(mic_scores, feature_names)

#### 線性回歸和正則化
# 回歸系數(shù): 根據(jù)線性回歸的系數(shù)判斷特征的重要性
lr = LinearRegression(normalize=True)
lr.fit(X, Y)
ranks["LinearRegression"] = rank_to_dict(np.abs(lr.coef_), feature_names)

# l1正則: Lasso的參數(shù)
lasso = Lasso(alpha=.05)
lasso.fit(X, Y)
ranks["Lasso"] = rank_to_dict(np.abs(lasso.coef_), feature_names)

# l2正則: 嶺回歸的參數(shù)
ridge = Ridge(alpha=7)
ridge.fit(X, Y)
ranks["Ridge"] = rank_to_dict(np.abs(ridge.coef_), feature_names)

#### 隨機森林特征選擇
# 平均不純度減少(Mean Decrease Impurity): 隨機森林建樹的過程中 根據(jù)不純度選擇特征的過程
rf = RandomForestRegressor()
rf.fit(X, Y)
ranks["RandomForestRegressor"] = rank_to_dict(rf.feature_importances_, feature_names)

#### 頂層特征選擇基于基礎(chǔ)模型的特征選擇
# 穩(wěn)定特征選擇(Stability Selection): 隨機lasso算法中實現(xiàn)穩(wěn)定特征選擇
rlasso = RandomizedLasso(alpha=0.04)
rlasso.fit(X, Y)
ranks["Stability"] = rank_to_dict(np.abs(rlasso.scores_), feature_names)

# 遞歸特征消除(Recursive Feature Elimination): 普通線性回歸(lr)實現(xiàn)遞歸特征消除
# 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)

#### 根據(jù)前面每個特征選擇的方式的得到每個特征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 0.66

x2 0.44 0.61 0.56 0.83 0.75 0.54 0.67 1 0.68

x3 0 0.34 0.5 0 0.05 0.1 0 1 0.25

x4 1 1 0.57 1 1 0.64 1 1 0.9

x5 0.1 0.2 0.27 0.51 0.88 0.25 0.63 0.78 0.45

x6 0 0 0.02 0 0.05 0 0 0.44 0.06

x7 0.01 0.07 0 0 0.01 0.01 0 0 0.01

x8 0.02 0.05 0.03 0 0.09 0.01 0 0.56 0.1

x9 0.01 0.09 0 0 0 0 0 0.11 0.03

x10 0 0.04 0.01 0 0.01 0 0 0.33 0.05

x11 0.29 0.43 0.6 0 0.59 0.61 0.32 1 0.48

x12 0.44 0.71 0.14 0 0.68 0.53 0.44 0.67 0.45

x13 0 0.23 0.48 0 0.02 0.09 0 0.89 0.21

x14 0.99 1 0 0.16 0.95 1 0.51 0.22 0.6

從以上結(jié)果中可以找到一些有趣的發(fā)現(xiàn)：

　　MIC對特征一視同仁，這一點上和關(guān)聯(lián)系數(shù)有點像,，另外,，它能夠找出 $x_{3}$ 和響應變量之間的非線性關(guān)系。

　　從LinearRegression的結(jié)果中,，可發(fā)現(xiàn)：特征之間存在線性相關(guān)關(guān)系,，每個特征都是獨立評價的，因此 $x_{1}, \dots, x_{4}$ 的得分和 $x_{11}, \dots x_{14}$ 的得分非常接近,，而噪音特征 $x_{6}, \dots, x_{10}$ 正如預期的那樣和響應變量之間幾乎沒有關(guān)系,。由于變量 $x_{3}$ 是二次的，因此 $x_{3}$ 和響應變量之間看不出有關(guān)系（除了MIC之外,，其他方法都找不到直接關(guān)系）,。這種方法能夠衡量出特征和響應變量之間的線性關(guān)系，但若想選出優(yōu)質(zhì)特征來提升模型的泛化能力,，這種方法就不是特別給力了,，因為所有的優(yōu)質(zhì)特征的相關(guān)特征也會被挑出來。

　　Lasso能夠挑出一些優(yōu)質(zhì)特征,，同時讓其他特征的系數(shù)趨于0,。當如需要減少特征數(shù)的時候它很有用，但是對于數(shù)據(jù)理解來說不是很好用,。例如在結(jié)果表中,， $x_{11}, x_{12}, x_{13}$ 的得分都是0，好像他們跟輸出變量之間沒有很強的聯(lián)系，但實際上不是這樣的,。

　　Ridge將回歸系數(shù)均勻的分攤到各個關(guān)聯(lián)變量上，從表中可以看出,， $x_{11}, \dots x_{14}$ 和 $x_{1}, \dots, x_{4}$ 的得分非常接近,。

　　隨機森林基于不純度的排序結(jié)果非常鮮明，在得分最高的幾個特征之后的特征,，得分急劇的下降,。而其他的特征選擇算法就沒有下降的這么劇烈。

　　穩(wěn)定性選擇常常是一種既能夠有助于理解數(shù)據(jù)又能夠挑出優(yōu)質(zhì)特征的這種選擇,，在結(jié)果表中就能很好的看出,。像Lasso一樣，它能找到那些性能比較好的特征 $x_{1}, x_{2}, x_{4}, x_{5}$ ,，同時,，與這些特征關(guān)聯(lián)度很強的變量也得到了較高的得分。

5. 總結(jié)

1. 對于理解數(shù)據(jù),、數(shù)據(jù)的結(jié)構(gòu),、特點來說，單變量特征選擇是個非常好的選擇,。盡管可以用它對特征進行排序來優(yōu)化模型,，但由于它不能發(fā)現(xiàn)冗余（例如假如一個特征子集，其中的特征之間具有很強的關(guān)聯(lián),，那么從中選擇最優(yōu)的特征時就很難考慮到冗余的問題）,。
1. 正則化的線性模型對于特征理解和特征選擇來說是非常強大的工具。L1正則化能夠生成稀疏的模型,，對于選擇特征子集來說非常有用,；相比起L1正則化，L2正則化的表現(xiàn)更加穩(wěn)定,，由于有用的特征往往對應系數(shù)非零,，因此L2正則化對于數(shù)據(jù)的理解來說很合適。由于響應變量和特征之間往往是非線性關(guān)系,，可以采用basis expansion(基展開)的方式將特征轉(zhuǎn)換到一個更加合適的空間當中,，在此基礎(chǔ)上再考慮運用簡單的線性模型。
1. 隨機森林是一種非常流行的特征選擇方法,，它易于使用,，一般不需要feature engineering、調(diào)參等繁瑣的步驟,，并且很多工具包都提供了平均不純度下降方法,。它的兩個主要問題，1是重要的特征有可能得分很低（關(guān)聯(lián)特征問題），2是這種方法對特征變量類別多的特征越有利（偏向問題）,。盡管如此,，這種方法仍然非常值得在你的應用中試一試。
1. 特征選擇在很多機器學習和數(shù)據(jù)挖掘場景中都是非常有用的,。在使用的時候要弄清楚自己的目標是什么,，然后找到哪種方法適用于自己的任務(wù)；當選擇最優(yōu)特征以提升模型性能的時候,，可以采用交叉驗證的方法來驗證某種方法是否比其他方法要好,；當用特征選擇的方法來理解數(shù)據(jù)的時候要留心，特征選擇模型的穩(wěn)定性非常重要,，穩(wěn)定性差的模型很容易就會導致錯誤的結(jié)論,；數(shù)據(jù)進行二次采樣然后在子集上運行特征選擇算法能夠有所幫助，如果在各個子集上的結(jié)果是一致的,，那就可以說在這個數(shù)據(jù)集上得出來的結(jié)論是可信的,，可以用這種特征選擇模型的結(jié)果來理解數(shù)據(jù)。

6. Tips

什么是卡方檢驗 ,？
　　用方差來衡量某個觀測頻率和理論頻率之間差異性的方法,。

什么是皮爾森卡方檢驗？
　　這是一種最常用的卡方檢驗方法,，它有兩個用途：1是計算某個變量對某種分布的擬合程度,，2是根據(jù)兩個觀測變量的 Contingency table 來計算這兩個變量是否是獨立的。主要有三個步驟：第一步用方差和的方式來計算觀測頻率和理論頻率之間卡方值,；第二步算出卡方檢驗的自由度（行數(shù)-1乘以列數(shù)-1）,；第三步比較卡方值和對應自由度的卡方分布，判斷顯著性,。

什么是 p-value ,？
　　簡單地說，p-value就是為了驗證假設(shè)和實際之間一致性的統(tǒng)計學意義的值,，即假設(shè)檢驗,。有些地方叫右尾概率，根據(jù)卡方值和自由度可以算出一個固定的p-value,，

什么是統(tǒng)計能力(statistical power) ?

什么是零假設(shè)(null hypothesis) ?
　　在相關(guān)性檢驗中,，一般會取“兩者之間無關(guān)聯(lián)”作為零假設(shè)，而在獨立性檢驗中,，一般會取“兩者之間是獨立”作為零假設(shè),。與零假設(shè)相對的是備擇假設(shè)（對立假設(shè)），即希望證明是正確的另一種可能,。

什么是多重共線性 ,？