Моделирование 10000 наблюдений из этой модели

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);

Перекрестная проверка модели линейной регрессии с использованием учеников SVM.

rng(1); % For reproducibility 
CVMdl = fitrlinear(X,Y,'CrossVal','on');

CVMdl является RegressionPartitionedLinear модель. По умолчанию программное обеспечение реализует 10-кратную перекрестную проверку. Можно изменить количество сгибов с помощью 'KFold' аргумент пары имя-значение.

Оцените среднее значение MSE тестовой выборки.

mse = kfoldLoss(CVMdl)

mse = 0.1735

Кроме того, можно получить кратные MSE, указав пару имя-значение. 'Mode','individual' в kfoldLoss.

Смоделировать данные, как в оценке k-кратной среднеквадратической ошибки.

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 

Перекрестная проверка модели линейной регрессии с использованием 10-кратной перекрестной проверки. Оптимизируйте целевую функцию с помощью SpaRSA.

CVMdl = fitrlinear(X,Y,'CrossVal','on','ObservationsIn','columns',...
    'Solver','sparsa');

CVMdl является RegressionPartitionedLinear модель. Он содержит свойство Trained, который является удержанием массива ячеек 10 на 1 RegressionLinear модели, которые программное обеспечение обучило с использованием обучающего набора.

Создайте анонимную функцию, измеряющую потери Хубера ($$\delta$$ = 1), то есть

где

$\underset{}{\overset{ejˆ}{{}_{}}}$ - остаток для наблюдения j. Пользовательские функции потери должны быть написаны в определенной форме. Правила записи пользовательской функции потери см. в разделе 'LossFun' аргумент пары имя-значение.

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);

Оцените среднюю потерю Huber по складкам. Также получите потери Huber для каждой складки.

mseAve = kfoldLoss(CVMdl,'LossFun',huberloss)

mseAve = -0.4447

mseFold = kfoldLoss(CVMdl,'LossFun',huberloss,'Mode','individual')

mseFold = 10×1

   -0.4454
   -0.4473
   -0.4452
   -0.4469
   -0.4434
   -0.4427
   -0.4471
   -0.4430
   -0.4438
   -0.4426

Чтобы определить хорошую силу лассо-штрафа для модели линейной регрессии, использующей наименьшие квадраты, реализуйте пятикратную перекрестную проверку.

Моделирование 10000 наблюдений из этой модели

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);

Создайте набор из 15 логарифмически разнесенных уровней регуляции от $$\mathrm{}{\mathrm{}}^{\mathrm{10-5}}$$ до $$\mathrm{}{\mathrm{}}^{\mathrm{10-1}}$$.

Lambda = logspace(-5,-1,15);

Выполните перекрестную проверку моделей. Чтобы увеличить скорость выполнения, транспонируйте данные предиктора и укажите, что наблюдения находятся в столбцах. Оптимизируйте целевую функцию с помощью SpaRSA.

X = X'; 
CVMdl = fitrlinear(X,Y,'ObservationsIn','columns','KFold',5,'Lambda',Lambda,...
    'Learner','leastsquares','Solver','sparsa','Regularization','lasso');

numCLModels = numel(CVMdl.Trained)

numCLModels = 5

CVMdl является RegressionPartitionedLinear модель. Поскольку fitrlinear реализует пятикратную перекрестную проверку, CVMdl содержит 5 RegressionLinear модели, которые программное обеспечение обучает на каждой складке.

Отображение первой обученной модели линейной регрессии.

Mdl1 = CVMdl.Trained{1}

Mdl1 = 
  RegressionLinear
         ResponseName: 'Y'
    ResponseTransform: 'none'
                 Beta: [1000x15 double]
                 Bias: [1x15 double]
               Lambda: [1x15 double]
              Learner: 'leastsquares'


  Properties, Methods

Mdl1 является RegressionLinear объект модели. fitrlinear построенный Mdl1 путем обучения на первых четырех складках. Поскольку Lambda - это последовательность сильных сторон регуляризации, вы можете думать о Mdl1 как 15 моделей, по одной для каждой силы регуляризации в Lambda.

Оценка перекрестно проверенного MSE.

mse = kfoldLoss(CVMdl);

Более высокие значения Lambda привести к предикторной переменной разреженности, которая является хорошим качеством регрессионной модели. Для каждой силы регуляризации следует обучить модель линейной регрессии, используя весь набор данных и те же опции, что и при перекрестной проверке моделей. Определите количество ненулевых коэффициентов на модель.

Mdl = fitrlinear(X,Y,'ObservationsIn','columns','Lambda',Lambda,...
    'Learner','leastsquares','Solver','sparsa','Regularization','lasso');
numNZCoeff = sum(Mdl.Beta~=0);

На том же рисунке постройте график перекрестно подтвержденного MSE и частоты ненулевых коэффициентов для каждой силы регуляции. Постройте график всех переменных на шкале журнала.

figure
[h,hL1,hL2] = plotyy(log10(Lambda),log10(mse),...
    log10(Lambda),log10(numNZCoeff)); 
hL1.Marker = 'o';
hL2.Marker = 'o';
ylabel(h(1),'log_{10} MSE')
ylabel(h(2),'log_{10} nonzero-coefficient frequency')
xlabel('log_{10} Lambda')
hold off

Выберите индекс силы регуляризации, который уравновешивает предикторную переменную разреженность и низкий MSE (например, Lambda(10)).

idxFinal = 10;

Извлеките модель с минимальным значением MSE.

MdlFinal = selectModels(Mdl,idxFinal)

MdlFinal = 
  RegressionLinear
         ResponseName: 'Y'
    ResponseTransform: 'none'
                 Beta: [1000x1 double]
                 Bias: -0.0050
               Lambda: 0.0037
              Learner: 'leastsquares'


  Properties, Methods

idxNZCoeff = find(MdlFinal.Beta~=0)

idxNZCoeff = 2×1

   100
   200

EstCoeff = Mdl.Beta(idxNZCoeff)

EstCoeff = 2×1

    1.0051
    1.9965

MdlFinal является RegressionLinear модель с одной прочностью регуляризации. Ненулевые коэффициенты EstCoeff близки к коэффициентам, моделирующим данные.