softmax regression分类相关

摘要：softmax regression是对logistic regression的二元分类问题的扩展，进行多分类。输出是每个样本相对某个类别的概率值，取概率值最大的即为对应的类别。

1、算法解释

（1）假设函数如下：

其中， 是模型的参数，
，这一项是为了归一化概率分布，使各项概率分布的和为1.

为了方便起见，将模型的参数写成k*(n+1)矩阵形式，该矩阵是将按行罗列起来得到的：
。

（2）损失函数
损失函数与logistic regression的类似，只是加入了类别：

其中，，表示1{表达式的值为真}=1。

（3）梯度函数
关于对损失函数的最小化问题，目前可采用梯度下降、L-BFGS方法，对其求导，得到梯度公式为：

即：

（4）参数更新

2、关于冗余问题的解决

在softmax regression有一个不寻常的特点，假设给减去，这时每个都变成了，于是有：

由此可见，几时对结果也丝毫没有影响。假设是损失函数的最小值，那么也会是它的最小值，这就会造成最小值不唯一。对于这个问题，可以采用增加一个权值衰减项来解决。

此时的损失函数变为：

梯度函数变为：

有了这个权重衰减项，就会使得损失函数收敛到唯一的最小值。

3、与logistics regression的关系

当类别k=2时，softmax regression就是logistic regression。softmax regression预测其中一个类别为：

另一个类别的概率为：

这与logistic regression的输出是相符的。

4、与k个二元分类器的关系

假设每个样本之间是互斥的，即每个样本只能单独的属于一个类别时，采样softmax regression进行分类。

当一个样本可以同时属于多个类别时，采用k个二元分类器。

5、gradient checking 梯度检测

使用梯度下降法解决参数更新问题时，很容易出错，所以最好对计算的损失函数进行梯度检测，以保证后续计算能较好的收敛。

梯度检测可以采用如下方法：当求出了损失函数的偏导数后，取一个参数值，计算出该参数值处的偏导数值；然后在该参数值附近取2个很近很近的点，利用损失函数在这个两个点值的差除以这2个点的距离，比较这两次计算出的结果是否相等，如果接近相等的话，则说明很大程度上，偏导数没有计算出错，这时候一定要记住不要再运行gradient checking，这样很耗时。

6、代码要点

（1）groundTruth=full(sparse(label,1:samNum,1));

生成一个10*60000的矩阵，表示每个样本对应的类别，并使得对应类别的值为1（通过sparse实现，然后将sparse用full变为全矩阵。）

7、代码实现

在本次试验中，没有进行特征提取，使用原始像素特征。

clear all;clc;addpath ..\kmeans_demo\minFunc%% load data% 使用minst数据，由于存储方式是ubyte,所以使用如下读取方式% 训练数据包括60000个28*28的图及其标签% 测试数据包括10000个28*28的图及其标签f=fopen('E:\MATLAB\Workspace\data\train-images.idx3-ubyte','r');train_data=fread(f,inf,'uint8');train_data=double(reshape(train_data(17:end),28*28,60000));train_data=[ones(size(train_data,2),1)';train_data];%给训练数据加截距fclose(f);f=fopen('E:\MATLAB\Workspace\data\train-labels.idx1-ubyte','r');train_label=fread(f,inf,'uint8');train_label=double(train_label(9:end))+1;fclose(f);f=fopen('E:\MATLAB\Workspace\data\t10k-images.idx3-ubyte','r');test_data=fread(f,inf,'uint8');test_data=double(reshape(test_data(17:end),28*28,10000));fclose(f);f=fopen('E:\MATLAB\Workspace\data\t10k-labels.idx1-ubyte','r');test_label=fread(f,inf,'uint8');test_label=double(test_label(9:end))+1;fclose(f);%% initilize these parameter% 输入的特征长度 feaSize ，在本次试验中，没有进行特征提取，使用原始像素特征% 参数向量 theta，随机初始化，长度是类别数和特征数的乘积,本人理解为每个特征在每个类别占的权重% 类别数 ClassNum% 初始化权重衰减因子 lambda% debug：用于梯度检测feaSize=28*28+1;ClassNum=10;lambda=1e-4;debug=false;if debug    feaSize = 8;    train_data = randn(8, 100);    train_label = randi(10, 100, 1);endtheta=0.005*randn(feaSize*ClassNum,1);%% compute cost function% 没有提取特征，直接进行损失函数的计算% 损失函数的参数应包括：权重衰减因子lambda、参数theta、样本种类 ClassNum% 训练样本train_data,样本标签train_label、特征数feaSize% grad:梯度函数，每个特征关于每个类别的梯度[cost,grad]=softmaxCost(theta,lambda,ClassNum,feaSize,train_data,train_label);disp(cost);%% gradient chccking% 梯度检测,检测参数的正确性，方便后面的模型训练if debug    numGrad=computeNumericalGradient(@(x) softmaxCost(x,...        lambda,ClassNum,feaSize,train_data,train_label),theta);    fprintf('the numGrad and grad is:\n');    disp([numGrad grad]);    diff=norm(numGrad-grad)/norm(numGrad+grad);   fprintf('the different is:%f\n',diff); end%% train the model % 训练分类模型，使用minFunc优化，参数包括：% 学习到的最优参数以及输入的特征个个数尺寸以及类别信息% opts：可选参数，包括max_iter等等% model的输出包括：if ~exist('opts','var')    opts=struct;endmodel=softmaxTrain(feaSize,lambda,ClassNum,train_data,train_label,opts);disp(model);%% prediction% 对训练的模型进行检测其正确率[pred]=softmaxPred(model,train_data);test_acc=mean(train_label(:)==pred(:));%若相等，test_label(:)==pred(:);的值就是1，不等则为0,fprintf('the accuacy is %f%%\n',100*test_acc);figure(1); hold on;  plot(find(pred(:)==1),'-'); plot(find(pred(:)==2),'--'); plot(find(pred(:)==3),'-'); plot(find(pred(:)==4),'--'); plot(find(pred(:)==5),'-'); plot(find(pred(:)==6),'--'); plot(find(pred(:)==7),'-'); plot(find(pred(:)==8),'--'); plot(find(pred(:)==9),''); plot(find(pred(:)==10),'--'); legend('1','2','3','4','5','6','7','8','9','10'); xlabel('class num'); ylabel('all num');figure(2); hold on;  plot(find(train_label(:)==1),'-'); plot(find(train_label(:)==2),'--'); plot(find(train_label(:)==3),'-'); plot(find(train_label(:)==4),'--'); plot(find(train_label(:)==5),'-'); plot(find(train_label(:)==6),'--'); plot(find(train_label(:)==7),'-'); plot(find(train_label(:)==8),'--'); plot(find(train_label(:)==9),'-'); plot(find(train_label(:)==10),'--'); legend('1','2','3','4','5','6','7','8','9','10'); xlabel('class num'); ylabel('all num');

softmaxCost.m

function [cost,grad]=softmaxCost(theta,lambda,ClassNum,feaSize,data,label)% groundTruth:补全的稀疏矩阵，内容是每个样本对应的类别值为1% z是theta*x，后对z进行归一化处理，在输入到exp函数中% cost：损失函数% p:每个样本对对应类别的概率值% grad：梯度函数,后面必须转成行向量，因为在minFunc中计算用行向量cost=0;samNum=size(data,2);theta=reshape(theta,ClassNum,feaSize);thetaGard=zeros(ClassNum,feaSize);groundTruth=full(sparse(label,1:samNum,1));z=theta*data;m=bsxfun(@minus,z,max(z,[],1));m=exp(m);p=bsxfun(@rdivide,m,sum(m));cost=-(1/samNum)*groundTruth(:)'*log(p(:))+(lambda/2)*sum((theta(:).^2));thetaGard=-(1/samNum)*(groundTruth-p)*data'+lambda*theta;grad=[thetaGard(:)];

softmaxTrain.m

function [model]=softmaxTrain(feaSize,lambda,ClassNum,data,label,opts)% 设置opts的各种取值if ~isfield(opts,'max_iter')    opts.max_iter=100endtheta=0.005*randn(ClassNum*feaSize,1);opts.method='lbfgs';minFuncOptions.display='on';[optsTheta,cost]=minFunc(@(p) softmaxCost(p,lambda,ClassNum,feaSize,data,label),theta,opts);model.optsTheta=reshape(optsTheta,ClassNum,feaSize);model.feaSize=feaSize;model.ClassNum=ClassNum;model.cost=[cost];

softmaxPred.m

function [pred]=softmaxPred(model,data)theta=model.optsTheta;pred=zeros(1,size(data,2));[nop,pred]=max(theta*data);%nop最大值，pred最大值的行标，也就是种类数

computeNumericalGraident.m

function numGrad=computeNumericalGradient(J,theta)numGrad=zeros(size(theta));epsilon=1e-4;n=size(theta,1)E=eye(n);for i=1:n    delta=E(:,i)*epsilon;    numGrad(i)=(J(theta+delta)-J(theta-delta))/(epsilon*2);  end

实验结果：使用训练集数据测试结果为：94.09%，使用测试集数据结果为92.12%。

8、参考文献

http://www.cnblogs.com/tornadomeet/archive/2013/03/23/2977621.html

http://www.cnblogs.com/tornadomeet/archive/2013/03/22/2975978.html

http://deeplearning.stanford.edu/wiki/index.php/Softmax%E5%9B%9E%E5%BD%92

http://www.cnblogs.com/tornadomeet/archive/2013/03/20/2970724.html

http://deeplearning.stanford.edu/wiki/index.php/Exercise:Sparse_Autoencoder#Step_3:_Gradient_checking

ps：有问题请及时指出，共同学习，共同进步！