1 项目简介
1.1 问题描述
淡水是我们最重要和最稀缺的自然资源之一，仅占地球总水量的 3%。它几乎触及我们日常生活的方方面面，从饮用、游泳和沐浴到生产食物、电力和我们每天使用的产品。获得安全卫生的供水不仅对人类生活至关重要，而且对正在遭受干旱、污染和气温升高影响的周边生态系统的生存也至关重要。题目来源于2023年intel的oneapi校园黑客松竞赛。
1.2 预期解决方案
通过参考英特尔的类似实现方案，预测淡水是否可以安全饮用和被依赖淡水的生态系统所使用，从而可以帮助全球水安全和环境可持续性发展。这里分类准确度和推理时间将作为评分的主要依据。
2 数据探索
先导入需要使用到的oneapi组件。
import modin.pandas as pd
import pandas
import os
os.environ["MODIN_ENGINE"] = "dask"
from modin.config import Engine
Engine.put("dask")
os.environ['NLSLANG']='SIMPLIFIED CHINESE_CHINA.UTF8'
pd.set_option( 'display.max_columns', 100)
warnings.filterwarnings('ignore')
pio.renderers.default='notebook'
intel_pal, color=['#0071C5','#FCBB13'], ['#7AB5E1','#FCE7B2']
temp=dict(layout=go.Layout(font=dict(family="Franklin Gothic", size=12),
height=500, width=1000))
Modin是一个Python第三方库，可以通过并行来处理大数据集。其中语法与pandas相近，拥有出色的性能弥补了pandas处理大型数据集的缺陷。
而Dask是一个用于分析计算的灵活的并行计算库，实现大型多维数据集分析的更快执行以及加速和扩展数据科学制作流程或工作流程的强大工具。
2.1 导入数据集
# Read data
data = pandas.read_csv('./data/dataset.csv')
test_data = pandas.read_csv('data/test_data.csv')
print("Data shape: {}\n".format(data.shape))
display(data.head())
data.info()
如下图所示，表格中共有24列的数据，19列为float64，2列为int64类型，剩余3列为非数值型数据。
2.2 因子化操作
非数值列不利于接下来的数据探索，将其因子化，转为数值型。
# 因子化操作，将非数值型转为数值型
display(data.head())
factor = pd.factorize(data['Color'])
print(factor)
data.Color = factor[0]
factor = pd.factorize(data['Source'])
print(factor)
data.Source = factor[0]
factor = pd.factorize(data['Month'])
print(factor)
data.Month = factor[0]
factor = pd.factorize(test_data['Color'])
test_data.Color = factor[0]
factor = pd.factorize(test_data['Source'])
test_data.Source = factor[0]
factor = pd.factorize(test_data['Month'])
test_data.Month = factor[0]
分别转换了Color、Source以及Month三列非数值型数据，转换结果如图所示。
2.3 数据统计性描述
为了方便对前述各项指标进行质量检查,包括对数据的分布情况、数据类型、最大值、最小值、均值、标准差等,另外重点查看变量的有效记录数的比例、完整百分比、有效记录数等作为字段筛选的依据
def display_stats(df):
"""
Function to display descriptive statistics of numerical variables,
includes skewness & kurtosis.
"""
df=data.describe()
skewness=data.skew()
kurtosis=data.kurtosis()
df=df.append([skewness, kurtosis],ignore_index=True)
idx=pd.Series(['count','mean','std','min','25%','50%','75%','max','skewness','kurtosis'],name='Summary Statistic')
df=pd.concat([df,idx],1).set_index('Summary Statistic')
display(df.style.format('{:,.3f}').
background_gradient(subset=(df.index[1:],df.columns[:]),
cmap='GnBu'))
display_stats(data)
分析得到的图像可知，部分变量的偏差值较大，需要进行处理，比如Iron，Lead等。
2.4 相关性分析
使用pandas库中的corr()函数计算数据框df中各个列（除了Target列）与Target列的相关系数，并使用abs()函数将其转换为绝对值，然后使用sort_values()函数按照相关系数的大小进行排序，得到一个Series类型的bar变量，其中每个元素代表一个特征与Target的相关系数的绝对值。
bar = data.corr()['Target'].abs().sort_values(ascending=False)[1:]
plt.bar(bar.index, bar, width=0.5)
# 设置figsize的大小
pos_list = np.arange(len(data.columns))
params = {
'figure.figsize': '20, 10'
}
plt.rcParams.update(params)
plt.xticks(bar.index, bar.index, rotation=-60, fontsize=10)
plt.show()
分析得到图像，从day开始的后几个相关性极低，接近于零。
将几个相关性低的变量删除，同时index明显是无关变量，也删除。
data = data.drop(
columns=['Index', 'Day', 'Time of Day', 'Month', 'Water Temperature', 'Source', 'Conductivity', 'Air Temperature'])
test_data = test_data.drop(
columns=['Index', 'Day', 'Time of Day', 'Month', 'Water Temperature', 'Source', 'Conductivity', 'Air Temperature'])
3 数据预处理
3.1 处理重复值与缺失值
先查看重复值和缺失值
display(data.isna().sum())
missing = data.isna().sum().sum()
duplicates = data.duplicated().sum()
print("\nThere are {:,.0f} missing values in the data.".format(missing))
print("There are {:,.0f} duplicate records in the data.".format(duplicates))
重复值与缺失值都比较多。
删除重复值并用众数填充缺失值。
from sklearn.impute import SimpleImputer
# 创建一个简单填充器，使用众数进行填充
imputer = SimpleImputer(strategy='most_frequent')
# 查找所有含有缺失值的列
columns_with_missing_values = data.columns[data.isnull().any()]
# 对于所有含有缺失值的列使用众数进行填充
for column in columns_with_missing_values:
data[column] = imputer.fit_transform(data[[column]])
data = data.drop_duplicates()
columns_with_missing_values = test_data.columns[test_data.isnull().any()]
for column in columns_with_missing_values:
test_data[column] = imputer.fit_transform(test_data[[column]])
test_data = test_data.drop_duplicates()
3.2 偏差值处理
首先，使用scipy.stats库中的pearsonr()函数计算数据框df中各个列与最后一列的皮尔逊相关系数，并检查每个特征与最后一列的相关性是否显著（p值小于0.05），如果不显著，则将该特征从数据框中删除。
接下来，检查每个数值型特征的方差是否小于等于0.1，如果小于等于0.1，则将该特征从数据框中删除。
from scipy.stats import pearsonr
variables = data.columns
var = data.var()
numeric = data.columns
for i in range(0, len(var) - 1):
if var[i] <= 0.1: # 方差大于10％
print(variables[i])
data = data.drop(numeric[i], 1)
test_data = test_data.drop(numeric[i], 1)
variables = data.columns
for i in range(0, len(variables)):
x = data[variables[i]]
y = data[variables[-1]]
if pearsonr(x, y)[1] > 0.05:
print(variables[i])
data = data.drop(variables[i], 1)
test_data = test_data.drop(variables[i], 1)
variables = data.columns
print(variables)
print(len(variables))
删除了lead列，剩余15列数据。
最后查看结果，没有缺失值和重复值了。
data = data.drop_duplicates()
test_data = test_data.drop_duplicates()
display(data.isna().sum())
missing = data.isna().sum().sum()
duplicates = data.duplicated().sum()
print("\nThere are {:,.0f} missing values in the data.".format(missing))
print("There are {:,.0f} duplicate records in the data.".format(duplicates))
3.3 数据不平衡
3.3.1 数据平衡检查
先查看正例与负例的比例，看是否平衡
print(data.Target.value_counts())
target = data.Target.value_counts()
target.rename(index={1: 'Target 1', 0: 'Target 0'}, inplace=True)
plt.pie(target, [0, 0.05], target.index, autopct='%1.1f%%')
plt.show()
结果如下图
算是比较平衡的，但是还可以进一步优化。
3.3.2 SMOTE算法
由于之前有些参数进行了因子化操作，所以不能直接使用smote方法。首先使用 One-Hot 编码对因子化的参数进行处理，然后再将编码后的数据与原始数据合并。这样就可以确保因子化的参数也能够被正确处理并参与到 SMOTE 方法中。
from imblearn.over_sampling import SMOTE
from sklearn.preprocessing import OneHotEncoder
# 对因子化的参数进行 One-Hot 编码
encoder = OneHotEncoder(sparse=False)
factorized_param = encoder.fit_transform(data[['Factorized_Param']])
data_encoded = pd.concat([data.drop(['Factorized_Param'], axis=1), pd.DataFrame(factorized_param)], axis=1)
smote_model = SMOTE(k_neighbors=2, random_state=42, sampling_strategy=0.45)
x_smote, y_smote = smote_model.fit_resample(data_encoded.iloc[:, :-1], data_encoded['Target'])
df_smote = pd.concat([x_smote, y_smote], axis=1)
# 对测试数据进行相同的处理
test_factorized_param = encoder.transform(test_data[['Factorized_Param']])
test_data_encoded = pd.concat([test_data.drop(['Factorized_Param'], axis=1), pd.DataFrame(test_factorized_param)], axis=1)
x_smote_test, y_smote_test = smote_model.fit_resample(test_data_encoded.iloc[:, :-1], test_data_encoded['Target'])
test_smote = pd.concat([x_smote_test, y_smote_test], axis=1)
df_smote.groupby('Target').count()
此时两类的比例为2:1
3.3.3 过拟合与欠拟合
先对插值后异常数据进行过拟合然后在单独对正常数据进行欠拟合。
from imblearn.over_sampling import RandomOverSampler
from imblearn.under_sampling import RandomUnderSampler
rus = RandomUnderSampler(random_state=0,sampling_strategy=0.45)
ros = RandomOverSampler(random_state=0,sampling_strategy=0.45)
X_resampled, y_resampled = ros.fit_resample(data.iloc[:,:-1],data['Target'])
df_resampled = pd.concat([X_resampled, y_resampled], axis=1)
X_uresampled, y_uresampled = rus.fit_resample(data.iloc[:,:-1],data['Target'])
df_uresampled = pd.concat([X_uresampled, y_uresampled], axis=1)
df_new = pd.concat([df_smote, df_resampled[df_resampled['Target']==1.0]], axis=0)
df_new = pd.concat([df_new[df_new['Target']==1.0],df_uresampled[df_uresampled['Target']==0.0]])
X_resampled, y_resampled = ros.fit_resample(test_data.iloc[:,:-1],test_data['Target'])
df_resampled = pd.concat([X_resampled, y_resampled], axis=1)
X_uresampled, y_uresampled = rus.fit_resample(test_data.iloc[:,:-1],test_data['Target'])
df_uresampled = pd.concat([X_uresampled, y_uresampled], axis=1)
test_new = pd.concat([test_smote, df_resampled[df_resampled['Target']==1.0]], axis=0)
test_new = pd.concat([test_new[test_new['Target']==1.0],df_uresampled[df_uresampled['Target']==0.0]])
df_new.groupby('Target').count()
两类数据比例已经接近1:1。
4 模型拟合
4.1 数据可视化函数
下面是需要用的几个函数，提前定义，分别用于绘制绘制ROC/PR曲线和预测目标分布、绘制预测数据分布
import plotly.io as pio
pio.renderers.default = "iframe"
def plot_model_res(model_name, y_test, y_prob):
intel_pal=['#0071C5','#FCBB13']
color=['#7AB5E1','#FCE7B2']
fpr, tpr, _ = roc_curve(y_test, y_prob)
roc_auc = auc(fpr,tpr)
precision, recall, _ = precision_recall_curve(y_test, y_prob)
auprc = average_precision_score(y_test, y_prob)
fig = make_subplots(rows=1, cols=2,
shared_yaxes=True,
subplot_titles=['Receiver Operating Characteristic<br>(ROC) Curve',
'Precision-Recall Curve<br>AUPRC = {:.3f}'.format(auprc)])
fig.add_trace(go.Scatter(x=np.linspace(0,1,11), y=np.linspace(0,1,11),
name='Baseline',mode='lines',legendgroup=1,
line=dict(color="Black", width=1, dash="dot")), row=1,col=1)
fig.add_trace(go.Scatter(x=fpr, y=tpr, line=dict(color=intel_pal[0], width=3),
hovertemplate = 'True positive rate = %{y:.3f}, False positive rate = %{x:.3f}',
name='AUC = {:.4f}'.format(roc_auc),legendgroup=1), row=1,col=1)
fig.add_trace(go.Scatter(x=recall, y=precision, line=dict(color=intel_pal[0], width=3),
hovertemplate = 'Precision = %{y:.3f}, Recall = %{x:.3f}',
name='AUPRC = {:.4f}'.format(auprc),showlegend=False), row=1,col=2)
fig.update_layout(template=temp, title="{} ROC and Precision-Recall Curves".format(model_name),
hovermode="x unified", width=900,height=500,
xaxis1_title='False Positive Rate (1 - Specificity)',
yaxis1_title='True Positive Rate (Sensitivity)',
xaxis2_title='Recall (Sensitivity)',yaxis2_title='Precision (PPV)',
legend=dict(orientation='v', y=.07, x=.45, xanchor="right",
bordercolor="black", borderwidth=.5))
fig.show()
def plot_distribution(y_prob):
plot_df=pd.DataFrame.from_dict({'Target 0':(len(y_prob[y_prob<=0.5])/len(y_prob))*100,
'Target 1':(len(y_prob[y_prob>0.5])/len(y_prob))*100},
orient='index', columns=['pct'])
fig=go.Figure()
fig.add_trace(go.Pie(labels=plot_df.index, values=plot_df.pct, hole=.45,
text=plot_df.index, sort=False, showlegend=False,
marker=dict(colors=color,line=dict(color=intel_pal,width=2.5)),
hovertemplate = "%{label}: <b>%{value:.2f}%</b><extra></extra>"))
fig.update_layout(template=temp, title='Predicted Target Distribution',width=700,height=450,
uniformtext_minsize=15, uniformtext_mode='hide')
fig.show()
4.2 XGBClassifier拟合
import time
## Prepare Train and Test datasets ##
scaler = RobustScaler()
X_train = scaler.fit_transform(df_new.drop(['Target'], axis=1))
y_train = df_new['Target']
print("Train Shape: {}".format(df_new.shape))
## Initialize XGBoost model ##
ratio = float(np.sum(y_train == 0)) / np.sum(y_train == 1)
parameters = {'scale_pos_weight': ratio.round(2),
'tree_method': 'hist',
'random_state': 21}
xgb_model = XGBClassifier(**parameters)
## Tune hyperparameters ##
strat_kfold = StratifiedKFold(n_splits=3, shuffle=True, random_state=21)
print("\nTuning hyperparameters..")
grid = {'min_child_weight': [1, 5, 10],
'gamma': [0.5, 1, 1.5, 2, 5],
'max_depth': [3, 4, 5],
}
grid_search = GridSearchCV(xgb_model, param_grid=grid,
cv=strat_kfold, scoring='roc_auc',
verbose=1, n_jobs=-1)
grid_search.fit(X_train, y_train)
print("Done!\nBest hyperparameters:", grid_search.best_params_)
print("Best cross-validation AUC: {:.4f}".format(grid_search.best_score_))
## Convert XGB model to daal4py ##
xgb = grid_search.best_estimator_
daal_model = d4p.get_gbt_model_from_xgboost(xgb.get_booster())
拟合结果
4.3 测试集测试
X_test = scaler.transform(test_new.drop(['Target'], axis=1))
y_test = test_new['Target']
print("Test Shape: {}".format(test_new.shape))
start_time = time.time()
daal_prob = d4p.gbt_classification_prediction(nClasses=2,
resultsToEvaluate="computeClassLabels|computeClassProbabilities",
fptype='float').compute(X_test, daal_model).probabilities # or .predictions
end_time = time.time()
inference_time_test = end_time - start_time
xgb_pred = pd.Series(np.where(daal_prob[:,1]>.5, 1, 0), name='Target')
xgb_auc = roc_auc_score(y_test, daal_prob[:,1])
xgb_f1 = f1_score(y_test, xgb_pred)
## Plot model results ##
print("\nTest F1 Accuracy: {:.2f}%, AUC: {:.5f}".format(xgb_f1*100, xgb_auc))
plot_model_res(model_name='XGBoost', y_test=y_test, y_prob=daal_prob[:,1])
plot_distribution(daal_prob[:,1])
print('预测类别为正常的数量为:{:d}'.format(len(daal_prob[:,1][daal_prob[:,1]<=0.5])))
print('预测类别为异常的数量为:{:d}'.format(len(daal_prob[:,1][daal_prob[:,1]>=0.5])))
print("Inference Time on Test Set: {:.4f} seconds".format(inference_time_test))
最后的结果是F1分数为0.8923，推理速度为0.1558秒
扫一扫
205
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
