【SAS】hpsvmプロシジャ_サポートベクターマシーン

今回はサポートベクターマシーンを試してみる。


Pythonでの実装コード


Pythonで簡単に実装できる。

#-- 線形SVM
from sklearn.svm import LinearSVC
mod1 = LinearSVC(C=1, random_state=999)
mod1.fit(X_train, y_train)
y_pred1 = mod1.predict(X_test)
 
#-- 非線形SVM
from sklearn.svm import SVC
klist=['poly','rbf','sigmoid']
for k in klist:
    mod1 = SVC(C=1, kernel=k, random_state=999)
    mod1.fit(X_train, y_train)
    y_pred1 = mod1.predict(X_test)


以前、次の記事でも触れている。
【機械学習_Python】サポートベクターマシーン - こちにぃるの日記

SASでの実装コード


ものすごくシンプルな記載。
ここにカーネルやらペナルティやらクロスバリデーションやらを設定できる。

proc hpsvm data = 訓練データ;
    target 目的変数 / order = 順序 level = 変数のタイプ;
    input 特徴量1, ... , 特徴量p;
    code file="hpsvm_fit.sas"; /* 将来データ用に予測アルゴリズム出力 */
run;
 
data _Pred_&cd.;
    set alldata;
    %inc hpsvm_fit; *-- 新しいデータで予測;
run;


SASでお試し(対象データ、データ加工)


今回もKaggleチュートリアルタイタニックデータを用いる。

データ加工は以前と同じなのでコードをたたんでいる。
【SAS】hpgenselect プロシジャ_罰則付きロジスティック回帰(LASSO) - こちにぃるの日記

★コードを見る場合ここをクリック★

/*---------------------------------------------------
    欠損値処理、変数作成
---------------------------------------------------*/
*-- 数値型:欠損値を中央値に置き換え;
proc stdize data=train out=train method=median reponly;
    var age;
run;

*-- カテゴリ型:欠損値を最頻値に置き換え;
data train;
    set train;
    Embarked = coalescec(Embarked,"S");
run;

/*---------------------------------------------------
    データ分割
---------------------------------------------------*/
proc surveyselect data=train rate=0.2 out=train outall method=srs seed=19;
run;
data train(drop=selected) test(drop=selected);
    set train;
    if selected=1 then output test; else output train;
run;

proc surveyselect data=train rate=0.2 out=train outall method=srs seed=19;
run;

proc format;
    value rollf 0="train-data" 1="validation-data" 2="test-data";
run;
data alldata;
   set train test(in=in1);
   format selected rollf.; *-- フォーマット適用;
   if in1 then selected=2; *-- テストデータに番号振る;
run;


SASでお試し(モデル作成)


線形SVMの実行。
kernel linear で線形SVM(デフォルトなので書かなくても良い)。
partition で訓練データとバリデーションデータを指定。
select cv=validateset でバリデーションデータを指定。クロスバリデーションで penalty のハイパーパラメータCを決定する。例では0.1~1の間で最適値を探してくれる。

proc hpsvm data=train;
    partition roleVar=selected(train='0' valid='1');
    performance details;
    
    target Survived / order=desc level=nominal;
    input Pclass Sex Age Fare Embarked SibSp Parch;
    
    kernel linear;
    penalty C=0.1 to 1;
    select cv=validateset seed=19;
    
    code file="hpsvm_linear_cv.sas";
run;


線形SVM(10-fold cross validation)の実行。
select cv=random fold=10 で10-fold クロスバリデーション。

proc hpsvm data=train;
    performance details;
    
    target Survived / order=desc level=nominal;
    input Pclass Sex Age Fare Embarked SibSp Parch;
    
    kernel linear;
    penalty C=0.1 to 1;
    select cv=random fold=10 seed=19;
    
    code file="hpsvm_linear_10fcv.sas";
run;


多項式カーネルでの非線形SVM(10-fold cross validation)の実行。
kernel polynom多項式カーネルを指定。
その他、rbf(Gaussカーネル)、sigmoid(シグモイドカーネル)がある。

proc hpsvm data=train;
    performance details;
    
    target Survived / order=desc level=nominal;
    input Pclass Sex Age Fare Embarked SibSp Parch;
    
    kernel polynom / degree=2;
    penalty C=0.1 to 1;
    select cv=random fold=10 seed=19;
    
    code file="hpsvm_polynom_10fcv.sas";
run;


作成モデルで予測する。

%macro _Pred(cd=);
    data _Pred_&cd.;
        set alldata;
        %inc hpsvm_&cd.;
        if P_Survived1 > 0.5 then U_Survived = 1; else U_Survived = 0;
    run;
%mend _Pred;
%_Pred(cd=linear_cv);
%_Pred(cd=linear_10fcv);
%_Pred(cd=polynom_10fcv);


SASでお試し(モデル評価)


評価指標を自作し、スコアを算出する。

%macro _Scores(inds=, tlabel=Survived, pred=U_Survived, group=selected);
    proc sql;
        create table &inds._Scores as
            select
                 %if &group.^="" %then %do; &group. , %end;
                 count(*) as N
                ,sum(case when &tlabel.=1 and &pred.=1 then 1 else 0 end) as TP
                ,sum(case when &tlabel.=1 and &pred.=0 then 1 else 0 end) as FN
                ,sum(case when &tlabel.=0 and &pred.=1 then 1 else 0 end) as FP
                ,sum(case when &tlabel.=0 and &pred.=0 then 1 else 0 end) as TN
                ,(calculated TP + calculated TN) / count(*)      as Accuracy    /*正診率(Accuracy)*/
                ,calculated TP / (calculated TP + calculated FN) as Recall      /*検出率(Recall) 、類語:感度*/
                ,calculated TN / (calculated FP + calculated TN) as Specificity /*特異度(Specificity)*/
                ,calculated TP / (calculated TP + calculated FP) as Precision   /*適合率(Precision)*/
                ,2 * calculated Recall * calculated Precision / (calculated Recall + calculated Precision) as f1_score /*F1スコア*/
            from &inds.
             %if &group.^="" %then %do; group by &group. %end;
        ;
    quit;
    
    title "&inds.";
    proc print data=&inds._Scores; run;
    title "";
%mend _Scores;
%_Scores(inds=_Pred_linear_cv);
%_Scores(inds=_Pred_linear_10fcv);
%_Scores(inds=_Pred_polynom_10fcv);


f:id:cochineal19:20210713230951p:plain

今回のいずれのモデルもハイパーパラメータ Cは 0.1 が最適値だったよう。
f1スコアは多項式カーネル0.70 くらい。

参考


SAS Help Center

本ブログは個人メモです。 本ブログの内容によって生じた損害等の一切の責任を負いかねますのでご了承ください。