第一步：读取数据

我这里是把文件存在当前文件夹，新建的data文件夹里。如果不知道把数据放在哪里可以先看一下path是什么，然后把data 文件夹存在path下的路径中。

file = 'bennett20tm.xlsx'
path = os.path.join(os.path.dirname(os.getcwd()),
                    'data',file)
df = pd.read_excel(path)

读完数据后可以看一下df的内容如下：读取到的数据有10行，其中Al,Fe,Mn,S,TOC,U_PPM,V_PPM,MO_PPM是岩石的各种化学物质的含量，也是岩性识别的特征。type则是岩石的类型也是岩性识别的标签。

---

---

第二步：对特征数据处理标准化

# standardization
from sklearn.preprocessing import StandardScaler
data_scaler = StandardScaler()
data_scaler.fit(df[df.columns[2:]])
scaled_dataframe = data_scaler.transform(df[df.columns[2:]])

这里使用的是z-core标准化，是将数据按期属性（按列进行）减去其均值，并除以其标准差。得到的结果是，对于每个属性/每列来说所有数据都聚集在0附近，方差为1。

其中df[df.columns[2:]]是df中从Al到MO_PPM所有的数据，标准化后的数据存在数组scaled_dataframe中，数据如下：

---

---

第三步：相关系数矩阵

new_cols = [i.split('_')[0] for i in df.columns[2:]]
df_scaled = pd.DataFrame(scaled_dataframe,
                         columns=new_cols)
cor_matrix = df_scaled.cov().round(2)
cor_mask = np.triu(np.ones_like(cor_matrix,dtype=bool))
import seaborn as sb
sb.heatmap(cor_matrix,annot=True,vmax=1,vmin=-1,center=0,
           cmap='vlag',mask=cor_mask)

其中new_cols是去除了df表头里后面的-%，只保留了Al,Fe,.....作为新表df_scaled的表头。然后cor_maxtrix是求df_scaled的相关系数，并且保留两位小数。后面的代码则是把相关系数矩阵可视化。后得到的结果为：

可以看到其中Mo和Mn的相关性比较低，值越小，相关性越低。

第四步：主成分分析（PCA）

首先介绍一下pca:主成分分析（Principal Component Analysis），是一种用于探索高维数据的技术。PCA通常用于高维数据集的探索与可视化。还可以用于数据压缩，数据预处理等。PCA可以把可能具有线性相关性的高维变量合成为线性无关的低维变量，称为主成分（principal components），新的低维数据集会尽可能的保留原始数据的变量，可以将高维数据集映射到低维空间的同时，尽可能的保留更多变量。

pca = PCA(n_components = 2) 
pca.fit(scaled_dataframe)
x_pca = pca.transform(scaled_dataframe)
np.unique(df.type)

mapping = {'BPOMZ':0,
           'CPOMZ':1,
           'P-Euxinic':2}
df5 = df['type'].replace(mapping,inplace=False)
df5=df5.to_frame()

plt.scatter(x_pca[:,0],x_pca[:,1],
            c=df5['type']) # c needs numerical data
plt.xlabel('first component')
plt.ylabel('second component')

比如这里设置了n_components=2,等于是把Al,Fe,Mn,S,TOC,U_PPM,V_PPM,MO_PPM这八个特征线性组合在一起，构成两个线性无关的特征。可以观察特征降维后的散点图，如下：

这两个降维后的特征基本上可以完全把三类岩石分隔开。（重合的地方较少）。但是考虑到两个特征还是比较少，所以我又设置了n_components=4，做了主成分分析。代码如下：

pca = PCA(n_components = 4) 
pca.fit(scaled_dataframe)
x_pca = pca.transform(scaled_dataframe)

df_out = pd.DataFrame(x_pca, 
                      columns = ['pc1',
                                 'pc2',
                                 'pc3',
                                 'pc4'])
df_out['type'] = df['type']

到这里数据预处理基本完成，下面要开始正式的分类了。

第五步：划分数据集

先把特征和标签用两个矩阵分别开，然后在整体中取30%的数据为测试集，70%的数据为训练集。

features = ['pc1','pc2','pc3','pc4']
label = 'type'

# read features and labels
#X = df_out[features].to_numpy()
#y = df_out[label].to_numpy()

X = df_out[features].values
y = df_out[label].as_matrix()



# split
X_train, X_test, y_train, y_test = \
    train_test_split(X,y,test_size=0.3,
                     random_state=0)
    

第六步：开始建模

# random forest
# let pipeline do the procedures automatically
pipe = Pipeline([('sr',StandardScaler()),
                 ('clf',RandomForestClassifier(
                     random_state=0))
    ])
                 
# hyperparameter optimization
param_range = [2,3,4,5,6,8,10,12,15,18,20]
param_grid = [{'clf__n_estimators': param_range,
               'clf__criterion':['entropy','gini'],
               'clf__max_depth':[2,4,6,8,10,12,14]
    }]

gs = GridSearchCV(
    estimator=pipe,
    param_grid=param_grid,
    scoring='accuracy',
    cv = 10,
    n_jobs = -1
    )

gs.fit(X_train, y_train)

print('Best params: ', gs.best_params_)

# assess performance using 
# test data that model has nevern seen

lr = gs.best_estimator_ # get best classifier
lr.fit(X_train, y_train)

y_pred = lr.predict(X_test)

# check accuracy
cm = confusion_matrix(y_test,y_pred)
display = ConfusionMatrixDisplay(
    cm,display_labels=lr.classes_).plot()
score = lr.score(X_test,y_test)

# KNN
pipe = Pipeline([('sr',StandardScaler()),
                 ('clf',KNeighborsClassifier())
    ])
                 
# hyperparameter optimization
param_range = [2,3,4,5,6,8,10,12,15,18,20]
param_grid = [{'clf__n_neighbors': param_range,
               'clf__weights':['uniform','distance']
    }] # uniform:距离大小不决定决策，distance：距离影响决策

gs = GridSearchCV(
    estimator=pipe,
    param_grid=param_grid,
    scoring='accuracy',
    cv = 10,
    n_jobs = -1
    )

gs.fit(X_train, y_train)
print('Best params: ', gs.best_params_)

# assess performance using 
# test data that model has nevern seen

lr = gs.best_estimator_ # get best classifier
lr.fit(X_train, y_train)

y_pred = lr.predict(X_test)

# check accuracy
cm = confusion_matrix(y_test,y_pred)
display = ConfusionMatrixDisplay(
    cm,display_labels=lr.classes_).plot()
score = lr.score(X_test,y_test)

这里使用了随机森林和knn两种比较简单的分类算法，关于参数调节是直接用了GridSearchCV这个库可以直接给出最佳的模型参数。其中随机森林给出的最佳参数是Best params:  {'clf__criterion': 'entropy', 'clf__max_depth': 4, 'clf__n_estimators': 20}，分类的得分如下图：

knn给出的最佳参数是Best params:  {'clf__n_neighbors': 2, 'clf__weights': 'distance'}，分类的得分如下图：

对比两个模型可以看出KNN的效果更好。