www.gusucode.com > stats 源码程序 matlab案例代码 > stats/SpecifyCustomRegressionLoss2Example.m

    %% Specify Custom Regression Loss
%%
% Simulate 10000 observations from this model
%
% $$y = x_{100} + 2x_{200} + e.$$
%
%
% * $X = {x_1,...,x_{1000}}$ is a 10000-by-1000 sparse matrix with 10%
% nonzero standard normal elements.
% * _e_ is random normal error with mean 0 and standard deviation
% 0.3.
%
rng(1) % For reproducibility
n = 1e4;
d = 1e3;
nz = 0.1;
X = sprandn(n,d,nz); 
Y = X(:,100) + 2*X(:,200) + 0.3*randn(n,1);
X = X'; % Put observations in columns for faster training
%%
% Train a linear regression model.  Reserve 30% of the observations as a
% holdout sample.
CVMdl = fitrlinear(X,Y,'Holdout',0.3,'ObservationsIn','columns'); 
Mdl = CVMdl.Trained{1}
%%
% |CVMdl| is a |RegressionPartitionedLinear| model. It contains
% the property |Trained|, which is a 1-by-1 cell array holding a
% |RegressionLinear| model that the software trained using the
% training set.
%%
% Extract the training and test data from the partition definition.
trainIdx = training(CVMdl.Partition);
testIdx = test(CVMdl.Partition);
%%
% Create an anonymous function that measures Huber loss ($\delta$ = 1),
% that is,
%
% $$\begin{array}{*{20}{c}}
% {L = \frac{1}{{\sum {{w_j}} }}\sum\limits_{j = 1}^n {{w_j}{\ell _j}} ,\;\;{\rm{where}}}\\
% {{\ell _j} = \left\{ {\begin{array}{*{20}{c}}
% {0.5{{\hat e}^2}}\\
% {\left| {\hat e} \right| - 0.5\;\;}
% \end{array}\begin{array}{*{20}{c}}
% {{\rm{for}}\;\;{\left| {\hat e} \right|} \le 1}\\
% {{\rm{otherwise}}}
% \end{array}} \right..}
% \end{array}$$
%
% $\hat e_j$ is the residual for observation _j_.  Custom loss functions must
% be written in a particular form.  For rules on writing a custom loss
% function, see the |LossFun| name-value pair argument.  
huberloss = @(Y,Yhat,W)sum(W.*((0.5*(abs(Y-Yhat)<=1).*(Y-Yhat).^2) + ...
    ((abs(Y-Yhat)>1).*abs(Y-Yhat)-0.5)))/sum(W);
%%
% Estimate the training- and test-sample regression loss using the
% Huber loss function.
eTrain = loss(Mdl,X(:,trainIdx),Y(trainIdx),'LossFun',huberloss,...
    'ObservationsIn','columns')
eTest = loss(Mdl,X(:,testIdx),Y(testIdx),'LossFun',huberloss,...
    'ObservationsIn','columns')