Statistical programming, matrix languages, and more

writing formula

Reply
New Contributor
Posts: 4

writing formula

Hello,

I wrote a formula in SAS IML. In that formula there are two type of parameters. One of them is from a existing file and the other one is created by SAS. When I run the program it works but the results are weird. I want to be sure that my code is correct. Could you check my code please. The formula which I work is attached.

data parms;
INFILE "C:\mirt\parms.txt";
INPUT a1 a2 a3 d1;   

run;

proc iml;

  call randseed(0);

  use parms;

  read  all var {a1 a2 a3} into a;

  read all var {d1} into d;

   mean = {0,0,0};  *change these values to transfrom theta1 and theta2 to a different mean;

   cov = {1 0.45 0.45,

            0.45 1 0.45,

          0.45 0.45 1};  *change the covariance to modify the correlation between dimensions;

  theta = RANDNORMAL(5000,mean,cov);

T=theta`;

z=j(5000,30);

p=j(5000,30);

u=j(5000,30);

  do j=1 to 5000;

        do i=1 to 30;

        

p[j,i]= exp((a[i,]*T[,j])+d[i,])/(1+(exp((a[i,]*T[,j])+d[i,])));
u[j,i]=ranuni(0);

                if p[j,i]<u[j,i] then z[j,i]=0;

                else if p[j,i]>=u[j,i]then z[j,i]=1;

        end;

end;

create response from z;

append from z;

names1={theta1 theta2 theta3};

create theta from theta [colname=names1];

append from theta;

quit;

data keep2;

  set response;

  file "C:\mirt\mirt_data.txt";

  put (col1-col30)(1.0);

run;

proc corr data=keep2 alpha;

  var col1-col30;

  ods output cronbachalpha=alpha;

  run;

data alpha_need;

  set alpha;

  if _n_=1;

  keep alpha;

  file "c:\mirt\alpha\alpha.txt" mod;

  put alpha;

run;

data theta;

  set theta;

  file "C:\mirt\theta.txt";

  put theta1 theta2 theta3;

run;

Attachment
Attachment
Attachment
SAS Super FREQ
Posts: 3,222

Re: writing formula

I don't know what "the results are weird" means.  What are you seeing that is different than you expect?

Your loop appears to reconstruct the formula in the paper.  I don't see anything obviously wrong, other than the code being slightly inefficient.

I suggest that you debug and check the code on a smaller example, like 5x3 instead of 5000x30.  That will enable you to check the computations by hand to make sure they are correct.

SAS Super FREQ
Posts: 3,222

Re: writing formula

The program is much easier to understand (and more efficient) if you vectorize the computations instead of using two DO loops.  In vector form, the main computation requires two calls:

p = logistic( theta*A` + d` );

call randgen(z,"Bernoulli", p);

The computation of P is equivalent to the longer computation that you have:

do j=1 to 5000;

   do i=1 to 30;

      p[j,i]= exp((a[i,]*T[,j])+d[i,])/(1+(exp((a[i,]*T[,j])+d[i,])));

   end;

end;

New Contributor
Posts: 4

Re: writing formula

Dr. Wicklin,

Thank you SO much for your answer. I really appreciate it. I will work on both method. By the way "weird results" means different results from literature. Now, I'm pretty sure these results are  because of my design.

Post a Question
Discussion Stats
  • 3 replies
  • 289 views
  • 0 likes
  • 2 in conversation