python做回归

逻辑回归

https://github/lawlite19/MachineLearning_Python/tree/master/LogisticRegression

全部代码

https://github/lawlite19/MachineLearning_Python/blob/master/LogisticRegression/LogisticRegression.py

代价函数

可以综合起来为：

其中：

为什么不用线性回归的代价函数表示，因为线性回归的代价函数可能是非凸的，对于分类问题，使用梯度下降很难得到最小值，上面的代价函数是凸函数

的图像如下，即y=1时：

可以看出，当

趋于1，y=1,与预测值一致，此时付出的代价cost趋于0，若

趋于0，y=1,此时的代价cost值非常大，我们最终的目的是最小化代价值

同理

的图像如下（y=0）：

梯度同样对代价函数求偏导：

可以看出与线性回归的偏导数一致

推导过程

正则化目的是为了防止过拟合

在代价函数中加上一项

注意j是重1开始的，因为theta(0)为一个常数项，X中最前面一列会加上1列1，所以乘积还是theta(0),feature没有关系，没有必要正则化

正则化后的代价：

# 代价函数def costFunction(initial_theta,X,y,inital_lambda): m = len(y) J = 0 h = sigmoid(np.dot(X,initial_theta)) # 计算h(z) theta1 = initial_theta.copy() # 因为正则化j=1从1开始，不包含0，所以复制一份，前theta(0)值为0 theta1[0] = 0 temp = np.dot(np.transpose(theta1),theta1) J = (-np.dot(np.transpose(y),np.log(h))-np.dot(np.transpose(1-y),np.log(1-h))+temp*inital_lambda/2)/m # 正则化的代价方程 return J

正则化后的代价的梯度

# 计算梯度def gradient(initial_theta,X,y,inital_lambda): m = len(y) grad = np.zeros((initial_theta.shape[0])) h = sigmoid(np.dot(X,initial_theta))# 计算h(z) theta1 = initial_theta.copy() theta1[0] = 0 grad = np.dot(np.transpose(X),h-y)/m+inital_lambda/m*theta1 #正则化的梯度 return grad

S型函数（即

）实现代码

# S型函数def sigmoid(z): h = np.zeros((len(z),1)) # 初始化，与z的长度一置 h = 1.0/(1.0+np.exp(-z)) return h

映射为多项式因为数据的feture可能很少，导致偏差大，所以创造出一些feture结合

eg:映射为2次方的形式:

实现代码：

# 映射为多项式def mapFeature(X1,X2): degree = 3; # 映射的最高次方 out = np.ones((X1.shape[0],1)) # 映射后的结果数组（取代X） ''' 这里以degree=2为例，映射为1,x1,x2,x1^2,x1,x2,x2^2 ''' for i in np.arange(1,degree+1): for j in range(i+1): temp = X1**(i-j)*(X2**j) #矩阵直接乘相当于matlab中的点乘.* out = np.hstack((out, temp.reshape(-1,1))) return out

使用scipy的优化方法梯度下降使用scipy中optimize中的fmin_bfgs函数

调用scipy中的优化算法fmin_bfgs（拟牛顿法Broyden-Fletcher-Goldfarb-Shanno

costFunction是自己实现的一个求代价的函数，

initial_theta表示初始化的值,

fprime指定costFunction的梯度

args是其余测参数，以元组的形式传入，最后会将最小化costFunction的theta返回

result = optimize.fmin_bfgs(costFunction, initial_theta, fprime=gradient, args=(X,y,initial_lambda))

运行结果data1决策边界和准确度

data2决策边界和准确度

使用scikit-learn库中的逻辑回归模型实现

https://github/lawlite19/MachineLearning_Python/blob/master/LogisticRegression/LogisticRegression_scikit-learn.py

导入包

from sklearn.linear_model import LogisticRegressionfrom sklearn.preprocessing import StandardScalerfrom sklearn.cross_validation import train_test_splitimport numpy as np

划分训练集和测试集

x_train,x_test,y_train,y_test = train_test_split(X,y,test_size=0.2)

归一化

scaler = StandardScaler() scaler.fit(x_train) x_train = scaler.fit_transform(x_train) x_test = scaler.fit_transform(x_test)

逻辑回归

model = LogisticRegression() model.fit(x_train,y_train)

预测

predict = model.predict(x_test) right = sum(predict == y_test) predict = np.hstack((predict.reshape(-1,1),y_test.reshape(-1,1))) # 将预测值和真实值放在一块，好观察 print predict print ('测试集准确率：%f%%'%(right*100.0/predict.shape[0])) #计算在测试集上的准确度

逻辑回归_手写数字识别_OneVsAll

https://github/lawlite19/MachineLearning_Python/blob/master/LogisticRegression

全部代码

https://github/lawlite19/MachineLearning_Python/blob/master/LogisticRegression/LogisticRegression_OneVsAll.py

随机显示100个数字我没有使用scikit-learn中的数据集，像素是20*20px，彩色图如下

灰度图：

实现代码：

# 显示100个数字def display_data(imgData): sum = 0 ''' 显示100个数（若是一个一个绘制将会非常慢，可以将要画的数字整理好，放到一个矩阵中，显示这个矩阵即可） - 初始化一个二维数组 - 将每行的数据调整成图像的矩阵，放进二维数组 - 显示即可 ''' pad = 1 display_array = -np.ones((pad+10*(20+pad),pad+10*(20+pad))) for i in range(10): for j in range(10): display_array[pad+i*(20+pad):pad+i*(20+pad)+20,pad+j*(20+pad):pad+j*(20+pad)+20] = (imgData[sum,:].reshape(20,20,order="F")) # order=F指定以列优先，在matlab中是这样的，python中需要指定，默认以行 sum += 1 plt.imshow(display_array,cmap='gray') #显示灰度图像 plt.axis('off') plt.show()

OneVsAll如何利用逻辑回归解决多分类的问题，OneVsAll就是把当前某一类看成一类，其他所有类别看作一类，这样有成了二分类的问题了

如下图，把途中的数据分成三类，先把红色的看成一类，把其他的看作另外一类，进行逻辑回归，然后把蓝色的看成一类，其他的再看成一类，以此类推...

可以看出大于2类的情况下，有多少类就要进行多少次的逻辑回归分类

手写数字识别共有0-9，10个数字，需要10次分类

由于数据集y给出的是0,1,2...9的数字，而进行逻辑回归需要0/1的label标记，所以需要对y处理

说一下数据集，前500个是0,500-1000是1,...,所以如下图，处理后的y，前500行的第一列是1，其余都是0,500-1000行第二列是1，其余都是0....

然后调用梯度下降算法求解theta

实现代码：

# 求每个分类的theta，最后返回所有的all_theta def oneVsAll(X,y,num_labels,Lambda): # 初始化变量 m,n = X.shape all_theta = np.zeros((n+1,num_labels)) # 每一列对应相应分类的theta,共10列 X = np.hstack((np.ones((m,1)),X)) # X前补上一列1的偏置bias class_y = np.zeros((m,num_labels)) # 数据的y对应0-9，需要映射为0/1的关系 initial_theta = np.zeros((n+1,1)) # 初始化一个分类的theta # 映射y for i in range(num_labels): class_y[:,i] = np.int32(y==i).reshape(1,-1) # 注意reshape(1,-1)才可以赋值 #np.savetxt("class_y.csv", class_y[0:600,:], delimiter=',') '''遍历每个分类，计算对应的theta值''' for i in range(num_labels): result = optimize.fmin_bfgs(costFunction, initial_theta, fprime=gradient, args=(X,class_y[:,i],Lambda)) # 调用梯度下降的优化方法 all_theta[:,i] = result.reshape(1,-1) # 放入all_theta中 all_theta = np.transpose(all_theta) return all_theta

预测之前说过，预测的结果是一个概率值，利用学习出来的theta代入预测的S型函数中，每行的最大值就是是某个数字的最大概率，所在的列号就是预测的数字的真实值,因为在分类时，所有为0的将y映射在第一列，为1的映射在第二列，依次类推

实现代码：

# 预测def predict_oneVsAll(all_theta,X): m = X.shape[0] num_labels = all_theta.shape[0] p = np.zeros((m,1)) X = np.hstack((np.ones((m,1)),X)) #在X最前面加一列1 h = sigmoid(np.dot(X,np.transpose(all_theta))) #预测 ''' 返回h中每一行最大值所在的列号 - np.max(h, axis=1)返回h中每一行的最大值（是某个数字的最大概率） - 最后where找到的最大概率所在的列号（列号即是对应的数字） ''' p = np.array(np.where(h[0,:] == np.max(h, axis=1)[0])) for i in np.arange(1, m): t = np.array(np.where(h[i,:] == np.max(h, axis=1)[i])) p = np.vstack((p,t)) return p

运行结果10次分类，在训练集上的准确度：

使用scikit-learn库中的逻辑回归模型实现https://github/lawlite19/MachineLearning_Python/blob/master/LogisticRegression/LogisticRegression_OneVsAll_scikit-learn.py

导入包

from scipy import io as spioimport numpy as npfrom sklearn import svmfrom sklearn.linear_model import LogisticRegression

加载数据

data = loadmat_data("data_digits.mat") X = data['X'] # 获取X数据，每一行对应一个数字20x20px y = data['y'] # 这里读取mat文件y的shape=(5000, 1) y = np.ravel(y) # 调用sklearn需要转化成一维的(5000,)

拟合模型

model = LogisticRegression() model.fit(X, y) # 拟合

预测

predict = model.predict(X) #预测 print u"预测准确度为：%f%%"%np.mean(np.float64(predict == y)*100)

输出结果（在训练集上的准确度）

（未完待续）

线性回归

（1）线性回归的主要内容：

有兴趣可以玩一下这个游戏：是猜相关系数的，给你一些散点图，猜相关系数，很难猜对，说明看图说明相关性不靠谱！

（2）线性回归怎么做？数学公式

一个简单线性回归的例子：

ols：

拟合优度：

（3） 假设检验

线性回归这位老师用的是statsmodels做的。

这就是用线性回归拟合出来的模型！

使用这个函数可以直接得到模型的一些结果，会得到这些表！

这个是预测的函数：predict（x）。这是这个包的函数！

二、多元线性回归

也可以用向量的方式。

参数是迭代计算的，此处不讲了。

模型的检验

只用与 模型的选择 调整的r方。避免的样本量的影响，预防过拟合！

python中的实现 试用stats model 多个变量用+连一起！

模型的一些值

变量的筛选，除了常规的，就是使用相前，向后还有逐步，还有全子集的方法。

不多介绍，上篇也讲过。

python需要自己写来实现，当然我不会，这里可以用spss去实现

这是python代码

跳过。。。

回归的假设。

线性

强影响点筛选，然后直接删除

强烈建议用spss进行多元线性回归建模。不过为了熟悉python明天打算用python更着去实现。

自己整理一下，多元线性回归的检验方法：

整体显著性检验：回归模型出来时候，可以之间看整体的F检验的p值

拟合优度检验：最简单的就看R方和调整的R方

系数检验：T检验结果，可以直接看P值

自相关检验：dw检验用得最多

异方差检验：Goldfeld - Quandt 检验法，这个方法挺好的的。其实还有其他检验方法，解决办法是建议取对数

共线性检验：vif值

不够细，明日具体操作的时候会细点！

翻译：猫儿爷/小兵/JuneXie审核：JuneXie

在2010年创建的Kaggle是一个数据科学平台，在此平台上用户可以共享、协作和竞争。Kaggle平台的一个主要特点是“竞争”，它为用户提供了在真实世界实践数据的能力，并在国际舞台中展现技能。

本指南将教你如何参加Kaggle比赛，包括探索数据、特性的创建和工程化（engineeringfeatures）、构建模型和提交预测。我们将用Python3和jupyter笔记本做。

用先进的回归技术做房价数据预测竞赛。

按照以下步骤取得成功的Kaggle参赛作品:

·获得数据

·探索数据

·特性的创建和工程化（engineeringfeatures）以及目标变量

·建立模型

·制作并提交预测

01

获取数据并创建我们的环境

我们需要为此次竞争找到数据。对特性和一些其他有用信息的描述包含在一个名称为data_description.txt的文件中。（附后）

下载数据并将其保存到一个文件夹中，在这个文件夹中，保存所有你需要的文件。

我们先来看看train.csv数据文件。在我们训练了模型之后，我们将预测结果存入test.csv数据文件。

首先，importPandas,一个很好的用Python语言函数库。接下来,importnumpy。

我们可以用Pandas读入csv文件。用pd.read_csv方法从csv文件中创建DataFrame。

让我们看看数据的大小。

Test文件有80列，train文件有81列。当然，这是由于前者测试数据不包括最终的销售价格信息。

接下来，我们将使用theDataFrame.head方法查看。

我们应该在我们的比赛文件夹中为比赛提供数据字典。你也可以在这里找到。

下面是在数据描述文件中的一些字段:

·SalesPrice——以美元计算的房地产的售价。这是你要预测的目标变量。

·MSSubClass——房地产建筑类别

·MSZoning——一般的地区区分类

·LotFrontage——与房地产相距的街道英尺距离（？）

·LotArea-用平方英尺表示的大小

·Stret——所在道路的类型

·Alley-所在巷通道的类型

·LotShape-----房地产形状的分类

·LandContour-----土地的平整性

·Utilities-----可用的实用工具类型

·LotConfig——大小配置

等等。

比赛要求你预测每个家庭的最终价格。在这一点上，我们应该开始考虑我们对艾姆斯、爱荷华州(Ames,Iowa,)住房价格的了解，以及我们在这个数据集中可能看到的东西。

从这些数据，我们看到了我们所期望的特性，YrSold(房子的最后出售时间)和SalePrice。还有一些我们可能没有预料到的，比如，LandSlope(建造在屋顶上的房屋的坡度),RoofMatl(用来建造屋顶的材料)。稍后，我们将不得不决定如何处理这些字段和其他特性字段。

可以在项目的探索阶段进行一些可视化预览工作，也需要将这些功能导入到我们的环境中。预览允许我们可视化数据的分布，检查离群值，查看可能忽略的其他模式。我们将使用一个流行的可视化库Matplotlib来做这件事。

步骤二：挖掘数据，特征工程化

挑战在于预测最终的房屋出售价格，房屋出售价格信息存在SalePrice列中。通常将要预测的值称为目标变量（targetvariable）。

可以应用Series.describe获取更多的信息。

Series.describe提供了系列（series）的信息，count：系列（series）包含的行数，mean：均值,std：标准差，min：最小值，max：最大值。

在我们的数据集中房屋售价均值接近$180,000，大部分房屋价格分布在$130,000~$215,000范围内。

接下来，检查数据偏斜度（skewness），并以此来衡量数据分布的形状。

当使用回归方法时，如果目标变量出现偏斜，则有必要对目标变量进行对数变换（log-transform）。通过对数变换，可以改善数据的线性度。判断何时进行对数变换超出了本文的范围

重要的是，由最终模型生成的预测数据也要经过对数变换，这样便可将预测结果转换为原始形式。np.log对数据进行变换，np.exp进行逆变换。

plt.hist可画出SalePrice的直方图，值得注意的是在分布右侧有较长的尾部。这个分布是正向偏差。

现在应用np.log对train.SalePric进行变换，然后重新计算偏斜率，并重新画图。通过变换，斜率更靠近0值，分布更接近正态分布。

已经对目标变量进行了对数变换，现在我们再考虑下数据特征。首先，选择出数值特征并进行绘图。select_dtypes方法返回与给定数据类型匹配的列的集合Idint64

MSSubClassint64

LotFrontagefloat64

LotAreaint64

OverallQualint64

OverallCondint64

YearBuiltint64

YearRemodAddint64

MasVnrAreafloat64

BsmtFinSF1int64

BsmtFinSF2int64

BsmtUnfSFint64

TotalBsmtSFint64

1stFlrSFint64

2ndFlrSFint64

LowQualFinSFint64

GrLivAreaint64

BsmtFullBathint64

BsmtHalfBathint64

FullBathint64

HalfBathint64

BedroomAbvGrint64

KitchenAbvGrint64

TotRmsAbvGrdint64

Fireplacesint64

GarageYrBltfloat64

GarageCarsint64

GarageAreaint64

WoodDeckSFint64

OpenPorchSFint64

EnclosedPorchint64

3SsnPorchint64

ScreenPorchint64

PoolAreaint64

MiscValint64

MoSoldint64

YrSoldint64

SalePriceint64

dtype:object

DataFrame.corr方法可以给出两个列之间的相关性，我们利用该方法用于检测特征和目标变量之间的相关性。

前5个特征与SalePrice正相关，后5个特征与SalePrice负相关。为了进一步了解Over11Qual，可以利用unique方法获取唯一值（uniquevalues）。

Over11Qual数据是1~10的整数。建立一个数据透视表（pivottable，），以便深入分析Over11Qual和SalePrice之间的关系。Pandas文档（Pandasdocs）演示了如何完成分析工作。先设定“index='OverallQual'”、“values='SalePrice'”，进而查看median的值。

为了更容易地可视化这个数据透视表(pivottable)，可以使用Series.plot方法创建一个条形图。

透过上图可以看出销售价格的中位数随着数量的增加而增加。接下来应用plt.scatter生成散点图，对居住面积GrLivArea和SalePrice的关系进行可视化分析。

通过上图可以看出居住面积与价格呈现正相关关系。接下来对GarageArea进行同样的分析。

通过上图，可以看出许多房屋的GarageArea是0，表明这些房屋没有配备车库，稍后我们会利用其他特征反映这个猜测。图中显示存在一些奇异值（outliners,https://en.wikipedia.org/wiki/Outlier），这些奇异值会使预测回归线远离真实的回归线，进而影响回归模型。为保证回归模型的准确性，我们需要剔除这些奇异值。剔除野值既是技巧也是科学，有很多技术可以处理野值。

将剔除野值后的数据重新组合成新的数据帧（dataframe）。

利用剔除离群值后的数据进行画图。

处理空白值（HandlingNullValues）

在这部分，检查空白值或者缺失值。

我们将创建一个数据帧DataFrame来查看顶部的空列(topnullcolumns)。将train.isnull.sum方法链接在一起，我们返回一系列每列中的空值的计数。

该文档可以帮助我们了解缺少的值。PoolQC表示泳池质量（PoolQuality）。当PoolArea为0时，泳池质量是NaN，或者是没有泳池。

我们可以在与车库相关的列之间找到相似的关系。

看看其他列，比如对MiscFeature列进行分析，利用Series.unique方法给出唯一值列表。

我们可以使用文档来了解这些值是什么意思：

这些值描述了房屋是否有超过100平方英尺的棚子、第二个车库等设施，稍后我们会利用这些信息。在处理缺失值时，需要更多的专业知识才能做出最好的决定。

处理非数字特征（Wranglingthenon-numericFeatures）

下面我们考虑非数字特征。

Count：非空观测值的计数，

Unique：唯一值的计数

Top：出现次数最多的值，

Freq：出现次数最多值的频率。

对于这些特征，我们希望用one-hotencoding（https://quora/What-is-one-hot-encoding-and-when-is-it-used-in-data-science）法建模。one-hotencoding可以将分类数据（categoricaldata）转化为数值。通过这种转化，模型可以理解特定观测属于哪一类别。

转换和工程化特征（transformingandengineeringfeatures）

在转换特征时，需要记住在适应模型前对训练数据进行的变换，并对测试数据进行同样的变换。

我们的模型预期train数据集的特征形状与test数据集相匹配，即在train数据集工程化的特征也可以应用在test数据集。下面将以Street数据集为例来演示这个过程，Street表示到达房屋的路是Gravel还是Paved。

Street列的唯一值（uniquevalues）是Pave和Grvl，用来描述到房屋道路的类型。在训练集中，5所房屋有碎石路。我们的模型需要数字数据，因此应用one-hotencoding法将Street数据转化为布尔列。创建一个新的列enc_street，pd.get_dummies方法可以完成将数据转换成布尔量的工作。如前文所述，需要对train数据和test数据进行同样的处理。

处理之后的值可以用来训练、测试模型，这是工程化的第一个特征。特征工程化是对数据特征进行处理，使其成为适合于机器学习和建立模型。将Street特征转变为布尔列，这就是工程化了一个特征。

让我们试着对其他特征工程化。下面将对SaleCondition进行处理，就像刚刚对OverallQual所做的一样。

值得注意的是Partial的价格中位数要比其他量高，因此将Partial作为新的特征。当SaleCondition等于Patrial时，赋值1，否则赋值0。采用处理Street的方法，完成赋值操作。

下面利用图像来分析新的特征。

你可以继续对更多的特征进行处理，来改善你建立的模型。

在建模前，需要处理缺失值。我们用均值对缺失值赋值，并用data保存，这个过程称为插值。DataFrame.interpolate方法可以简单方便的实现这个过程。但是应用此方法建立的模型，可能在处理新数据时达不到最优性能。处理缺失值是建模过程中重要的一步，创造性和洞察力在这个操作中大有裨益。这是本教材中另一个可以拓展的领域。

确认data所有列值均为0。

0

步骤三：构建线性模型

让我们来完成建模前数据准备的最后步骤：分离建模的特征变量和目标变量。将特征变量分配给x，将目标变量分配给y，如上所述，使用np.log来改变模型的y变量。

使用data.drop([features],axis=1这个语句来告诉pandas想要排除的列。很明显SalePrice不会包含在变量内，同样id只是SalePrice并没有关联关系的一个索引。

y=np.log(train.SalePrice)

X=data.drop(['SalePrice''Id'],axis=1)

让我们来将数据进行划分并开始建模。使用scikit-learn里的train_test_split函数生成一份训练集和一份测试集。

使用这种方式划分数据，使我们能评估当模型面对以前从未遇到的数据时，它的执行情况会怎样。如果将测试数据都用来训练模型，那将会很难判断模型是不是过度拟合了。

train_test_split会返回4个对象

X_train：用作训练数据的一个子集

X_test：用作测试数据的一个子集--会被用来测试训练好的模型

y_train：与X_train对应的目标变量SalesPrice

y_test：与X_test对应的目标变量SalesPrice

第一个参数值X表示一组预测数据集，y是目标变量。接下来设置random_state=42，它提供了可重复的结果，因为sci-kitlearn的train_test_split函数会随机的划分数据。test_size参数告诉函数样本的占比。在这个例子中，约33%的数据被用作测试集。

开始建模

我们先创建一个线性回归模型。首先实例化模型

接下来，训练这个模型。模型拟合是一个针对不同类型模型而变化的过程。简而言之，评估预测因子和目标变量之间的关联程度，以便对新的数据集做出准确的预测。

用X_train和y_train训练模型，用X_test和y_test衡量模型，lr.fit方法会对指定的特征变量和目标变量进行线性回归训练。

model=lr.fit(X_train,y_train)

评估模型的表现并展示结果

现在，我们要评估模型的表现。每个比赛可能会对提交作品有不同的评估方式。在这场比赛中，Kaggle将使用均方根误差（RMSE）评估我们的提交作品。

我们还会考虑r平方的值，r平方值是数据与拟合回归线的接近程度的一种度量，它的取值在0到1之间，1意味着目标中的所有方差都由数据解释。一般来说，较高的r平方值意味着更好的拟合。

model.score方法默认会返回r平方值

这标志着我们的特征变量能解释目标变量中约89％的差异(todo?是否意译为准确度89%)。点击上面的链接可以了解更多。

接下来，我们来关注下RMSE。为此，我们用建立好的模型对测试数据集进行预测。

predictions=model.predict(X_test)

model.predict在训练模型之后使用这个方法，它会返回给定一组预测变量的预测列表。

函数mean_squared_error需要两个数组来计算RMSE。

fromsklearn.metricsimport

mean_squared_errorprint('RMSEis:\n'mean_squared_error(y_testpredictions))

RMSEis:0.0178417945196

RMSE这个值比r平方值更直观，它表示我们的预测值和实际值之间的距离。

我们可以用散点图来查看这种关系。

如果预测值与实际值相同，则该图将为直线y=x，因为每个预测值x将等于每个实际值y。

尝试改进模型

接下来我们尝试使用岭正则化(RidgeRegularization)来减少非关键特征变量的影响，岭正则化(RidgeRegularization)是收缩非关键特征的回归系数的过程。

现在再次实例化模型，岭正则化(RidgeRegularization)模型使用一个参数α，用来控制正则化的强度。

我们将通过循环赋于几个不同的α值来进行实验，并�...