# SAS 相關筆記

### proc MI 多重插補法(補遺漏值)
1.透過直方圖，先觀察補值的變數狀況

proc univariate data= test ;
 histogram v2_creatine / cbarline=grey cfill=ligr;
     inset n nmiss mean std min Q1 median Q3 max /  
     header='Descriptive Statistics' pos=ne noframe;
   histogram v2_al / cbarline=grey cfill=ligr;
     inset n nmiss mean std min Q1 median Q3 max /  
     header='Descriptive Statistics' pos=ne noframe;
run;

2.尋找每個變數最適合的Box-cox轉換，將其轉換為常態分佈

*於資料中增加const變數，並設定為常數1;
data test.lab_base_data_a3_2; set test.lab_base_data_a3;
   const=1;
run;

3.利用proc transreg來選擇線性回歸中變數的最佳Box-cox轉換
如2.所新增的常數變數，在線性回歸時會自動去掉，線性回歸相當於對Box-cox轉換後去擬合的常態分佈

*以v2_esr為應變數，const為自變數做線性回歸;
*Box-cox變數之參數從-3~3，間隔為0.1中挑選;
proc transreg data=test.lab_base_data_a3_2;
   model BoxCox(v2_esr / lambda=-3 to 3 by .1) = identity(const);
run;   

從下圖中可觀察到最適合

4.將最大/最小值存入巨集中

proc sql noprint;
  select min(v2_esr), max(v2_esr)
    into :minv2_esr, :maxv2_esr 
  from test.lab_base_data_a3;
quit;

5.開始進行多重插補法Proc MI
輸出的資料檔即為以插補完的數據

/*nimpute為指定插補次數，
將4.存入巨集中的變數在這裡呼叫使用*/
proc mi
  data=test.lab_base_data_a3 out=test.mi_lab nimpute=5 
  minimum=&minv2_esr. maximum=&maxv2_esr.;
  *lambda的參數為3.中觀察到的參數
  Transform boxcox(v2_esr/lambda=0.2);
*放入要補值的變數以及相關變數;
  var v2_esr;
run;

6.透過proc univariate來觀察插補完的數據狀況

/*_Imputation_代表插補的順序編號
因一開始設定nimpute=5，所以會產生五組插補結果，
觀察這五組的平均值及標準誤差*/
proc univariate data=test.mi_lab;
  var v2_esr ;  by _Imputation_;
  output out=test.mi_lab_results mean=mean stderr=std;
  run;

7.觀察變數各次插補的分析結果

/*modeleffects參數估計的狀況
stderr 參數估計的標準誤差*/
proc mianalyze data=test.mi_lab_results;
  modeleffects mean;
  stderr std;
run;

### Trajectory 軌跡分析
code 參考網址:https://www.andrew.cmu.edu/user/bjones/cnorm.htm
https://zhuanlan.zhihu.com/p/653507660
1.整理資料格式，將long form轉換成wide form
long form格式為

| ID | date | value |
| -------- | -------- | -------- |
| aaa | 2003-05-09| 100 |
| aaa | 2005-09-09| 80 |
| aba | 2003-03-20| 75 |
| aba | 2004-05-19| 120 |
| aba | 2005-10-06| 90 |

wide form須為一人一筆，可透過轉置來調整

proc transpose data=case.lab_info out=case.lab_info_a;
var median_eos;id lab_yr;by idcode;run;

資料格式會轉換成
| ID | value_1 | value_2 | value_3 | t_1 | t_2 | t_3 |
| -------- | -------- | -------- | -------- |-------- | -------- | -------- |
| aaa | 100| 80 | |1|3|5|
| aba | 75 | 120 | 90 |1 |3| 5|

2.整理完格式，直接執行官方的四個macro sas code，就接續分析

/*var 要分析的變數，indep 為時間變數，要設定變數的最大值MAX，ngroups 要分組的組數
  order 選取分組要擬合的類型，線性(1)、二元quadratic(2)、立方cubic(3)*/
PROC TRAJ DATA=OPPOSITN OUTPLOT=OP OUTSTAT=OS OUT=OF OUTEST=OE ITDETAIL;
    ID ID; VAR value_1-value_3; INDEP t_1-t_3;
    MODEL CNORM; MAX 10; NGROUPS 3; ORDER 3 3 3;
    /*model -> continuous normal model（連續常態模型), logit (類別變項)，zip ， BETA */
    RUN;
%TRAJPLOT(OP,OS,'Opposition vs. Scaled Age','Cnorm Model','Opposition','Scaled Age')
/*將擬合的組別成圖，也可觀察組別比例*/


### Forest plot
用HR當範例

/*匯出HR的結果*/
ods output ParameterEstimates=PE;/*ods output OddsRatios=orci;*/
proc phreg data=case.ana_data2;
   class group(ref='0') sex(ref='M') / param=ref;
   model follow_month*death(0) = group sex age;
run;
/*繪圖*/
proc sgplot data=PE;
   scatter y=Parameter x=HazardRatio / xerrorlower=HRLowerCL xerrorupper=HRUpperCL markerattrs=(symbol=CircleFilled);
   refline 1 / axis=x lineattrs=(pattern=shortdash);
   xaxis type=log label="Hazard Ratio (95% CI)";
   yaxis discreteorder=data;
run;