***************************************************************************; ** PROGRAM: Gange_2MRCV.sas **; ** AUTHOR: Christopher R. Bilder **; ** Department of Statistics **; ** University of Nebraska-Lincoln **; ** chris@chrisbilder.com **; ** DATE: 5-2-05 **; ** PURPOSE: Finds a vector of multinomial probabilities subject to the **; ** constraints of the Ho model. The fit of the Ho model **; ** provides configuations for the iterative proportional fitting**; ** algorithm. **; ** NOTES: **; ** 1) Do not distribute copies of this program without the permission of**; ** the author. **; ** 2) Copyright 2005 Christopher R. Bilder **; ** 3) See control_2MRCV.sas for more information **; ***************************************************************************; **********************************************************************; ** **; ** NAME: MACRO MAKE_SUBTABLE_SAME_IML **; ** PURPOSE: MAKE IML SUB-TABLES for W and Y pairs in PROC IML **; ** **; **********************************************************************; %MACRO MAKE_SUBTABLE_SAME_IML(var, numb, add); k = 0; %do i=(1+&add) %to (&numb+&add-1); %do j=(&i+1) %to (&numb+&add); m_&i._&j = repeat(0, 2, 2); m_&i._&j[1, 1] = &var._m[1+k*4 ,1]/&n; m_&i._&j[1, 2] = &var._m[2+k*4 ,1]/&n; m_&i._&j[2, 1] = &var._m[3+k*4 ,1]/&n; m_&i._&j[2, 2] = &var._m[4+k*4 ,1]/&n; k = k + 1; %end; %end; %MEND MAKE_SUBTABLE_SAME_IML; **********************************************************************; ** **; ** NAME: MACRO MAKE_SUBTABLE_SAME **; ** PURPOSE: MAKE SUB-TABLES for W and Y pairs **; ** **; **********************************************************************; %MACRO MAKE_SUBTABLE_SAME(var, numb);; %do i = 1 %to (&numb-1); %do j = (&i + 1) %to &numb; proc freq data=data2x2 noprint; tables &var.&i*&var.&j / out=out_subtable sparse; run; *Fix data set up to save; data out_subtable2; set out_subtable; var1 = "&var.&i"; var2 = "&var.&j"; if count = 0 then count = 0.5; keep var1 var2 count; run; data out_subtable_&var; set out_subtable_&var out_subtable2; run; %end; %end; %MEND MAKE_SUBTABLE_SAME; **********************************************************************; ** **; ** NAME: MACRO MAKE_SUBTABLE **; ** PURPOSE: MAKE SUB-TABLES for W and Y pairs **; ** **; **********************************************************************; %MACRO MAKE_SUBTABLE; k = 0; %do i=1 %to &nrow2; %do j=(&nrow2+1) %to (&nrow2+&ncol2); m_&i._&j = repeat(0, 2, 2); m_&i._&j[1, 1] = m[1+k*4 ,1]/&n; m_&i._&j[1, 2] = m[2+k*4 ,1]/&n; m_&i._&j[2, 1] = m[3+k*4 ,1]/&n; m_&i._&j[2, 2] = m[4+k*4 ,1]/&n; k = k + 1; %end; %end; %MEND MAKE_SUBTABLE; **********************************************************************; ** **; ** NAME: MACRO P_MAT **; ** PURPOSE: Initial multinomial probabilities **; ** **; **********************************************************************; %MACRO P_MAT; %let X = do item1 = 0 to 1%str(;); %let endit = end%str(;); *Generalize for number of items; %do j = 2 %to &ncol; %let X = &X.do item&j=0 to 1%str(;); %let endit = &endit end%str(;); %end; *Create initial matrix of all items and their corresponding joint probabilities; data set1_init; &X; p = 0.5; output; &endit; run; %MEND P_MAT; **********************************************************************; ** **; ** NAME: MACRO MARGINAL **; ** PURPOSE: Find marginal probabilities from the multinomial **; ** probabilities **; ** VARIABLES: Y1=general item marginal of interest **; ** **; **********************************************************************; %MACRO MARGINAL(Y1); do i =1 to nrow(p_matrix) ; if p_matrix[i,&Y1] = 1 then marginal[1,&Y1] = p_matrix[i,ncol(p_matrix)] + marginal[1,&Y1]; end; %MEND MARGINAL; **********************************************************************; ** **; ** NAME: MACRO OR **; ** PURPOSE: Find OR from the joint probabilities **; ** VARIABLES: Y1 and Y2 = general items marginal of interest **; ** **; **********************************************************************; %MACRO OR_CHECK(Y1, Y2); OR[ &Y1, &Y2] = p_&Y1._&Y2[1,1]*p_&Y1._&Y2[2,2]/(p_&Y1._&Y2[1,2]*p_&Y1._&Y2[2,1]); %MEND OR_CHECK; **********************************************************************; ** **; ** NAME: MACRO PAIRWISE **; ** PURPOSE: Find the pairwise probabilities from the current **; ** multinomial probabilities **; ** VARIABLES: Y1 and Y2 = general items marginal of interest **; ** **; **********************************************************************; %MACRO PAIRWISE(Y1, Y2); *Intitialize pairwise matrix; p_&Y1._&Y2 = repeat(0,2,2); do i = 1 to nrow(p_matrix); *Sum over only items with 0,0 or 0,1 or 1,0, or 1,1 and put results in a matrix; if p_matrix[i,&Y1] = 0 & p_matrix[i,&Y2] = 0 then p_&Y1._&Y2[1,1] = p_matrix[i,&ncol+1] + p_&Y1._&Y2[1,1]; if p_matrix[i,&Y1] = 0 & p_matrix[i,&Y2] = 1 then p_&Y1._&Y2[1,2] = p_matrix[i,&ncol+1] + p_&Y1._&Y2[1,2]; if p_matrix[i,&Y1] = 1 & p_matrix[i,&Y2] = 0 then p_&Y1._&Y2[2,1] = p_matrix[i,&ncol+1] + p_&Y1._&Y2[2,1]; if p_matrix[i,&Y1] = 1 & p_matrix[i,&Y2] = 1 then p_&Y1._&Y2[2,2] = p_matrix[i,&ncol+1] + p_&Y1._&Y2[2,2]; end; %MEND PAIRWISE; **********************************************************************; ** **; ** NAME: MACRO ADJUST **; ** PURPOSE: Adjust the multinomial probabilities for the current **; ** pairwise probabilities and the desired pairwise **; ** probabilities **; ** NOTES: This is the Iterative Proportional Fitting Method **; ** VARIABLES: Y1 and Y2 = general items marginal of interest **; ** **; **********************************************************************; %MACRO ADJUST(Y1, Y2); do i =1 to nrow(p_matrix); *Adjusts probabilties for variables with specified 0 or 1; if p_matrix[i, &Y1]=0 & p_matrix[i,&Y2]=0 then p_matrix[i,ncol(p_matrix)]=p_matrix[i,ncol(p_matrix)]/p_&Y1._&Y2[1,1]*m_&Y1._&Y2[1,1]; if p_matrix[i, &Y1]=0 & p_matrix[i,&Y2]=1 then p_matrix[i,ncol(p_matrix)]=p_matrix[i,ncol(p_matrix)]/p_&Y1._&Y2[1,2]*m_&Y1._&Y2[1,2]; if p_matrix[i, &Y1]=1 & p_matrix[i,&Y2]=0 then p_matrix[i,ncol(p_matrix)]=p_matrix[i,ncol(p_matrix)]/p_&Y1._&Y2[2,1]*m_&Y1._&Y2[2,1]; if p_matrix[i, &Y1]=1 & p_matrix[i,&Y2]=1 then p_matrix[i,ncol(p_matrix)]=p_matrix[i,ncol(p_matrix)]/p_&Y1._&Y2[2,2]*m_&Y1._&Y2[2,2]; end; %MEND ADJUST; **********************************************************************; ** **; ** NAME: MACRO ALL **; ** PURPOSE: Runs the above MACROS **; ** **; **********************************************************************; %MACRO ALL; proc iml; edit set1_init; read all into p_matrix; edit obstats_orig; read all var{pred} into m; *print m; %MAKE_SUBTABLE; *print m_1_3 m_1_4 m_2_3 m_2_4; edit out_subtable_W; read all var{count} into W_m; edit out_subtable_Y; read all var{count} into Y_m; *print W_m Y_m; %MAKE_SUBTABLE_SAME_IML(W, &nrow2, 0); %MAKE_SUBTABLE_SAME_IML(Y, &ncol2, &nrow2); *print m_1_2; *print m_3_4; *Set initial values; tol = &tol; save_p=1; *initial set of multinomial probabilities; iter = 1; *Loop to do IPF; do j=1 to &gange_iter while(max(abs(save_p-p_matrix[,&ncol+1]))>tol); *Used to determine convergence; save_p = p_matrix[,&ncol+1]; %do a = 1 %to (&ncol-1); %do b = (&a+1) %to &ncol; %PAIRWISE(&a, &b); %ADJUST(&a, &b); %end; %end; error = max(abs(save_p-p_matrix[,&ncol+1])); iter = iter + 1; end; p = T(p_matrix[,&ncol+1]); *print 'Final joint probabilities:' p; *print 'Iter:' iter; *Check ORs and marginals; OR = repeat(0,&ncol,&ncol); marginal = repeat(0,1,&ncol); %do a = 1 %to (&ncol-1); %do b = (&a+1) %to &ncol; %OR_CHECK(&a, &b); %end; %MARGINAL(&a); %end; %MARGINAL(&ncol); *print 'Check Marginals:' marginal; *print 'Gange program check ORs:' OR; *print 'Upper triangle part contains ORs for items: (1,2) element is OR_item1_item2'; *print ' The items are numbered as W1=item1, W2=item2, ..., Y1=item(r+1),...'; *print 'Lower triangle and diagonal will be 0''s'; *Find ending error; error = max(abs(save_p-T(p))); *print error; *print 'data set #' &sim_numb; create save_p from p_matrix; append from p_matrix; *save information about the convergence; info = &sim_numb || error || iter; col = {'sim_numb' 'error' 'iter_conv'}; create info_set from info[colname=col]; append from info; quit; %MEND ALL; **********************************************************************; ** **; ** NAME: MACRO DO_IT_GANGE **; ** PURPOSE: wrapper MACRO **; ** **; **********************************************************************; %MACRO DO_IT_GANGE; *Set initial set of multinomial probabilities; %P_MAT; *Set save data sets for W and Y subtables; data out_subtable_W; set _null_; run; data out_subtable_Y; set _null_; run; %MAKE_SUBTABLE_SAME(W, &nrow2); %MAKE_SUBTABLE_SAME(Y, &ncol2); *Perform the calculations to find the multinomial probabilities; %ALL; *Can sabe to data set outside of SAS; data set1_tau; set save_p; *file "&dir./&file..txt"; *put col1-col%eval(&ncol+1); run; *save information about the convergence; data save_info_set; set save_info_set info_set; run; %MEND DO_IT_GANGE; /* *NEED in boot_2MRCV.sas: data save_info_set; set _null_; run; *saves information about convergence; */