%macro clacah(data,id,nomvar,noeud=-1,class=5,arbre=arbre,classes=classes,chi2=yes,graph=no); /* Macro écrite par Hélène Mathian (équipe PARIS) et Patrick Brossier (GIP RECLUS) data tableau des données à classer id liste des identificateurs des observations nomvar liste des variables noeud nombre de noeud à lister class nombre de classes à générer arbre tableau résultant décrivant l'arbre classes tableau des observations avec la classe qui leur a été attribuée chi2 yes permet de recoder les données avec la métrique du chi2 no si les données ne nécessitent pas de recodage. graph yes pour obtenir les graphiques des classes no ne sort pas de graphiques Cette macro génère le graphique de l'arbre, le trait rouge indique la coupure de l'arbre avec le nombre de classes désiré. Prévoir l'instruction goptions pour récupérer ce graphique. Modif du 29/04/99 les valeurs manquantes sont mises à 0. Modif du 10/06/99 si toutes les variables sont manquantes ou nulles l'observation est rejetée, arrêt si toutes les observations d'une variable est nulle ou manquante. Modif du 18/06/99 écriture des noms de classes sur la ligne rouge. Modif du 01/07/99 contrôle de l'unicité des identificateurs. Modif du 16/01/01 changer sum=&nomvar par sum=&listevar dans proc means. (HM) Modif du 24/06/02 ajout des graphiques Modif du 29/11/02 support de liste de variables de plus de 200 caractères Modif du 21/01/06 changeùent de nom des variables de travail dans tableau SAS, ex. _i, _v Modif du 06/02/06 graphique pour chi2=no Modif du 25/04/06 déplacement du format xdf Modif du 26/04/06 Ajout du format clacah faisant la liaison entre le n° et le nom du cluster Macros utilisées: %obsnvars (attention variables globales nvars nobs) %listid (attention variables globales listevar nid) %listevar (attention variable globale listevar) %forme (attention variables globales temp typetemp) Version du 26/04/06 */ run; %global nvars nobs nid listevar; %obsnvars(&data(keep=&nomvar)); %if &nvars<=0 or &nobs<=0 %then %goto fin; %let noeudn=%eval(&nobs*2-1); %if &noeud=-1 %then %let noeud=&noeudn; %listid(&data,&id); %listevar(&data,&nomvar); %let chi2=%upcase(&chi2); %let graph=%upcase(&graph); %if &SYSSCP=MAC %then %let tempo=~'temp':; %if &SYSSCP=WIN %then %let tempo=!temp\; proc printto log=%quote("&tempo.goptions") new; proc goptions nolist pattern title; proc printto log=log; run; goptions reset=pattern; data _null_; infile %quote("&tempo.goptions") dlm=";"; file %quote("&tempo.clacah"); length ligne $ 132; input @1 ligne; if substr(ligne,1,5) in('TITLE' 'PATTE' 'FOOTN') then put ligne ";"; run; %if &chi2=YES %then %do; /* Calcul du poids des colonnes */ proc means data=&data noprint; var &nomvar; output out=tempo1b(drop=_type_ _freq_) sum=fj1-fj&nvars; /* Création du tableau des profils*/ data tempo2(keep=_ID &nomvar fi) tempo8(keep=_id &id);set &data; array nj(_i) fj1-fj&nvars; array nij(_i) &nomvar; retain ntot mult; if _N_=1 then do; set tempo1b; ntot=sum(of fj1-fj&nvars); mult=10**(int(log10(ntot))+1); /* 29/11/02 */ %do i=1 %to &nvars; _i=&i; /* A10/06/99 *4 */ if nj<=0 then do; _va=trim(%quote("&&lstv&i")); put "ERROR: la colonne " _va "est nulle"; stop; end; %end; end; ni=sum(of &nomvar); /* M10/06/99 */ if ni<=0 then do;put "WARNING: Individu " &id "est supprimé";delete;end; fi=ni*mult/ntot; _ID=&nid; do over nij; if nij=. then nij=0; nij=(nij/ni)/sqrt(nj/ntot); end; run; %let met=; %end; %else %do; data tempo2(keep=_ID &nomvar) tempo8(keep=_id &id);set &data; ni=sum(of &nomvar); array nij &nomvar; do over nij; if nij=. then nij=0; end; _ID=&nid; run; %let met=std; %end; proc sort data=tempo2;by _ID; /* A01/07/99 */ data _null_;set tempo2 end=eof;by _ID; retain _erreur_ 0; if ^(first._ID & last._ID) then do; put "ERROR: l'identificateur" _ID " n'est pas unique"; _erreur_=1; end; if eof & _erreur_ then stop; run; %obsnvars(tempo2(keep=&nomvar)); %if &nvars<=0 or &nobs<=0 %then %goto fin; %let noeudn=%eval(&nobs*2-1); %if &noeud=-1 %then %let noeud=&noeudn; proc cluster data=tempo2 method=ward &met print=&noeud outtree=&arbre; var &nomvar; id _ID; %if &chi2=YES %then freq fi;; run; data tempo1;set &arbre(where=(_sprsq_>0)); _sprsq_=round(_sprsq_*1000); proc chart data=tempo1;hbar _NCL_/freq=_sprsq_ discrete freq cfreq; run; proc tree data=&arbre noprint ncl=&class out=&classes; copy &nomvar; id _id; run; proc freq data=&classes;tables cluster*clusname/list out=tempo1; %forme; data _null_; set tempo1; file forme; if _N_=1 then put "proc format; value clacah"; put cluster +(-1) "='" cluster clusname +(-1) "'"; run; %inc forme;; run; proc sort data=&classes;by _id &nomvar; /* début construction arbre */ data tempo1;set &arbre; if _id=' '; keep _name_ _sprsq_ _NCL_; rename _name_=_parent_ _sprsq_=psprsq_ _NCL_=pncl_; label _name_='Parent of Observation or Cluster' _sprsq_='Semi-Partial R-Squared of parent' _NCL_='Number of Clusters of parent'; proc sort data=tempo1;by _parent_; proc sort data=&arbre out=tempo3;by _parent_; data tempo2;merge tempo1 tempo3;by _parent_; retain aj 0; if pncl_ =. then pncl_=0; if _ncl_=&nobs then do;_ncl_=_ncl_+aj;aj=aj+1;end; keep _name_ _parent_ PSPRSQ_ _sprsq_ _freq_ _NCL_ pncl_; run; proc sort data=tempo2;by pncl_ descending _sprsq_; data tempo3;set tempo2 end=eof;by pncl_; array _v _v1-_v&noeudn; array _d _d1-_d&noeudn; array _f _f1-_f&noeudn; array _l _l1-_l&noeudn; retain _v1-_v&noeudn; if first.pncl_ then _sens=1; else _sens=-1; _v(_ncl_)=pncl_*_sens; if eof then do;do _i=1 to &noeudn;_l(_i)=0;end; do _i=&nobs to &noeudn; _l(_i)=1; _k=abs(_v(_i)); bouc:if _k>0 then do; _l(_k)=_l(_k)+1; _k=abs(_v(_k)); goto bouc; end; end; _d(1)=1;_f(1)=_l(1); do _i=2 to &noeudn; _k=_v(_i); if _k>0 then do; _f(_i)=_f(_k); _d(_i)=_f(_i)-_l(_i)+1;end; else if _k<0 then do; _d(_i)=_d(-_k); _f(_i)=_d(_i)+_l(_i)-1;end; end; do _i=&noeudn to 1 by -1; _k=_v(_i); m=abs(_k); if _k>0 then _f(_k)=(_d(_i)+_f(_i))/2; else if _k<0 then _d(-_k)=(_d(_i)+_f(_i))/2; end; do _i=1 to &noeudn; debut=_d(_i);fin=_f(_i); _ncl_=_i; output; end; end; keep _ncl_ debut fin; run; proc sort data=tempo3;by _ncl_; proc sort data=tempo2;by _ncl_; data tempo4;merge tempo3 tempo2;by _ncl_; proc means data=tempo4 noprint;var debut _sprsq_; output out=tempo6 min=mx my max=px py; run; data tempo5; h=1; id=1; x=0;y=0;output; y=100;output; x=100*sqrt(2);output; y=0;output; x=0; output; keep id h x y; run; proc sort data=&classes(keep=cluster clusname) out=tempo1;by cluster; data tempo1;set tempo1;by cluster; if first.cluster; rename clusname=_name_; keep clusname; proc sort data=&arbre(keep=_name_ _sprsq_) out=tempo2;by _name_; proc sort data=tempo1;by _name_; data tempo1;merge tempo1(in=check) tempo2;by _name_; if check; proc means data=tempo1 noprint;var _sprsq_; output out=tempo3 max=max; proc sort data=tempo2;by descending _sprsq_; data tempo2;set tempo2; retain min; if _N_=1 then set tempo3(keep=max); if _sprsq_>max then min=_sprsq_;else do;output;stop;end; keep min max; run; data tempo7;set tempo4 end=eof; /* Modifié le 18/06/99 */ if _n_=1 then do; set tempo6; /* mx my px py */ set tempo2; /* min max */ yl=((((min+max)/2)-my)*95/(py-my))+5; c2=sqrt(2)*100; end; length function color $ 8 position $ 1; retain xsys ysys '2' color 'black' hsys '4' c2 yl; if PSPRSQ_^=. then do; function='move'; x=(((debut+fin)/2)-mx)*c2/(px-mx); y=((PSPRSQ_-my)*95/(py-my))+5; yp=y; output; function='draw'; y=((_SPRSQ_-my)*95/(py-my))+5; output; if yyl then do; function='label'; y=yl; angle=0; position='2'; size=.75; output; end; end; if _ncl_>=&nobs then do; function='label'; x=(fin-mx)*c2/(px-mx); y=0; angle=90; position='6'; size=1; output; end; else do; function='move'; x=(debut-mx)*c2/(px-mx); y=((_SPRSQ_-my)*95/(py-my))+5; output; function='draw'; x=(fin-mx)*c2/(px-mx); output; end; if eof then do; color='red'; function='move'; x=0; y=yl; output; function='draw'; x=c2; output; end; keep x y xsys ysys hsys size function color _name_ position angle; rename _name_=text; run; proc gmap data=tempo5 map=tempo5; id id; choro h/discrete nolegend annotate=tempo7 coutline=white; pattern c=black v=e; run; quit; goptions reset=pattern; /* fin construction arbre */ /* proc sort data=&classes;by _id; déjà classé */ proc sort data=tempo8;by _id; data &classes;merge tempo8 &classes;by _id; drop _id; proc sort data=&arbre;by _id; data &arbre;merge tempo8 &arbre;by _id; drop _id; proc sort data=&data(keep=&id &nomvar) out=tempo8;by &id; /* M22/06/99 */ proc sort data=&classes;by &id; data &classes;merge tempo8 &classes(in=check keep=&id cluster clusname);by &id; /* M10/06/99 */ if check; run; %if &chi2=YES %then %do; /* CALCUL DES POIDS COLONNES Fj*/ data tempo1;set tempo1b; ntot=sum(of fj1-fj&nvars); array fj fj1-fj&nvars; do over fj; fj=fj/ntot; end; run; /* CALCUL DES POIDS lignes Fi et des poids Fij */ data tempo2;set &classes(keep=&nomvar cluster clusname); if _N_=1 then set tempo1; array fij &nomvar; do over fij; if fij=. then fij=0; fij=fij/ntot; end; fi=sum(of &nomvar); run; /* Création du fichier des profils par classes */ proc sort data=tempo2;by cluster clusname; proc means data=tempo2 noprint; var fi &nomvar; by cluster clusname; output out=tempo3(drop=_type_ _freq_) sum=; run; proc means data=tempo2 noprint; var fi &nomvar; output out=tempo3b(drop=_type_ _freq_) sum=; run; /* Calcul des inerties colonnes */ data tempo3;set tempo3 tempo3b(in=check); rename fi=poids; if check then cluster=99; run; data tempo2;set tempo2; array fij &nomvar; array fj fj1-fj&nvars; do over fij; fij=(fij-fi*fj)**2/(fi*fj); end; INRi=sum(of &nomvar); run; proc means data=tempo2 noprint; var &nomvar INRi; output out=tempo4(drop=_type_ _freq_) sum=INR1-INR&nvars INRtot; run; /*concatenation des fij, fj et inr(j) et calcul des aides*/ data tempo3;set tempo3; if _N_=1 then do;set tempo1(drop=ntot); set tempo4; end; array fnj &nomvar; array fj fj1-fj&nvars; array inrj inr1-inr&nvars; array corj cor1-cor&nvars; array ctrj ctr1-ctr&nvars; do over fnj; fnj=fnj/POIDS; corj=(fnj-fj)**2/fj; ctrj=(POIDS*(fnj-fj)**2/fj)/inrj; end; RHO2=sum(of cor1-cor&nvars); INR=POIDS*RHO2/INRtot; do over corj; corj=sign(fnj-fj)*corj/rho2; end; run; %let mult=1000; %end; %else %do; Proc standard data=&classes(keep=&nomvar cluster clusname) out=tempo2 mean=0 std=1; var &nomvar; run; /* Création du fichier des classes*/ proc sort data=tempo2;by cluster clusname; proc means data=tempo2 noprint; var &nomvar; by cluster clusname; output out=tempo3(drop=_type_ _freq_) sum=&listevar n=POIDS; run; /* Calcul des aides*/ data tempo3;set tempo3; array xnj &nomvar; array corj cor1-cor&nvars; array ctrj ctr1-ctr&nvars; do over xnj; xnj=xnj/POIDS; corj=xnj**2; ctrj=POIDS*xnj**2/&nobs; end; RHO2=sum(of cor1-cor&nvars); POIDS=POIDS/&nobs; INR=POIDS*RHO2/&nvars; do over corj; corj=sign(xnj)*corj/rho2; end; run; %let mult=100; %end; proc sort data=tempo3;by cluster;run; data tempo3b;set tempo3; array xnj &nomvar; array corj cor1-cor&nvars; array ctrj ctr1-ctr&nvars; do over xnj; xnj=round(xnj*1000); corj=round(corj*1000); ctrj=round(ctrj*1000); end; POIDS=round(POIDS*1000); INR=round(INR*1000); RHO2=round(RHO2*&mult); run; /* impression modifiée le 18/06/99 pour impression plus près d'ADDAD */ goptions reset=title; %if &chi2=YES %then %do; title1 TABLEAU DES PROFILS DES CLASSES; title2 cluster=99 pour le profil moyen; title3 Toutes les valeurs sont multipliées par 1000; %end; %else %do; title1 TABLEAU DES VALEURS CENTREES-REDUITES DES CLASSES; title2 Toutes les valeurs sont multipliées par 1000,; title3 sauf RHO2 qui est multiplié par 100; %end; title4 TABLEAU DES CONTRIBUTIONS DES VARIABLES A L''ELOIGNEMENT DE LA; title5 CLASSE (COSINUS**2 SIGNE DES ANGLES CLASSES-CDG et VARIABLE); title6 TABLEAU DES CONTRIBUTIONS RELATIVES DES CLASSES; title7 A L''INERTIE SUR L''AXE; proc print data=tempo3b noobs; id cluster clusname; var POIDS INR RHO2 %do i=1 %to &nvars; %scan(&listevar,&i) cor&i ctr&i %end; ;run; goptions reset=(pattern title); %inc %quote("&tempo.clacah"); %if &graph=YES %then %do; data tempo1 tempo2;set tempo3; if cluster=99 then output tempo2; else output tempo1; data tempo2;set tempo2; array _v_(_i) &listevar; array _g_(_i) _v1-_v&nvars; do over _v_; if _v_=. then _v_=0; _g_=_v_; end; keep _v1-_v&nvars; axis2 label=none; proc format;value xdf %do i=1 %to &nvars; &i="%scan(&listevar,&i)" %end; ; %if &chi2=YES %then %do; data tempo4;set tempo2; array _g_(_i) _v1-_v&nvars; do over _g_; _moy_=_g_*100; output; end; keep _moy_ _i; goptions reset=pattern; proc gchart data=tempo4; vbar _i/ sumvar=_moy_ discrete coutline=black raxis=axis2 maxis=axis2; format _i xdf.; pattern v=s c=gray r=&nvars; run; quit; goptions reset=pattern; data tempo3;set tempo1; array _v_(_i) &listevar; array _g_(_i) _v1-_v&nvars; if _N_=1 then set tempo2; do over _v_; if _v_=. then _v_=0; _moy_=(_v_-_g_)*100; output; end; keep _moy_ _i cluster; %end; %else %do; data tempo3;set tempo1; array _v_(_i) &listevar; do over _v_; if _v_=. then _v_=0; _moy_=_v_*100; output; end; keep _moy_ _i cluster; %end; proc means data=tempo3 noprint; var _moy_; output out=tempo5(drop=_type_ _freq_) min=minx max=maxx range=ranx; data _null_; set tempo5; bax=10**floor(log10(max(abs(minx),abs(maxx))/10)); aax=floor(minx/(bax*5))*(bax*5); abx=ceil(maxx/(bax*5))*(bax*5); acx=round(ranx/10,bax*5); if acx=0 then acx=round(ranx/10,bax*.5); abx=aax+ceil((abx-aax)/acx)*acx; gg=put(aax,best.); call symput('aax',compress(gg)); gg=put(abx,best.); call symput('abx',compress(gg)); gg=put(acx,best.); call symput('acx',compress(gg)); run; axis1 order=&aax to &abx by &acx label=none; proc gchart data=tempo3;by cluster; vbar _i/ sumvar=_moy_ discrete patternid=by coutline=black raxis=axis1 maxis=axis2; format _i xdf.; %inc %quote("&tempo.clacah"); ; format cluster clacah.; run; quit; %end; %else %do; goptions reset=(title pattern); %inc %quote("&tempo.clacah"); %end; PROC DATASETS DD=WORK NOLIST NOWARN FORCE NOFS; DELETE TEMPO1 TEMPO2 TEMPO3 TEMPO4 TEMPO5 TEMPO6 TEMPO7 tempo8 %if &chi2=YES %then tempo1b TEMPO3b;; run; quit; %fin: %mend;