BookmarkSubscribeRSS Feed
thanoon
Calcite | Level 5

hi all

i want to simulate mv ordinal variables in sas/iml and i have errors in this code is there any one help me please .

proc iml;

load module=_all_;                     /* load the modules */

    /* P1   P2    P3   p4   p5   p6   p7   p8   p9   p10 */

P = {0.25  0.50  0.20  0.10  0.20  0.30  0.25  0.15  0.20  0.25,

     0.75  0.20  0.15  0.25  0.15  0.10  0.15  0.25  0.35  0.45,

      .    0.30  0.25  0.15  0.25  0.30  0.35  0.65  0.55  0.15,

      .     .    0.40  0.25  0.35  0.10  0.25  0.10  0.15  0.20,

      .     .     .    0.15  0.25  0.20  0.30  0.15  0.35  0.45,   

      .     .     .     .    0.45  0.25  0.55  0.15  0.65  0.25,

      .     .     .     .     .    0.35  0.25  0.25  0.75  0.15,

      .     .     .     .     .     .    0.15  0.50  0.55  0.35,

      .     .     .     .     .     .     .    0.25  0.15  0.45,

      .     .     .     .     .     .     .    0.30  0.25  0.25};

/* expected values and variance for each ordinal variable */

Expected = OrdMean(P) // OrdVar(P);

varNames = "X1":"X10";

print Expected[r={"Mean" "Var"} c=varNames];

/* test the RandMVOrd function */

Delta = {1.0  0.4  0.3  0.5  0.3  0.6  0.4  0.1  0.5  0.2,

         0.4  1.0  0.4  0.2  0.4  0.6  0.8  0.2  0.4  0.1,

         0.3  0.4  1.0  0.5  0.4  0.6  0.1  0.3  0.2  0.4,

         0.5  0.2  0.5  1.0  0.2  0.1  0.4  0.5  0.6  0.5,

  0.3  0.4  0.4  0.2  1.0  0.3  0.2  0.4  0.1  0.6,

  0.6  0.6  0.6  0.1  0.3  1.0  0.5  0.6  0.7  0.3,

  0.4  0.8  0.1  0.4  0.2  0.5  1.0  0.2  0.8  0.1,

  0.1  0.2  0.3  0.5  0.4  0.6  0.2  1.0  0.2  0.4,

  0.5  0.4  0.2  0.6  0.1  0.7  0.8  0.2  1.0  0.5,

  0.2  0.1  0.4  0.5  0.6  0.3  0.1  0.4  0.5  1.0};

call randseed(54321);

X = RandMVOrdinal(1000, P, Delta);

print first[label="First 1000 Obs: Multivariate Ordinal"];

23 REPLIES 23
thanoon
Calcite | Level 5

after some correction and still have some errors .

proc iml;

load module=_all_;                     /* load the modules */

    /* P1   P2    P3    p4    p5    p6    p7    p8    p9   p10 */

P = {0.25  0.50  0.20  0.20  0.20  0.30  0.10  0.10  0.20  0.10,

     0.75  0.20  0.15  0.25  0.15  0.10  0.15  0.10  0.10  0.10,

      .    0.30  0.25  0.15  0.25  0.10  0.15  0.10  0.10  0.10,

      .     .    0.40  0.25  0.20  0.10  0.10  0.10  0.15  0.10,

      .     .     .    0.15  0.10  0.20  0.10  0.15  0.10  0.10,   

      .     .     .     .    0.10  0.10  0.10  0.15  0.10  0.10,

      .     .     .     .     .    0.10  0.15  0.10  0.15  0.10,

      .     .     .     .     .     .    0.15  0.10  0.10  0.10,

      .     .     .     .     .     .     .    0.10  0.15  0.10,

      .     .     .     .     .     .     .     .     .    0.10};

/* expected values and variance for each ordinal variable */

Expected = OrdMean(P) // OrdVar(P);

varNames = "X1":"X10";

print Expected[r={"Mean" "Var"} c=varNames];

/* test the RandMVOrd function */

Delta = {1.0  0.4  0.3  0.5  0.3  0.6  0.4  0.1  0.5  0.2,

         0.4  1.0  0.4  0.2  0.4  0.6  0.8  0.2  0.4  0.1,

         0.3  0.4  1.0  0.5  0.4  0.6  0.1  0.3  0.2  0.4,

         0.5  0.2  0.5  1.0  0.2  0.1  0.4  0.5  0.6  0.5,

  0.3  0.4  0.4  0.2  1.0  0.3  0.2  0.4  0.1  0.6,

  0.6  0.6  0.6  0.1  0.3  1.0  0.5  0.6  0.7  0.3,

  0.4  0.8  0.1  0.4  0.2  0.5  1.0  0.2  0.8  0.1,

  0.1  0.2  0.3  0.5  0.4  0.6  0.2  1.0  0.2  0.4,

  0.5  0.4  0.2  0.6  0.1  0.7  0.8  0.2  1.0  0.5,

  0.2  0.1  0.4  0.5  0.6  0.3  0.1  0.4  0.5  1.0};

call randseed(54321);

X = RandMVOrdinal(1000, P, Delta);

print first[label="First 1000 Obs: Multivariate Ordinal"];

thanoon
Calcite | Level 5

hello everyone

i wait your answers regarding previous post .

yours sincerely

Rick_SAS
SAS Super FREQ

The matrix Delta is supposed to represent the correlations between the ordinal variables.  Your matrix is not positive definite, and therefore is not a valid correlation matrix.  Run this code

v = eigval(Delta);

print v;

and notice that the matrix has negative eigenvalues. Generating valid correlation matrices are discussed in Chap 10 of my book. See also Appendix B, which is available from the book's web site.

thanoon
Calcite | Level 5

thank you so much dr. rick for your answer

i changed the correlation matrix to positive definite and also i have some errors please help me to correct this command.

proc iml;

load module=_all_;                     /* load the modules */

    /* P1   P2    P3    p4    p5    p6    p7    p8    p9   p10 */

P = {0.25  0.50  0.20  0.20  0.20  0.30  0.10  0.10  0.20  0.10,

     0.75  0.20  0.15  0.25  0.15  0.10  0.15  0.10  0.10  0.10,

      .    0.30  0.25  0.15  0.25  0.10  0.15  0.10  0.10  0.10,

      .     .    0.40  0.25  0.20  0.10  0.10  0.10  0.15  0.10,

      .     .     .    0.15  0.10  0.20  0.10  0.15  0.10  0.10,   

      .     .     .     .    0.10  0.10  0.10  0.15  0.10  0.10,

      .     .     .     .     .    0.10  0.15  0.10  0.15  0.10,

      .     .     .     .     .     .    0.15  0.10  0.10  0.10,

      .     .     .     .     .     .     .    0.10  0.15  0.10,

      .     .     .     .     .     .     .     .     .    0.10};

/* expected values and variance for each ordinal variable */

Expected = OrdMean(P) // OrdVar(P);

varNames = "X1":"X10";

print Expected[r={"Mean" "Var"} c=varNames];

/* test the RandMVOrd function */

Delta = {1.00 0.01 0.26 0.27 0.22 -0.08 0.39 -0.26 0.96   -0.01,

         0.01 1.00 0.06 0.09 0.07 0.18 0.18 0.16 -0.09 0.79,

         0.26 0.06 1.00 -0.17 -0.19 -0.12 0.19 -0.11 0.28 0.14,

         0.27 0.09 -0.17 1.00 0.95 0.06 0.20 0.01 0.09 0.07,

         0.22 0.07 -0.19 0.95 1.00 0.11 0.17 0.04 0.06 0.05,

        -0.08 0.18 -0.12 0.06 0.11 1.00 0.08 0.27 -0.10 0.18,

         0.39 0.18 0.19 0.20 0.17 0.08 1.00 0.05 0.31 0.12,

        -0.26 0.16 -0.11 0.01 0.04 0.27 0.05 1.00 -0.31 0.20,

         0.96 -0.09 0.28 0.09 0.06 -0.10 0.31 -0.31 1.00 -0.11,

        -0.01 0.79 0.14 0.07 0.05 0.18 0.12 0.20 -0.11 1.00);

call randseed(54321);

X = RandMVOrdinal(1000, P, Delta);

print first[label="First 1000 Obs: Multivariate Ordinal"];

Rick_SAS
SAS Super FREQ

The SAS Log reports error messages that you an use to debug your programs. If you look at the SAS Log, the error message is

ERROR 79-322: Expecting a }.

which lets you know that you need to use a curly brace to close the definition of the matrix Delta.

Let me anticipate your next question. After you fix the curly brace and run the program again, the program will still report an error. Why? Because the covariance matrix is still not valid. The covariance between two ordinal variables is special, and not all matrices are allowed. This is shown in the MV binary case on p. 154-155.

Let me give a one-dimensional analogy. For a Bernoulli random variable with probability of success p, the variance is not arbitrary. It is constrained to be p(1-p).  If you ask statistical software to simulate data from a Bernoulli process with mean 0.5 and variance 0.6, the software will report that it cannot fulfill your request because such a distribution is mathematically impossible.  That is essentially what is happening here. The software is reporting an error because you are asking it to simulate from a distribution that does not exist.

So what can you do? I know it is tempting to use simulation to "invent" data when you have none, but this example shows how difficult it is to invent data that have specify properties. Instead, the 10 variables that you are simulating should reflect real variables for which you have data.  Each variable has some number of categories and you can estimate the probability of each category from the data. If you've measured the 10 variables on dozens of subjects, you can estimate the correlations between the variables.  You can then use those empirical estimates to simulate new data that have similar properties.

thanoon
Calcite | Level 5

dear dr. rick

thank you for your help i want to invent new ordinal variables without depending on real variables how can i do that please ?

i think you mean special correlation matrix not covarince matrix because delta is a correlation matrix .i am  trying to simulate 10 ordinal variables by using data step and conducting a spearman correlation (special correlation as you said) and put it in delta and run the analysis i see same errors .

regards

Rick_SAS
SAS Super FREQ

If you just want them correlated in an arbitrary manner, you can generate MV normal data and then bin each univariate variable into some number of categories, as shown in Figure 9.2 on p. 159.

If you want to specify the exact correlation matrix for the MV distribution, then you need to work through the math to figure out a valid correlation matrix, given the number of categories and the probabilities. This is not easy, as you have discovered.

thanoon
Calcite | Level 5

this is my problem i want data with high correlation to use it in structural equation models and as you know this is the difference between sas data step and sas/iml . my question how can i manpulate the sas data step method to be like sas/iml.

in previous comment - first part if i simulate variables using mvnormal and change it to sevral categories and apply spearman correlation matrix then put it in delta matrix you think i will simulate mvordinal data without errors, regarding second part i hope to show me how can i use the hath methods.

regards

thanoon
Calcite | Level 5

dear dr. Rick

can you give me example for a delta matrix to put in the commands above please because i cannot correct this command.

regards

Rick_SAS
SAS Super FREQ

I do not know how to "invent" a 10 x 10 covariance matrix for 10 ordinal variables that have specified means and number of categories. I cannot recall any published papers that describe how to solve this problem.  The papers that I cite in Chapter 9 of my book all assume that the correlation matrix is given.

thanoon
Calcite | Level 5

then i must depend on correlation matrix given to simulate mvordinal because example in ch9 with 3 variables only and i need 10 please give me any correlation matrix.

regards

Rick_SAS
SAS Super FREQ

If any correlation matrix and mean structure suffices, Kaiser, Trager, and Leisch (2011) give a 9x9 example, which I have modified into a 10x10 example. Good luck.

%include "RandMVOrd.sas";

proc iml;
load module=_all_;                     /* load the modules */
/* 10x10 example based on Kaiser, Trager, Leisch (2011) */
m = {0.05 0.05 0.05 0.05 0.05 0.10 0.10 0.15 0.20 0.20}`;
N = nrow(m);
P = j(N,N);
P[,1] = m;
do i = 2 to N;
   P[,i] = m[ranperm(N)];
end;
v = 1 || do(0.45, 0.05, -0.05);
Corr = toeplitz(v);
print P, Corr;

X = RandMVOrdinal(100, P, Corr);
print (X[1:10,]);

thanoon
Calcite | Level 5

i am so sorry dr. theree are some errors in sas log .

proc iml;

load module=_all_;                     /* load the modules */

    /* P1    P2     P3       p4     p5      p6       p7      p8     p9      p10 */

P = {0.05 0.05 0.05 0.1    0.05 0.05 0.1     0.1 0.1    0.2,

     0.05 0.05 0.2    0.1    0.2    0.05 0.15 0.05 0.05    0.1,

     0.05 0.1    0.05 0.2    0.05 0.05 0.05 0.2    0.2    0.1,

     0.05 0.2    0.2    0.15 0.1    0.05 0.05 0.2    0.05    0.15,

     0.05 0.1    0.05 0.05 0.1    0.15 0.2    0.05 0.05 0.05,

     0.1 0.15 0.1    0.05 0.05 0.2    0.2    0.05 0.2    0.05,

     0.1 0.05 0.1    0.05 0.15 0.1    0.05 0.05 0.05 0.05,

     0.15 0.2    0.15 0.2    0.05 0.1    0.1    0.05 0.1    0.2,

     0.2 0.05 0.05 0.05 0.2    0.2    0.05 0.1    0.15 0.05,

     0.2 0.05 0.05 0.05 0.05 0.05 0.05 0.15 0.05 0.05};

/* expected values and variance for each ordinal variable */

Expected = OrdMean(P) // OrdVar(P);

varNames = "X1":"X10";

print Expected[r={"Mean" "Var"} c=varNames];

/* test the RandMVOrd function */

Delta = {1  0.45 0.4   0.35 0.3    0.25 0.2    0.15 0.1    0.05,

        0.45 1 0.45 0.4    0.35 0.3    0.25 0.2    0.15 0.1,

        0.4  0.45 1   0.45    0.4    0.35 0.3    0.25 0.2    0.15,

        0.35  0.4 0.45 1    0.45 0.4    0.35 0.3    0.25 0.2,

        0.3  0.35 0.4   0.45     1    0.45 0.4    0.35 0.3    0.25,

        0.25  0.3 0.35 0.4    0.45   1    0.45 0.4    0.35 0.3,

        0.2  0.25 0.3   0.35    0.4    0.45 1    0.45 0.4    0.35,

        0.15  0.2 0.25 0.3   0.35 0.4    0.45   1    0.45 0.4,

        0.1  0.15 0.2   0.25   0.3    0.35 0.4    0.45   1    0.45,

        0.05  0.1 0.15 0.2 0.25 0.3    0.35 0.4   0.45     1};

call randseed(54321);

X = RandMVOrdinal(1000, P, Delta);

print first[label="First 1000 Obs: Multivariate Ordinal"];

Rick_SAS
SAS Super FREQ

The program that I gave you works without error. This problem might require further time and effort on your part, but I think I've given you a good start. Good luck.

sas-innovate-2024.png

Join us for SAS Innovate April 16-19 at the Aria in Las Vegas. Bring the team and save big with our group pricing for a limited time only.

Pre-conference courses and tutorials are filling up fast and are always a sellout. Register today to reserve your seat.

 

Register now!

Multiple Linear Regression in SAS

Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.

Find more tutorials on the SAS Users YouTube channel.

From The DO Loop
Want more? Visit our blog for more articles like these.
Discussion stats
  • 23 replies
  • 2887 views
  • 5 likes
  • 2 in conversation