/*construction of data set*/


data a;
do i = 1 to 90;
do j = 0+i to 90;
output;
end;
end;

run;


data a;
set a;
where i ne j;
run;



data i;
  keep sender_effect_i reciever_effect_i i age;
  mu1=0; mu2=0; var1=.75; var2=1.5; rho= .6; 
  do i = 1 to 90;
   sender_effect = mu1+sqrt(var1)*rannor(0);
    reciever_effect = (mu2+rho*(sqrt(var2)/sqrt(var1))*(sender_effect-mu1)) +
                        sqrt(var2*(1-rho**2))*rannor(0);
  age = 11+ 5*rannor(0);

  reciever_effect_i =.05*age +  reciever_effect;
  sender_effect_i = -.05*age +  sender_effect;
		output;
  end;
run;
/*
proc corr data = i cov;
var sender_effect_i reciever_effect_i;
run;
*/

/*
proc means data = i;
var sender_effect_i reciever_effect_i;
run;
*/


data j;
set i;
where i ne 1;  /*becuase 1 is never j but it doesn't matter becuase when you merge on j, there will never be a j = 1, even though there is in the data set work.j*/
j = i;
sender_effect_j = sender_effect_i;
reciever_effect_j = reciever_effect_i;
run;


data j;
set j (drop = i sender_effect_i reciever_effect_i);
run;



data a2;
merge a i;
by i;

run;


data a2;
set a2;
where j is not missing;
run;


proc sort data = a2;
by j;
run;



data a3;
merge j a2;
by j;
run;


data j_vars;
set i (keep = i age);
age_j = age;
j = i;
run;



data a3;
merge j_vars (drop = age) a3 (keep = i j sender_effect_i reciever_effect_i sender_effect_j reciever_effect_j);
by j;



run;



proc sort data = a3;
by i j;
run;




data a3b;
merge i (keep = i age) a3;
by i;
run;


data a4;
set a3b;
where sender_effect_i is not missing;

/*set value of density and reciprocity*/

age_diff = abs(age - age_j);


u = 0 + -.5*age_diff;
p = .5;

k = (1 + exp(u + sender_effect_i + reciever_effect_j) + exp(u + sender_effect_j + reciever_effect_i) 
+ exp(p + 2*u + sender_effect_i + sender_effect_j + reciever_effect_i + reciever_effect_j));


/*	p1 = Yij = 1, Yji = 0
	p2 = Yij = 0, Yji = 1
	p3 = Yij = 1, Yji = 1
	p4 = Yij = 0, Yji = 0*/


p1 = exp(u + sender_effect_i + reciever_effect_j)/k;

p2 = exp(u + sender_effect_j + reciever_effect_i)/k;

p3 = exp(p + 2*u + sender_effect_i + sender_effect_j + reciever_effect_i+ reciever_effect_j)/k;

p4 =  1/k;

sump = p1 + p2 + p3 + p4;


y10 = "1 0";
y01 = "0 1";
y11 = "1 1";
y00 = "0 0";

array m{4} y10 y01 y11 y00;
x=m{rantbl(0,p1,p2,p3,p4)};

c1 = substr(x,1,1);
c2 = substr(x,3,1);

cell1 = c1 + 0;
cell2 = c2 + 0;

age_diff2 = age_diff;

run;





proc transpose data = a4 out=trans;
var cell1 cell2;
by i j;
run;




data a5;
set trans;

if _name_ = "cell2" then ego_id = j;
if _name_ = "cell2" then alter_id = i;
if _name_ = "cell1" then ego_id = i;
if _name_ = "cell1" then alter_id = j;



run;





/*try the glimix macro*/


/*First creating reciprocity variable*/


proc sort data = a5;
by ego_id alter_id;
run;


data reciprocity;
set a4 (keep = i j cell1 cell2);
ego_id = i;

alter_id = j;

run;

proc sort data = reciprocity out = ego_recip;
by ego_id alter_id;
run;



proc sort data = a5;
by ego_id alter_id;
run;

data glimmix;
merge ego_recip (keep = ego_id alter_id cell2)
	  a5;
by ego_id alter_id;

run;

data alter_recip;
set ego_recip (drop = ego_id alter_id);

ego_id = j;
alter_id = i;
run;

proc sort data = alter_recip;
by ego_id alter_id;

run;


data glimmix2;
merge alter_recip (keep = alter_id ego_id cell1)
      glimmix;
by ego_id alter_id;

if cell1 = . then cell1 = 0;
if cell2 = . then cell2 = 0;

reciprocity = cell1 + cell2;



run;


data egovars;
set i (keep = i age);
ego_id = i;
age_ego = age;
run;

data altvars;
set i (keep = i age);
alter_id = i;
age_alter = age;

run;


proc sort data = glimmix2;
by alter_id;
run;

data glimmix3;
merge glimmix2 altvars;
by alter_id;
run;


proc sort data = glimmix3;
by ego_id alter_id;
run;


data glimmix4;
merge glimmix3 egovars;
by ego_id;
run;


data glimmix5;


set glimmix4 (keep = ego_id alter_id col1 reciprocity age_ego age_alter);

age_diff = abs(age_ego - age_alter);

friend = col1;

/*this is where I generate an ascii data file for tom*/

file "c:\sim_network.dat";

put ego_id 1-2 alter_id 4-5 friend 7 reciprocity 9 age_ego 11-12 age_alter 14-15 age_diff 17-18;

run;


/*run glimmix model*/


%include "c:/glimmix.sas";


    %glimmix(data=glimmix5, procopt=method=reml                                                                        
       covtest, stmts=%str(                                           
          class ego_id alter_id;                             
          model col1= reciprocity age_alter age_ego age_diff;                                       
          random ego_id alter_id ;                                 
       ),                                                     
       error=binomial,link=logit               
    ) ;


run;