设为首页收藏本站

EPS数据狗论坛

 找回密码
 立即注册

QQ登录

只需一步,快速开始

查看: 1495|回复: 0

SAS逻辑回归之二分类

[复制链接]

40

主题

378

金钱

580

积分

初级用户

发表于 2019-8-19 14:46:33 | 显示全部楼层 |阅读模式

数据集这里用的是australian,有14个自变量Xi,一个因变量Y,Y值只取0或1。

代码如下:
  1. /*逻辑回归数据集australian(690个观测值,每个含14个属性,目标变量y(0、1))*/  
  2.     /*导入数据集australian到逻辑库work中*/  
  3.     proc import out=aus  
  4.         datafile="\\vmware-host\Shared Folders\桌面\SAS\\data\australian.csv"      /*文件路径*/  
  5.         dbms=csv replace;                               /*文件类型指定*/  
  6.         delimiter=',';  
  7.         getnames=yes;                                   /*是否将第一列作为列名*/  
  8.     run;  
  9.       
  10.     /*查看数据集*/  
  11.     proc print data=aus;  
  12.     run;  
  13.       
  14.     /****************************  使用交叉验证法选择最优模型  *****************************************/  
  15.       
  16.     /*利用10-折交叉验证法计算测试集上的预测准确率*/  
  17.     %let k=10;                            /*定义宏变量-交叉验证的折数k*/  
  18.     %let rate=%sysevalf((&k-1)/&k);       /*给出交叉验证的样本抽样比率(因为宏变量k的本质是文本,不能直接参与运算,要将其视为数字计算要用%evalf or %sysevalf)*/  
  19.       
  20.     /*生成交叉验证的10个样例,保存在cv中*/  
  21.     proc surveyselect data=aus   
  22.                   out=cv            /*生成的样例全部放在数据集cv中*/  
  23.                   seed=158  
  24.                   samprate=&rate    /*抽样比率设定,宏变量rate的调用要加&*/  
  25.                   outall            /*输出全部数据*/  
  26.                   reps=10;          /*指定样本重复的次数*/  
  27.     run;  
  28.       
  29.     /*交叉验证的生成数据集中,selected列为1表示该行为训练集样本,0表示测试集样本,这里为new_y赋值,  
  30.       若selected=1,则可获得Y的值,若为0,该行的new_y为空。接下来给出new_y为空行的预测值。*/  
  31.     data cv;  
  32.       set cv;  
  33.        if selected then new_y=Y;  
  34.       run;  
  35.       
  36.     /*逻辑回归主程序 - 10折交叉验证*/  
  37.     ods output parameterestimates=paramest   /*输出交叉验证的参数估计值*/  
  38.                association=assoc;            /*输出交叉验证的C统计量*/  
  39.     proc logistic data=cv des;            /*des控制以Y=1来建模*/  
  40.         /* class new_y (param=ref ref='yes');  若new_y是分类变量,则用class对其参数化处理,这里选择处理方式为ref,以“yes”作为参考水平,以便于后续odds的计算*/  
  41.         model new_y=X1-X14 / SELECTION=STEPWISE SLE=0.1 SLS=0.1;  
  42.         by replicate;                         /*以交叉验证的组别来分组建模*/  
  43.         output out=out1(where=(new_y=.))    /*只给出测试集的预测结果(即new_y为空的样本)*/  
  44.                p=y_hat;  
  45.     run;  
  46.     ods output close;  
  47.       
  48.     data out1;   
  49.         set out1;  
  50.         if y_hat>0.5 then pred=_LEVEL_ ;     /* PHAT为logistic方程针对每个观察体计算的属于该组别的概率,若PHAT>0.5,则属于该组别(这里level为1),否则,属于另一组别 */  
  51.         else pred=0;                     /* 本例为二分类,概率依照level(1)计算,因此另一类为0 */  
  52.     run;  
  53.       
  54.     /*汇总交叉验证的结果*/  
  55.     /*计算预测准确率(测试集中预测准确的样本占预测总样本的概率)*/  
  56.     data out2;  
  57.         set out1;  
  58.         if Y=pred then d=1;  /*d为真实值和预测值的误差,这里设无误差为1,有误差为0*/  
  59.         else d=0;  
  60.     run;  
  61.       
  62.     proc summary data=out2;  
  63.         var d;  
  64.         by replicate;  
  65.         output out=out3 sum(d)=d1;   /*预测正确的个数*/  
  66.     run;  
  67.       
  68.     data out3;  
  69.         set out3;  
  70.         acc=d1/_freq_;   /*预测准确率*/  
  71.         keep replicate acc;  
  72.     run;  
  73.       
  74.     /*结果中加入交叉验证的C统计量(度量观测值和预测值之间的一致性,越大越好)*/  
  75.     data assoc;  
  76.         set assoc;  
  77.         where label2="c";  
  78.         keep replicate cvalue2;  
  79.     run;  
  80.       
  81.     /*合并交叉验证的统计结果*/  
  82.     data cvresult;  
  83.     merge assoc(in=ina) out3(in=inb);  
  84.     keep replicate cvalue2 acc;  
  85.     run;  
  86.       
  87.     proc print data=cvresult;  
  88.     title'交叉验证组号、c统计量、预测准确率';  
  89.     run;  
  90.       
  91.     title '交叉验证最优模型选择:组号、预测准确率';  
  92.     ods output SQL_Results=cvparam;      /*保存最优模型结果在cvparam数据集中*/  
  93.     proc sql ;  
  94.         select replicate,acc from cvresult having acc=max(acc);  
  95.     quit;  
  96.     ods output close;  
  97.       
  98.       
  99.       
  100.     /***************** 以交叉验证的最优结果组进行建模  *************************************/  
  101.     /*以最优组合从cv的10个样例中拿出最优样例,作为训练集和测试集*/  
  102.     /*取出最优组号对应的selected=1的行,作为训练集train,其余的作为测试集test*/  
  103.     proc sql ;  
  104.         create table train as  
  105.         select * from cv where replicate in (select replicate from cvparam)  
  106.         having selected=1;  
  107.         create table test as  
  108.         select * from cv where replicate in (select replicate from cvparam)  
  109.         having selected=0;  
  110.     run;  
  111.       
  112.     TITLE '--------Logistic Regression - 数据集Neur - 建模方法 STEPWISE ---------------------------';  
  113.       
  114.     /* 逻辑回归主程序 - 通过训练集建立logistic模型*/  
  115.     proc logistic data=train  DES                    /*根据分类值从大到小选择建模组别,此处为yes*/  
  116.                         covout outest=Nout_step  /*输出建模参数估计值及变量间的协方差矩阵*/  
  117.                         outmodel=model            /*输出建模结果(若想要通过已有的建模结果来预测新数据集,这里可以用inmodel实现)*/  
  118.                         simple;                          /*输出变量的简单统计量*/   
  119.             /* class Y (param=ref ref='yes');  若Y是分类变量,则用class对其参数化处理,这里选择处理方式为ref,以“yes”作为参考水平,以便于后续odds的计算*/  
  120.             MODEL Y=X1-X14                             /*logistic回归模型:反应变量=自变量1 2 3...*/  
  121.                           / SELECTION=STEPWISE           /*选择建模方式 - 逐步排除法*/   
  122.                             SLE=0.1 SLS=0.1              /*变量在模型中的显著程度,默认为0.05*/   
  123.                             details                      /*输出模型界定的过程,包括自变量的检定和相关系数的值*/  
  124.                             lackfit                      /*输出HL拟合优度*/  
  125.                             RSQ                          /*模型解释度R方*/  
  126.                             STB                          /*输出标准化模型后的参数*/  
  127.                             CL                           /*参数估计和置信区间*/  
  128.                             itprint                      /*输出分析每个步骤的统计量*/  
  129.                             corrb                        /*输出变量的相关矩阵*/  
  130.                             covb                         /*输出变量的协方差矩阵*/  
  131.                             ctable                       /*输出不同阈值下的二分类变量的分组情况,类似于ROC曲线上的每个点的值*/  
  132.                             influence                    /*输出观察体中每个变量统计量,便于找出对分析结果影响力较大的观察体*/  
  133.                             IPLOTS ;                     /*针对influence的结果画出图形,影响力过高的观察体在图形上都会显得特别突出*/  
  134.      score data=train outroc=train_roc;            /*通过score语句得到训练集上一系列的sensitivity和specificity,画出ROC曲线*/  
  135.      score data=test   
  136.            out=test_pred   
  137.            outroc=test_roc;                      /*通过score来预测测试集,结果保存在test_pred中,画出ROC曲线*/  
  138.     OUTPUT out=train_pred                        /*保存模型预测结果在该数据集中,数据集中包含的列由以下添加的统计量给出*/  
  139.                 P=PHAT  lower=LCL upper=UCL              /*输出文件中包含每个观察体属于logistic方程预测组别的概率,用PHAT作列名,LCL和UCL为置信上下限的值*/  
  140.                 RESCHI=RESCHI  RESDEV=RESDEV             /*Pearson残差和偏差残差,找出与模型不太符合的观察体*/  
  141.                 DIFCHISQ=DIFCHISQ  DIFDEV=DIFDEV         /*检测观察体对对皮尔森卡方适合度和对偏激统计量的影响程度,越大说明与模型越不符*/  
  142.                                                          /* 还可加入的统计量:C、CBAR、DFBETAS、H、XBETA、STDXBETA */  
  143.                / ALPHA=0.1;                              /*界定P值的信赖度,默认为0.05,对应信赖度为95%,这里为90%*/  
  144.     run;      
  145.     quit;  
  146.       
  147.     /*   
  148.     逻辑回归主程序 - 根据logistic模型对测试集进行预测(有需要时可使用独立的logistic过程对新数据进行预测)  
  149.     proc logistic inmodel=model;                 
  150.         SCORE data=test                           
  151.               outroc=predict_roc;               
  152.     run;        
  153.     */  
  154.       
  155.     /* 训练集的预测结果中只给出了预测概率,接下来根据0.5分界将观察体归到具体的类中,加一列“pred”(预测组别)*/  
  156.     data train_pred;   
  157.         set train_pred;  
  158.         if PHAT>0.5 then pred=_LEVEL_ ;     /* PHAT为logistic方程针对每个观察体计算的属于该组别的概率,若PHAT>0.5,则属于该组别(这里level为1),否则,属于另一组别 */  
  159.         else pred=0;                     
  160.     run;  
  161.       
  162.     /* 输出混淆矩阵 - 训练集*/  
  163.     ods output CrossTabFreqs=ct_train;   /*保存混淆矩阵表(训练集)*/  
  164.     ods trace on;  
  165.     proc freq data=train_pred;  
  166.         tables Y*pred;  
  167.     run;  
  168.     ods trace off;  
  169.     ods output close;  
  170.       
  171.     proc sql;  
  172.         create table acc1 as  
  173.         select sum(percent) from ct_train where (Y=pred and Y ^=.);  
  174.     proc print data=acc1;  
  175.     title '训练集上的预测准确率';  
  176.     run;  
  177.       
  178.       
  179.     /* 输出混淆矩阵及准确率等指标 - 测试集*/  
  180.     ods output CrossTabFreqs=ct_test; /*保存混淆矩阵表(测试集)*/  
  181.     proc freq data=test_pred;  
  182.         tables F_Y*I_Y ;  
  183.     run;  
  184.     ods output close;  
  185.       
  186.     proc sql;  
  187.         create table acc2 as  
  188.         select sum(percent) from ct_test where (F_Y=I_Y and F_Y ^='');  
  189.     proc print data=acc2;  
  190.     title '测试集上的预测准确率';  
  191.     run;
复制代码
您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

关闭

站长推荐上一条 /1 下一条

客服中心
关闭
在线时间:
周一~周五
8:30-17:30
QQ群:
653541906
联系电话:
010-85786021-8017
在线咨询
客服中心

意见反馈|网站地图|手机版|小黑屋|EPS数据狗论坛 ( 京ICP备09019565号-3 )   

Powered by BFIT! X3.4

© 2008-2028 BFIT Inc.

快速回复 返回顶部 返回列表