回帰用の一般化加法モデルの学習
この例では、最適なパラメーターで回帰用の一般化加法モデル (GAM)に学習させる方法と、学習済みモデルの予測性能を評価する方法を示します。この例では、最初に一変量の GAM に最適なパラメーター値 (線形項のパラメーター) を特定し、次に二変量の GAM の値 (交互作用項のパラメーター) を特定します。また、この例では、特定の予測に対する項のローカル効果を調べて、予測子に対する予測の部分依存を計算することで、学習済みモデルを解釈する方法についても説明します。
標本データの読み込み
標本データ セット NYCHousing2015 を読み込みます。
load NYCHousing2015データ セットには、2015 年のニューヨーク市における不動産の売上に関する情報を持つ 10 の変数が含まれます。この例では、これらの変数を使用して売価 (SALEPRICE) を解析します。
データ セットを前処理します。1000 ドル以下の SALEPRICE は、現金対価なしの所有権移転を示すものと仮定します。このような SALEPRICE をもつ標本を削除します。さらに、関数 isoutlier で識別される外れ値も削除します。その後、datetime 配列 (SALEDATE) を月番号に変換して、応答変数 (SALEPRICE) を最後の列に移動します。LANDSQUAREFEET、GROSSSQUAREFEET、および YEARBUILT に含まれるゼロを NaN に変更します。
idx1 = NYCHousing2015.SALEPRICE <= 1000; idx2 = isoutlier(NYCHousing2015.SALEPRICE); NYCHousing2015(idx1|idx2,:) = []; NYCHousing2015.SALEDATE = month(NYCHousing2015.SALEDATE); NYCHousing2015 = movevars(NYCHousing2015,'SALEPRICE','After','SALEDATE'); NYCHousing2015.LANDSQUAREFEET(NYCHousing2015.LANDSQUAREFEET == 0) = NaN; NYCHousing2015.GROSSSQUAREFEET(NYCHousing2015.GROSSSQUAREFEET == 0) = NaN; NYCHousing2015.YEARBUILT(NYCHousing2015.YEARBUILT == 0) = NaN;
table の最初の 3 行を表示します。
head(NYCHousing2015,3)
BOROUGH NEIGHBORHOOD BUILDINGCLASSCATEGORY RESIDENTIALUNITS COMMERCIALUNITS LANDSQUAREFEET GROSSSQUAREFEET YEARBUILT SALEDATE SALEPRICE
_______ ____________ ____________________________ ________________ _______________ ______________ _______________ _________ ________ _________
2 {'BATHGATE'} {'01 ONE FAMILY DWELLINGS'} 1 0 1103 1290 1910 2 3e+05
2 {'BATHGATE'} {'01 ONE FAMILY DWELLINGS'} 1 1 2500 2452 1910 7 4e+05
2 {'BATHGATE'} {'01 ONE FAMILY DWELLINGS'} 1 2 1911 4080 1931 1 5.1e+05
関数datasampleを使用して 1000 個の標本を無作為に選択し、関数cvpartitionを使用して観測値を学習セットとテスト セットに分割します。テスト用に 10% のホールドアウト標本を指定します。
rng('default') % For reproducibility NumSamples = 1e3; NYCHousing2015 = datasample(NYCHousing2015,NumSamples,'Replace',false); cv = cvpartition(size(NYCHousing2015,1),'HoldOut',0.10);
学習インデックスとテスト インデックスを抽出し、学習データ セット用とテスト データ セット用の table を作成します。
tbl_training = NYCHousing2015(training(cv),:); tbl_test = NYCHousing2015(test(cv),:);
最適なハイパーパラメーターでの GAM の学習
名前と値の引数 OptimizeHyperparameters を使用して、交差検証損失を最小化するハイパーパラメーターで GAM に学習させます。
OptimizeHyperparameters を 'auto' または 'all' として指定し、一変量パラメーターおよび二変量パラメーター両方に最適なハイパーパラメーター値を計算できます。あるいは、'auto-univariate' または 'all-univariate' オプションを使用して一変量パラメーターに最適な値を求めてから、'auto-bivariate' または 'all-bivariate' オプションを使用して二変量パラメーターに最適な値を求めることもできます。この例では、'all-univariate' および 'all-bivariate' を使用します。
一変量の GAM に学習させます。FitStandardDeviation を true として指定し、応答変数の標準偏差にもモデルを当てはめます。標準偏差モデルを標準偏差推定の精度に当てはめる場合、最適なハイパーパラメーターの使用が推奨されます。名前と値の引数 InitialLearnRateForPredictors、MaxNumSplitsPerPredictor および NumTreesPerPredictor の最適な値を fitrgam で求めるため、OptimizeHyperparameters として 'all-univariate' を指定します。再現性を得るために、'expected-improvement-plus' の獲得関数を使用します。ShowPlots を false として、Verbose を 0 として指定し、プロットの表示とメッセージの表示をそれぞれ無効にします。
Mdl_univariate = fitrgam(tbl_training,'SALEPRICE','FitStandardDeviation',true, ... 'OptimizeHyperparameters','all-univariate', ... 'HyperparameterOptimizationOptions',struct('AcquisitionFunctionName','expected-improvement-plus', ... 'ShowPlots',false,'Verbose',0))
Mdl_univariate =
RegressionGAM
PredictorNames: {'BOROUGH' 'NEIGHBORHOOD' 'BUILDINGCLASSCATEGORY' 'RESIDENTIALUNITS' 'COMMERCIALUNITS' 'LANDSQUAREFEET' 'GROSSSQUAREFEET' 'YEARBUILT' 'SALEDATE'}
ResponseName: 'SALEPRICE'
CategoricalPredictors: [2 3]
ResponseTransform: 'none'
Intercept: 5.1868e+05
IsStandardDeviationFit: 1
NumObservations: 900
HyperparameterOptimizationResults: [1×1 BayesianOptimization]
Properties, Methods
fitrgam は、最適な推定実行可能点を使用する RegressionGAM モデル オブジェクトを返します。最適な推定実行可能点は、ベイズ最適化プロセスの基となる目的関数モデルに基づいて目的関数値の信頼限界の上限を最小化するハイパーパラメーターのセットを示します。HyperparameterOptimizationResults プロパティから、または関数 bestPoint を使用して、最適な点を取得できます。
x = Mdl_univariate.HyperparameterOptimizationResults.XAtMinEstimatedObjective
x=1×3 table
InitialLearnRateForPredictors MaxNumSplitsPerPredictor NumTreesPerPredictor
_____________________________ ________________________ ____________________
0.063687 1 61
bestPoint(Mdl_univariate.HyperparameterOptimizationResults)
ans=1×3 table
InitialLearnRateForPredictors MaxNumSplitsPerPredictor NumTreesPerPredictor
_____________________________ ________________________ ____________________
0.063687 1 61
最適化プロセスの詳細については、OptimizeHyperparameters を使用した GAM の最適化を参照してください。
二変量の GAM に学習させます。名前と値の引数 Interactions、InitialLearnRateForInteractions、MaxNumSplitsPerInteraction および NumTreesPerInteraction の最適な値を fitrgam で求めるため、OptimizeHyperparameters として 'all-bivariate' を指定します。x の一変量パラメーターの値を使用して、交互作用項に最適なパラメーター値が x の値に基づいて特定されるようにします。
Mdl = fitrgam(tbl_training,'SALEPRICE','FitStandardDeviation',true, ... 'InitialLearnRateForPredictors',x.InitialLearnRateForPredictors, ... 'MaxNumSplitsPerPredictor',x.MaxNumSplitsPerPredictor, ... 'NumTreesPerPredictor',x.NumTreesPerPredictor, ... 'OptimizeHyperparameters','all-bivariate', ... 'HyperparameterOptimizationOptions',struct('AcquisitionFunctionName','expected-improvement-plus', ... 'ShowPlots',false,'Verbose',0))
Mdl =
RegressionGAM
PredictorNames: {'BOROUGH' 'NEIGHBORHOOD' 'BUILDINGCLASSCATEGORY' 'RESIDENTIALUNITS' 'COMMERCIALUNITS' 'LANDSQUAREFEET' 'GROSSSQUAREFEET' 'YEARBUILT' 'SALEDATE'}
ResponseName: 'SALEPRICE'
CategoricalPredictors: [2 3]
ResponseTransform: 'none'
Intercept: 5.1679e+05
Interactions: [3×2 double]
IsStandardDeviationFit: 1
NumObservations: 900
HyperparameterOptimizationResults: [1×1 BayesianOptimization]
Properties, Methods
最適な二変量のハイパーパラメーターを表示します。
Mdl.HyperparameterOptimizationResults.XAtMinEstimatedObjective
ans=1×4 table
Interactions InitialLearnRateForInteractions MaxNumSplitsPerInteraction NumTreesPerInteraction
____________ _______________________________ __________________________ ______________________
3 0.0010182 21 302
Mdl のモデル表示には、モデルのプロパティの一部のみが表示されます。モデルのプロパティの完全な一覧を表示するには、ワークスペースで変数名 Mdl をダブルクリックします。Mdl の変数エディターが開きます。あるいは、コマンド ウィンドウでドット表記を使用してプロパティを表示できます。たとえば、ReasonForTermination プロパティを表示します。
Mdl.ReasonForTermination
ans = struct with fields:
PredictorTrees: 'Terminated after training the requested number of trees.'
InteractionTrees: 'Terminated after training the requested number of trees.'
ReasonForTermination プロパティを使用して、各線形項および各交互作用項に対して、学習済みモデルに指定した数の木が含まれているかどうかを確認できます。
Mdl の交互作用項を表示します。
Mdl.Interactions
ans = 3×2
3 6
4 6
5 8
Interactions の各行は 1 つの交互作用項を表し、交互作用項の予測子変数の列インデックスを格納します。Interactions プロパティを使用して、モデル内の交互作用項とそれらが fitrgam でモデルに追加された順序を確認できます。
予測子名を使用して Mdl の交互作用項を表示します。
Mdl.PredictorNames(Mdl.Interactions)
ans = 3×2 cell
{'BUILDINGCLASSCATEGORY'} {'LANDSQUAREFEET'}
{'RESIDENTIALUNITS' } {'LANDSQUAREFEET'}
{'COMMERCIALUNITS' } {'YEARBUILT' }
新しい観測値での予測性能の評価
テスト標本 tbl_test とオブジェクト関数 predict および loss を使用して、学習済みモデルの性能を評価します。これらの関数をもつ完全またはコンパクトなモデルを使用できます。
学習データ セットの性能を評価するには、再代入オブジェクト関数resubPredictおよびresubLossを使用します。これらの関数を使用するには、学習データを含む完全なモデルを使用しなければなりません。
コンパクトなモデルを作成して、学習済みモデルのサイズを縮小します。
CMdl = compact(Mdl); whos('Mdl','CMdl')
Name Size Bytes Class Attributes CMdl 1x1 11975596 classreg.learning.regr.CompactRegressionGAM Mdl 1x1 12170960 RegressionGAM
線形項と交互作用項の両方を含めることにより得た結果と線形項のみを含めることにより得た結果を比較します。
テスト データ セット tbl_test について、応答を予測し、平均二乗誤差を計算します。
[yFit,ySD,yInt] = predict(CMdl,tbl_test); L = loss(CMdl,tbl_test)
L = 1.2746e+11
学習済みモデルに交互作用項を含めずに予測された応答と誤差を調べます。
[yFit_nointeraction,ySD_nointeraction,yInt__nointeraction] = predict(CMdl,tbl_test,'IncludeInteractions',false); L_nointeractions = loss(CMdl,tbl_test,'IncludeInteractions',false)
L_nointeractions = 1.2531e+11
交互作用項が含まれていない場合の方がテスト データ セットの誤差が小さくなっています。
並べ替えられた真の応答を予測応答および予測区間と一緒にプロットします。
yTrue = tbl_test.SALEPRICE; [sortedYTrue,I] = sort(yTrue); figure ax = nexttile; plot(sortedYTrue,'o') hold on plot(yFit(I)) plot(yInt(I,1),'k:') plot(yInt(I,2),'k:') legend('True responses','Predicted responses', ... '95% Prediction interval limits','Location','best') title('Linear and interaction terms') hold off nexttile plot(sortedYTrue,'o') hold on plot(yFit_nointeraction(I)) plot(yInt__nointeraction(I,1),'k:') plot(yInt__nointeraction(I,2),'k:') ylim(ax.YLim) title('Linear terms only') hold off
2 つのプロットの予測区間は同じような幅になっています。
予測の解釈
関数plotLocalEffectsを使用して、最初のテスト観測値についての予測を解釈します。また、関数plotPartialDependenceを使用して、モデル内のいくつかの重要な項の部分依存プロットを作成します。
テスト データの最初の観測値に対する応答値を予測し、予測に対する CMdl 内の項のローカル効果をプロットします。プロットに切片項を含めるには、'IncludeIntercept',true を指定します。
yFit = predict(CMdl,tbl_test(1,:))
yFit = 5.3526e+05
figure
plotLocalEffects(CMdl,tbl_test(1,:),'IncludeIntercept',true)
関数 predict で、最初の観測値 tbl_test(1,:) の売価を返します。関数 plotLocalEffects で、予測に対する CMdl 内の項のローカル効果を示す横棒グラフを作成します。ローカル効果の各値は、tbl_test(1,:) の予測売価への各項の寄与を示します。
BUILDINGCLASSCATEGORY の部分従属の値を計算し、並べ替えられた値をプロットします。学習データ セットとテスト データ セットの両方を指定して、両方のセットの部分依存の値を計算します。
[pd,x] = partialDependence(CMdl,'BUILDINGCLASSCATEGORY',[tbl_training; tbl_test]); [pd_sorted,I] = sort(pd); x_sorted = x(I); x_sorted = reordercats(x_sorted,I); figure plot(x_sorted,pd_sorted,'o:') xlabel('BUILDINGCLASSCATEGORY') ylabel('SALEPRICE') title('Patial Dependence Plot')
プロットされたラインは、学習済みモデルにおける予測子 BUILDINGCLASSCATEGORY と応答 SALEPRICE の間の平均化された部分関係性を表します。
テスト データ セットを使用して、項 RESIDENTIALUNITS と LANDSQUAREFEET の部分依存プロットを作成します。
figure plotPartialDependence(CMdl,["RESIDENTIALUNITS","LANDSQUAREFEET"],tbl_test)
x 軸 (RESIDENTIALUNITS) と y 軸 (LANDSQUAREFEET) 小目盛りは、指定されたデータ内の一意の予測子値を表します。予測子値にはいくつかの外れ値があり、RESIDENTIALUNITS と LANDSQUAREFEET の値のほとんどはそれぞれ 5 未満と 5000 未満です。プロットは、RESIDENTIALUNITS の値が 5 を超えると、SALEPRICE の値が大きくは変わらないことを示しています。
参考
fitrgam | RegressionGAM | CompactRegressionGAM | plotLocalEffects | plotPartialDependence | bayesopt | optimizableVariable
