Guest User

Untitled

a guest
Jun 29th, 2018
100
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
SAS 24.67 KB | None | 0 0
  1. %macro bartlett(data=_last_, var=, group=, bygroup=_none);                                                                              
  2.                                                                                                                                        
  3.                                                                                                                                        
  4. /*****************************************************************/                                                                    
  5. /* The Bartlett macro fills in the template for calculations on  */                                                                    
  6. /* a number of samples. The variable containing the observations */                                                                    
  7. /* is "var". The samples are identified by the classification    */                                                                    
  8. /* variable "group". The calculation may be repeated on several  */                                                                    
  9. /* groups of samples. The groups are identified by the "bygroup" */                                                                    
  10. /* variable.                                                     */                                                                    
  11. /* Furthermore, the macro calculates Bartlett's test for         */                                                                    
  12. /* homogeniety of variance.                                      */                                                                    
  13. /* (Samples for which all the observations are equal and groups  */                                                                    
  14. /* for which the number of non-missing observations is less than */                                                                    
  15. /* 2 are excluded from Bartletts test althought they appear in   */                                                                    
  16. /* the template for calculations.                                */                                                                    
  17. /*                                                               */                                                                    
  18. /* Due to the layout of the output it is recommanded that the    */                                                                    
  19. /* "group" variable either is numerical or can be identified     */                                                                    
  20. /* from the first four characters.                               */                                                                    
  21. /*                                                               */                                                                    
  22. /* The lines in the template of calculation are 80 characters    */                                                                    
  23. /* long so it is recommended to have a linesize of 80 or larger  */                                                                    
  24. /* (OPTIONS LS=80;)                                              */                                                                    
  25. /*                                                               */                                                                    
  26. /* The data set "_bartlett" contains nearly all values written   */                                                                    
  27. /* on the output.                                                */                                                                    
  28. /*                                                               */                                                                    
  29. /* Reference:                                                    */                                                                    
  30. /*         Preben Blaesild & Jørgen Granfeldt:                   */                                                                    
  31. /*         STATISTICS with APPLICATIONS IN BIOLOGY AND GEOLOGY,  */                                                                    
  32. /*         2002.                                                 */                                                                    
  33. /*                                                               */                                                                    
  34. /* Danish versions:                                              */                                                                    
  35. /* J›rn Attermann 19-01-94, 14-04-94, 08-06-94                   */                                                                    
  36. /* First english version:                                        */                                                                    
  37. /* Preben Blaesild & J›rgen Granfeldt 03-01-01                   */                                                                    
  38. /*                                                               */                                                                    
  39. /*****************************************************************/                                                                    
  40.                                                                                                                                        
  41.                                                                                                                                        
  42.                                                                                                                                        
  43.                                                                                                                                        
  44. data _mact1;                                                                                                                            
  45.   set &data;                                                                                                                            
  46.   _none=1;                                                                                                                              
  47.                                                                                                                                        
  48. proc sort;                                                                                                                              
  49.   by &bygroup &group;                                                                                                                  
  50.                                                                                                                                        
  51. proc means noprint;                                                                                                                    
  52.   var &var;                                                                                                                            
  53.   output out=_mact2 n=_ni sum=_Si USS=_USSi mean=_meani;                                                                                
  54.   by &bygroup &group;                                                                                                                  
  55.                                                                                                                                        
  56. data _bartlett;                                                                                                                        
  57.   set _mact2;                                                                                                                          
  58.   by &bygroup;                                                                                                                          
  59.   retain _n _S _USS _CT _SSD _f _s_invfi _s_flns _k . _ekskl 0;                                                                        
  60.   if first.&bygroup then do;                                                                                                            
  61.     _n=.;                                                                                                                              
  62.     _S=.;                                                                                                                              
  63.     _USS=.;                                                                                                                            
  64.     _CT=.;                                                                                                                              
  65.     _SSD=.;                                                                                                                            
  66.     _f=.;                                                                                                                              
  67.     _s_invfi=.;                                                                                                                        
  68.     _s_flns=.;                                                                                                                          
  69.     _k=.;                                                                                                                              
  70.     _ekskl = 0;                                                                                                                        
  71.   end;                                                                                                                                  
  72.   if _ni > 0 then _CTi=_Si**2/_ni;                                                                                                      
  73.   else _CTi = .;                                                                                                                        
  74.   _SSDi=_USSi-_CTi;                                                                                                                    
  75.   _fi=_ni-1;                                                                                                                            
  76.   if _fi >= 1 & _SSDi > 0 then do;                                                                                                      
  77.     _vari=_SSDi/_fi;                                                                                                                    
  78.     _s_invfi+1/_fi;                                                                                                                    
  79.     _s_flns+_fi*log(_vari);                                                                                                            
  80.     _k+1;                                                                                                                              
  81.   end;                                                                                                                                  
  82.   else do;                                                                                                                              
  83.     _ekskl + 1;                                                                                                                        
  84.     _fi = 0;                                                                                                                            
  85.     _vari = .;                                                                                                                          
  86.   end;                                                                                                                                  
  87.   _n+_ni;                                                                                                                              
  88.   _S+_Si;                                                                                                                              
  89.   _USS+_USSi;                                                                                                                          
  90.   _CT+_CTi;                                                                                                                            
  91.   _SSD+_SSDi;                                                                                                                          
  92.   _f+_fi;                                                                                                                              
  93.   if last.&bygroup then do;                                                                                                            
  94.     _var=_SSD/_f;                                                                                                                      
  95.     if _k > 1 then do;                                                                                                                  
  96.       _C=1+1/3/(_k-1)*(_s_invfi-1/_f);                                                                                                  
  97.       _2lnQ=_f*log(_var)-_s_flns;                                                                                                      
  98.       _Ba=_2lnQ/_C;                                                                                                                    
  99.       _testss=1-probchi(_Ba,_k-1);                                                                                                      
  100.     end;                                                                                                                                
  101.   end;                                                                                                                                  
  102.                                                                                                                                        
  103. /* Print calculated quantities */                                                                                                      
  104.                                                                                                                                        
  105.   file print;                                                                                                                          
  106.   if _n_=1 then do;                                                                                                                    
  107.     put / @34 'Macro "BARTLETT"';                                                                                                      
  108.     put / @37 "Data set: %upcase(&data)";                                                                                              
  109.     put @28 "Response variable: %upcase(&var)";                                                                                        
  110.     put @31 "Group variable: %upcase(&group)";                                                                                          
  111.     %if &bygroup ^= _none %then %do;                                                                                                    
  112.       put @24 "Bygroup variable: %upcase(&bygroup)";                                                                                    
  113.     %end;                                                                                                                              
  114.   end;                                                                                                                                  
  115.   if first.&bygroup then do;                                                                                                            
  116.     %if &bygroup^=_none %then %do;                                                                                                      
  117.       _bychar=trim(left(&bygroup));                                                                                                    
  118.       put // @33 "&bygroup = " @44 _bychar $28.;                                                                                        
  119.       put '            ____________________________________________________';                                                          
  120.     %end;                                                                                                                              
  121.     put // @28 'Calculations in k samples:' /;                                                                                          
  122.     put  @63 'Sample' @74 'Sample' ;                                                                                                    
  123.     put  @1 'i' @7 'ni' @14 'Si' @25 'USSi' @37 'Si2/ni' @49 'SSDi'                                                                    
  124.           @58 'fi' @62 'Variance' @75 'Mean' /;                                                                                        
  125.   end;                                                                                                                                  
  126.   _grpchar=trim(left(&group));                                                                                                          
  127.   put @1 _grpchar $4. @6 _ni 3.0 @10 _Si 9.3 @20 _USSi 12.4 @33 _CTi 12.4                                                              
  128.       @46 _SSDi 10.4 @57 _fi 3.0 @61 _vari 10.5 @72 _meani 9.4;                                                                        
  129.   if last.&bygroup then do;                                                                                                            
  130.     put '------------------------------------------------------------====================';                                            
  131.     put @6 _n 3.0 @10 _S 9.3 @20 _USS 12.4 @33 _CT 12.4                                                                                
  132.       @46 _SSD 10.4 @57 _f 3.0 @61 _var 10.5;                                                                                          
  133.     put // @33 "Bartlett's test:";                                                                                                      
  134.     put / @23 'C' @32 '-2lnQ' @44 'Ba' @52 'Pr > ChiSq' /;                                                                              
  135.     put @ 17 _C 10.5 @28 _2lnQ 10.5 @39 _Ba 10.5 @50 _testss 10.5;                                                                      
  136.     if _ekskl > 0 then do;                                                                                                              
  137.       put /// @5 'NB: ' _ekskl ' group(s) excluded from the test.';                                                                    
  138.     end;                                                                                                                                
  139.     put _PAGE_;                                                                                                                        
  140.   end;                                                                                                                                  
  141.   keep                                                                                                                                  
  142.     %if &bygroup^=_none %then %do;                                                                                                      
  143.       &bygroup                                                                                                                          
  144.     %end;                                                                                                                              
  145.     &group _ni _n _Si _S _USSi _USS _CTi _CT _SSDi _SSD _fi _f _vari _var                                                              
  146.     _meani _C _Ba _2lnQ _testss;                                                                                                        
  147. run;                                                                                                                                    
  148.                                                                                                                                        
  149.                                                                                                                                        
  150. %mend bartlett;                                                                                                                        
  151.                                                                                                                                        
  152.                                                                                                                                        
  153. DATA kraft;                                                                                                                            
  154. DO i=1 TO 4;                                                                                                                            
  155.    INPUT person@@;                                                                                                                      
  156.    INPUT vaegt@@;                                                                                                                      
  157.    DO j=1 TO 4;                                                                                                                        
  158.       INPUT kraft@@;                                                                                                                    
  159.       OUTPUT;                                                                                                                          
  160.    END;                                                                                                                                
  161. END;                                                                                                                                    
  162. DATALINES;                                                                                                                              
  163. 1 73.2 3253 3136 3245 3321                                                                                                              
  164. 2 78.3 3526 3621 3548 3612                                                                                                              
  165. 3 75.3 3325 3417 3345 3442                                                                                                              
  166. 4 72.5 3123 3112 3201 3159                                                                                                              
  167. ;                                                                                                                                      
  168. RUN;                                                                                                                                    
  169.                                                                                                                                        
  170. /* Omregning til 1000 Newton med 2 decimaler*/                                                                                          
  171.                                                                                                                                        
  172. DATA kraft;                                                                                                                            
  173. SET kraft;                                                                                                                              
  174. kraft=round(kraft/1000,0.01);                                                                                                          
  175. RUN;                                                                                                                                    
  176.                                                                                                                                        
  177. PROC GLM DATA=kraft;                                                                                                                    
  178.                                                                                                                                        
  179. %bartlett(data=kraft, var=i, group=j);                                                                                                  
  180. RUN;                                                                                                                                    
  181.                                                                                                                                        
  182.                                                                                                                                        
  183. PROC PRINT;                                                                                                                            
  184. RUN;
Add Comment
Please, Sign In to add comment