第一步:读取数据
我这里是把文件存在当前文件夹,新建的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的效果更好。