文章目录
SVC真实数据案例：预测明天是否会下雨
一般来说很少能一次性导入所有需要的库，可以一边码一边添
import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import LabelEncoder
from math import radians, sin, cos, acos
import sys
import re
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import OrdinalEncoder
from sklearn.preprocessing import StandardScaler
from time import time
import datetime
from sklearn.svm import SVC
from sklearn.model_selection import cross_val_score
from sklearn.metrics import roc_auc_score, recall_score
from sklearn.metrics import confusion_matrix as CM
import matplotlib.pyplot as plt
from sklearn.metrics import roc_curve as ROC
from sklearn.metrics import accuracy_score as AC
weather = pd.read_csv("weatherAUS5000.csv", index_col=0)
print(weather.head())
print(weather.info())
#探索数据集特征
object是文本类型数据集，一共有22个特征（字段）大部分字段有缺失值
#定义特征值x和标签y
x = weather.iloc[:, :-1]
y = weather.iloc[:, -1]
print(x.shape)
#(5000, 21)
print(y.shape)
#(5000,)
print(np.unique(y))
#观察y的标签类型有几种
#array(['No', 'Yes'], dtype=object)
print(x.isnull().mean())
#查看每一列的缺失值占比
缺失率低的字段可以不进行填补
在实际情况下，通常很难得到测试集，所以我们在这个案例中提前将训练集和测试集进行拆分，拆分之后要更新索引
xtrain, xtest, ytrain, ytest = train_test_split(x, y, test_size=0.3, random_state=420)
for i in [xtrain, xtest, ytrain, ytest]:
i.index = range(i.shape[0])
查看是否存在样本不均衡的问题
y_label_count = ytrain.value_counts()
print(y_label_count)
可以看出来yes是少数类，no是多数类，少数类的标签是1，多数类的标签是0，多少有一点样本不均衡，比例基本上是在3:1左右
#对标签进行编码
encorder = LabelEncoder().fit(ytrain)
y_train = pd.DataFrame(encorder.transform(ytrain))
y_test = pd.DataFrame(encorder.transform(ytest))
观察特征中具体的值，是否存在不符合常理的情况
xtrain.describe([0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.99]).T
xtest.describe([0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.99]).T
其中温度有零下摄氏度的情况，有负数是正常表现。没有特别严重的偏态问题，就是量纲不统一，一会一起处理。
测试集和训练集的分布比较接近，没有特别严重的问题
处理困难特征：时间
接下来处理最难的两个特征：时间和地点，因为时间和地点在一定程度上影响降雨量，但是我们又没有办法让计算机直接了解到这一点，所以需要对这两个特征进行转换。首先，我们需要把时间转换成数字型；其次，还需要将时间转化为月份，因为不同的月份所处季节不同，降雨的概率也不同。地点就转换为气候类型。
#首先我们来看一下时间，同一个时间是否只出现一次，如果是，那有可能是时序数据，如果不是，考虑是否能进行分类
xtrain.iloc[:, 0].value_counts()
可以发现每个时间并不是只出现一次。如果我们把它当作分类型变量处理，类别太多，有2141类，如果换成数值型，会被直接当成连续型变量，如果做成哑变量，我们特征的维度会爆炸。那我们可以换一种思路，既然算法处理的是列与列之间的关系，我是否可以把”今天的天气会影响明天的天气“这个指标转换成一个特征呢？因此，我们可以将时间对气候的连续影响，转换为”今天是否下雨“这个特征，巧妙地将样本对应标签之间的联系，转换成是特征与标签之间的联系 了。
#增加特征值：今天是否会下雨，将降雨量大于1的定义为会，其余的定义为不会
xtrain.loc[xtrain["Rainfall"] >= 1, "Rainday"] = "yes"
xtrain.loc[xtrain["Rainfall"] < 1, "Rainday"] = "no"
xtrain.loc[xtrain["Rainfall"] == np.nan, "Rainday"] = np.nan
xtest.loc[xtrain["Rainfall"] >= 1, "Rainday"] = "yes"
xtest.loc[xtrain["Rainfall"] < 1, "Rainday"] = "no"
xtest.loc[xtrain["Rainfall"] == np.nan, "Rainday"] = np.nan
其次，将时间中的月份进行提取
xtrain.loc[0,"Date"].split("-")
#['2015', '08', '24']
xtrain.loc[0,"Date"].split("-")[1]
#'08'
int(xtrain.loc[0,"Date"].split("-")[1])
#08
xtrain["Date"] = xtrain["Date"].apply(lambda x:int(x.split("-")[1]))
xtest["Date"] = xtest["Date"].apply(lambda x:int(x.split("-")[1]))
xtrain = xtrain.rename(columns={"Date": "Month"})
xtest = xtest.rename(columns={"Date": "Month"})
处理困难特征：地点
常识上来说，我们认为地点肯定是对明天是否会下雨存在影响的。比如说，如果其他信息都不给出，我们只猜测，“伦敦明天是否会下雨”和”北京明天是否会下雨“，我一定会猜测伦敦会下雨，而北京不会，因为伦敦是常年下雨的城市，而北京的气候非常干燥。对澳大利亚这样面积巨大的国家来说，必然存在着不同的城市有着不同的下雨倾向的情况。但尴尬的是，和时间一样，我们输入地点的名字对于算法来说，就是一串字符，"London"和"Beijing"对算法来说，和0，1没有区别。同样，我们的样本中含有49个不同地点，如果做成分类型变量，算法就无法辨别它究竟是否是分类变量。也就是说，我们需要让算法意识到，不同的地点因为气候不同，所以对“明天是否会下雨”有着不同的影响。如果我们能够将地点转换为这个地方的气候的话，我们就可以将不同城市打包到同一个气候中，而同一个气候下反应的降雨情况应该是相似的。
这是由澳大利亚气象局和澳大利亚建筑规范委员会（ABCB）制作统计的，澳大利亚不同地区不同城市的所在的气候区域划分。总共划分为八个区域，非常适合我们用来做分类。如果能够把49个地点转换成八种不同的气候，这个信息应该会对是否下雨的判断比较有用。基于气象局和ABCB的数据，我为大家制作了澳大利亚主要城市所对应的气候类型数据，并保存在csv文件city_climate.csv当中。然后，我使用以下代码，在google上进行爬虫，爬出了每个城市所对应的经纬度，并保存在数据cityll.csv当中，大家可以自行导入，来查看这个数据。爬虫的笔记之后再补充！
city_climate = pd.read_csv("Cityclimate.csv")
city_climate.head()
这个表格里分别是城市与其对应的气候
接下来爬取澳大利亚每个城市的经纬度，之所以这么做，是因为数据的location是气候站的名称，而不是城市本身的名称，因此不是每一个城市都能直接获取到其气候，需要将气候站的经纬度和城市的经纬度进行距离计算，找到与其相近的城市地点，将其气候类型进行对应。
cityll = pd.read_csv("cityll.csv", index_col=0)
cityll
float(cityll.loc[0, "Latitude"][: -1])#去掉度数，变成浮点数
cityll["Latitude"] = cityll["Latitude"].apply(lambda x: float(x[: -1]))
cityll["Longitude"] = cityll["Longitude"].apply(lambda x: float(x[: -1]))
#观察一下所有的经纬度方向都是一致的，全部是南纬，东经，因为澳大利亚在南半球，东半球
#所以经纬度的方向我们可以舍弃了
citylld = cityll.iloc[:, [0, 1, 2]]
#将气候类型进行填补
citylld["Climate"] = city_climate.iloc[:, -1]
爬取我们原始数据的气候站名称即经纬度
samplecity = pd.read_csv("samplecity.csv", index_col=0)
samplecity["Latitude"] = samplecity["Latitude"].apply(lambda x: float(x[: -1]))
samplecity["Longitude"] = samplecity["Longitude"].apply(lambda x: float(x[: -1]))
samplecity.head()
samplecityd = samplecity.iloc[:, [0, 1, 2]] #舍去经纬度方向
#首先使用radians将角度转换成弧度
citylld.loc[:, "slat"] = citylld.iloc[:, 1].apply(lambda x : radians(x))
citylld.loc[:, "slon"] = citylld.iloc[:, 2].apply(lambda x : radians(x))
samplecityd.loc[:,"elat"] = samplecityd.iloc[:,1].apply(lambda x : radians(x))
samplecityd.loc[:,"elon"] = samplecityd.iloc[:,2].apply(lambda x : radians(x))
for i in range(samplecityd.shape[0]):
slat = citylld.loc[:, "slat"]
slon = citylld.loc[:, "slon"]
elat = samplecityd.loc[i, "elat"]
elon = samplecityd.loc[i, "elon"]
dist = 6371.01 * np.arccos(np.sin(slat)*np.sin(elat) +
np.cos(slat)*np.cos(elat)*np.cos(slon.values - elon))
city_index = np.argsort(dist)[0]
#每次计算后，取距离最近的城市，然后将最近的城市和城市对应的气候都匹配到samplecityd中
samplecityd.loc[i, "closest_city"] = citylld.loc[city_index, "City"]
samplecityd.loc[i, "climate"] = citylld.loc[city_index, "Climate"]
#确认无误后需要取出样本城市所对应的气候并保存
locafinal = samplecityd.iloc[:, [0, -1]]
#在这里设定locafinal的索引为地点，是为了之后进行map的匹配
locafinal.columns = ["Location", "Climate"]
locafinal = locafinal.set_index(keys="Location")
locafinal.to_csv(r"D:\samplelocation2.csv")
print(locafinal.head())
#将location中的内容替换，并且确保匹配进入的气候字符串中不含有逗号，气候两边不含有空格
#我们使用re这个模块来消除逗号
#re.sub（希望替换的值，希望被替换成的值，要操作的字符串）
xtrain["Location"] = xtrain["Location"].map(locafinal.iloc[:, 0]).apply(lambda x: re.sub(",", "", x.strip()))
xtest["Location"] = xtest["Location"].map(locafinal.iloc[:, 0]).apply(lambda x: re.sub(",", "", x.strip()))
#修改特征内容之后，我们使用新列名“Climate”来替换之前的列名“Location”
#注意这个命令一旦执行之后，就再没有列"Location"了，使用索引时要特别注意
xtrain = xtrain.rename(columns={"Location": "Climate"})
xtest = xtest.rename(columns={"Location": "Climate"})
处理分类变量
用众数填补缺失值
#首先找出，分类型特征有哪些
cate = xtrain.columns[xtrain.dtypes == "object"].tolist()
#除了特征类型为“object”的特征们，还有虽然用数字表示，但是本质为分类型特征的云层遮蔽程度
cloud = ["Cloud9am", "Cloud3pm"]
cate = cate + cloud
#对于分类型特征，我们使用众数进行填补
si = SimpleImputer(missing_values=np.nan, strategy="most_frequent")
si.fit(xtrain.loc[:, cate])
#然后我们用训练集中的众数来同时填补训练集和测试集
xtrain.loc[:, cate] = si.transform(xtrain.loc[:, cate])
xtest.loc[:, cate] = si.transform(xtest.loc[:, cate])
xtrain.loc[:, cate].isnull().mean()
xtest.loc[:, cate].isnull().mean()
进行编码
#处理分类变量：将分类型变量编码
oe = OrdinalEncoder().fit(xtrain.loc[:, cate])
#用训练集的编码结果来编码训练和测试矩阵
xtrain.loc[:, cate] = oe.transform(xtrain.loc[:, cate])
xtest.loc[:, cate] = oe.transform(xtest.loc[:, cate])
处理连续型变量
#处理连续型变量：填补缺失值
col = xtrain.columns.tolist()
for i in cate:
col.remove(i)
#实例化，填补策略为“mean”均值
impmean = SimpleImputer(missing_values = np.nan, strategy="mean").fit(xtrain.loc[:, col])
xtrain.loc[:, col] = impmean.transform(xtrain.loc[:, col])
xtest.loc[:, col] = impmean.transform(xtest.loc[:, col])
ss = StandardScaler().fit(xtrain.loc[:, col])
xtrain.loc[:, col] = ss.transform(xtrain.loc[:, col])
xtest.loc[:, col] = ss.transform(xtest.loc[:, col])
统一量纲：标准化
将数据均值调整为0，方差为1，但是数据并不一定服从正太分布
ss = StandardScaler().fit(xtrain.loc[:, col])
xtrain.loc[:, col] = ss.transform(xtrain.loc[:, col])
xtest.loc[:, col] = ss.transform(xtest.loc[:, col])
建模与模型评估
选择最佳的kernel
times = time()
for kernel in ["linear", "poly", "rbf", "sigmoid"]:
clf = SVC(kernel=kernel, gamma="auto", degree=1, cache_size=5000, class_weight = "balanced").fit(xtrain, y_train)
result = clf.predict(xtest)
score = clf.score(xtest, y_test)
recall = recall_score(y_test, result)
auc = roc_auc_score(y_test, clf.decision_function(xtest))
print("%s 's testing accuracy %f, recall is %f', auc is %f" %(kernel, score, recall, auc))
print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
可以看出rbf的效果最好，我仿佛眼瞎了，用的全都是linear，但是换汤不换药，我先把流程过一遍。选出最佳的kernel之后，选择其他最佳参数
模型调参
追求最高Recall
如果想要recall分数高，那就要加大少数类的权重，让模型更多得学习少数类
times = time()
for kernel in ["linear", "poly", "rbf", "sigmoid"]:
clf = SVC(kernel=kernel, gamma="auto", degree=1, cache_size=5000, class_weight={1: 10}#注意，这里写的其实是，类别1：10，隐藏了类别0：1这个比例).fit(xtrain, y_train)
result = clf.predict(xtest)
score = clf.score(xtest, y_test)
recall = recall_score(y_test, result)
auc = roc_auc_score(y_test, clf.decision_function(xtest))
print("%s 's testing accuracy %f, recall is %f', auc is %f" %(kernel, score, recall, auc))
print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
当把标签1的权重调整为10的时候，可以看出linear效果最佳
追求高准确率
先查看特异度，了解一下多数类被预测对的比例有多大
clf = SVC(kernel="linear",
gamma="auto",
cache_size=5000).fit(xtrain, y_train)
result = clf.predict(xtest)
cm = CM(y_test, result, labels=(1, 0))
#print(cm)
'''
[[ 149 194]
[ 47 1110]]
'''
specificity = cm[1, 1]/ cm[1, :].sum()
print(specificity)
#0.9593777009507347,特异度很高，说明大多数的多数类已经被准确预测到了
了解到特异度的值很大，说明大多数多数类已经预测对了，可以试试看使用class_weight将模型向少数类的方向稍微调整，已查看我们是否有更多的空间来捕捉少数类从而提升精确度
irange = np.linspace(0.01, 0.05, 10)
for i in irange:
times = time()
clf = SVC(kernel="linear",
gamma="auto",
cache_size=5000,
class_weight={1: 1+i}).fit(xtrain, y_train)
result = clf.predict(xtest)
score = clf.score(xtest, y_test)
recall = recall_score(y_test, result)
auc = roc_auc_score(y_test, clf.decision_function(xtest))
print("under ratio 1:%f testing accuracy %f, recall is %f', auc is %f" %
(1 + i, score, recall, auc))
print(datetime.datetime.fromtimestamp(time() - times).strftime("%M:%S:%f"))
整体来看，准确率是在0.84左右徘徊
追求平衡
我们前面经历了多种尝试，选定了线性核，并发现调节class_weight并不能够使我们模型有较大的改善。现在我们来试试看调节线性核函数的C值能否有效果：
c_range = np.linspace(0.01, 20, 20)
recallall = []
aucall = []
scoreall = []
for c in c_range:
times = time()
clf = SVC(kernel="linear", C=c, cache_size=5000, class_weight="balanced").fit(xtrain, y_train)
result = clf.predict(xtest)
score = clf.score(xtest, y_test)
recall = recall_score(y_test, result)
auc = roc_auc_score(y_test, clf.decision_function(xtest))
recallall.append(recall)
aucall.append(auc)
scoreall.append(score)
print("under C %f, testing accuracy is %f,recall is %f', auc is %f" %
(c, score, recall, auc))
print(datetime.datetime.fromtimestamp(time() - times).strftime("%M:%S:%f"))
print(max(aucall), c_range[aucall.index(max(aucall))])
plt.figure()
plt.plot(c_range, recallall, c="red", label="recall")
plt.plot(c_range, aucall, c="black", label="auc")
plt.plot(c_range, scoreall, c="orange", label="accuracy")
plt.legend()
plt.show()
可以看出调c值还不如不调的效果好嘞
画ROC曲线，找到最佳阈值，通过调整阈值来判断准确率和召回率
clf = SVC(kernel="linear", gamma="auto", degree=1, cache_size=5000, class_weight="balanced").fit(xtrain, y_train)
FPR, Recall, thresholds = ROC(y_test, clf.decision_function(xtest), pos_label=1)
area = roc_auc_score(y_test, clf.decision_function(xtest))
plt.figure()
plt.plot(FPR, Recall, color='red',
label='ROC curve (area = %0.2f)' % area)
plt.plot([0, 1], [0, 1], color='black', linestyle='--')
plt.xlim([-0.05, 1.05])
plt.ylim([-0.05, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('Recall')
plt.title('Receiver operating characteristic example')
plt.legend(loc="lower right")
plt.show()
maxindex = (Recall - FPR).tolist().index(max(Recall - FPR))
best_threshold = thresholds[maxindex]
print(best_threshold)
#-0.12138496714357605
times = time()
clf = SVC(kernel="linear",
cache_size=5000,
class_weight="balanced").fit(xtrain, y_train)
prob = pd.DataFrame(clf.decision_function(xtest))
prob.loc[prob.iloc[:, 0] >= thresholds[maxindex], "y_pred"] = 1
prob.loc[prob.iloc[:, 0] < thresholds[maxindex], "y_pred"] = 0
prob.loc[:, "y_pred"].isnull().sum()
#检查模型本身的准确度
score = AC(y_test, prob.loc[:, "y_pred"].values)
recall = recall_score(y_test, prob.loc[:, "y_pred"])
print("testing accuracy %f,recall is %f" % (score, recall))
print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
testing accuracy 0.790667,recall is 0.807580
00:04:920161
调整了阈值之后，达到一个相对平衡的状态
总之，这个结果并不让人满意，如果是打比赛的话，远远不够
扫一扫
8857
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
