topic Re: Matrix functions in Proc MCMC in Statistical Procedures

Matrix functions in Proc MCMC

Mehdi_R — Sun, 05 Jul 2015 20:44:07 GMT

Hello everyone,

I have noticed that PROC MCMC is limited to apply functions on matrices. We can add, subtract, multiply matrices, or take determinant and interse of them. But if I want exponential or log of each entry of a matrix, I need to program it. For example:

array MAT[1000,1000] ;

array cov[1000,1000] ;

do i = 1 to 1000;

do j = i+1 to 1000;

cov[i, j] = exp( MAT[i, j]);

cov[j, i] = cov[i, j];

end;

cov[i,i] = exp(MAT[i, i]);

end;

This slows down the program. Don't you know any simpler or more efficient way of applying functions on matrices? How can I take the exponential of a matrix without a loop?

Re: Matrix functions in Proc MCMC

JacobSimonsen — Mon, 06 Jul 2015 21:38:37 GMT

You can make use of the matrixfunctions in fcmp. But unfortunately they can not be called immediately from a datastep or procedure. Therefore, what you can do is to make a matrix-function which call the fcmp's expmatrix-function. Your function can then be called from a datastep or from a procedure.

proc fcmp outlib = work.func.matrix;

subroutine matrixexp(m[*,*],t,y[*,*]);

outargs y;

call expmatrix(m, t, y);

endsub;

run;

option cmplib=(work.func);

data _NULL_;

array a{2,2} _temporary_;

array y{2,2} _temporary_;

call matrixexp(a,4,y);

run;

I haven't tested if the functions declared in this way can be called from PROC MCMP, but they work within some other procedures (nlmixed for instance).

I made an idea on this site about making the matrix functions directly available without first declare them with PROC FCMP. To my surprise it didnt get so positive votes :smileyshocked: https://communities.sas.com/ideas/1570

Re: Matrix functions in Proc MCMC

Mehdi_R — Tue, 07 Jul 2015 06:15:01 GMT

It was awesome! Thanks.

Re: Matrix functions in Proc MCMC

Mehdi_R — Tue, 07 Jul 2015 06:57:10 GMT

Thanks Jacob, but I think still there is problem. Please take a look at this page:

Base SAS(R) 9.2 Procedures Guide

Good news is that it works in PROC MCMC. But bad news is that the page above says that CALL EXPMATRIX does not exponentiate each element of the matrix! While I need to exponentiate each element of a matrix in PROC MCMC. Do you have another suggestion or comment?

Thanks:)

Re: Matrix functions in Proc MCMC

JacobSimonsen — Tue, 07 Jul 2015 07:19:55 GMT

Sorry, I did not read your question correct. If you want to exponentiate elementwise, then change the fcmp function to

proc fcmp outlib = work.func.matrix;

subroutine matrixexp(m[*,*],y[*,*]);

outargs y;

do i=1 to dim(m,1);

do j=1 to dim(m,2);

y[i,j]=exp(m[i,j]);

end;

endsub;

run;

Btw, you can not avoid the loop. Therefore hiding the matrix exponential within a function will not speed up the program.

Re: Matrix functions in Proc MCMC

Mehdi_R — Tue, 07 Jul 2015 08:12:17 GMT

Thanks Jacob for your quick reply.

I thought that there might be a more efficient way to calculate exponential of a matrix (matrix of exponential of entries) in PROC MCMC. Because MCMC involves heavy calculations and I wanted to avoid such loop... It seems there is no built-in function in SAS which calculates exponential of a matrix...

Do you think writing loops in proc mcmc is faster or writing loops in proc fcmp and then call the subroutine to proc mcmc?

Thanks again Jacob:)

Re: Matrix functions in Proc MCMC

JacobSimonsen — Tue, 07 Jul 2015 08:37:07 GMT

No, I think it will be almost same speed. Also, even if there was a built in function to make elementwise exponential, then there still is the loop;

Though, it could be made such that the elements are exponentiated in parallel, but don't ask me how that should be programmed

If you want to gain a bit, you can add m to outargs. I think some copying of values then can be avoided. It will not change much.

Re: Matrix functions in Proc MCMC

Mehdi_R — Tue, 07 Jul 2015 08:47:53 GMT

Thank you Jacob. It was really helpful:)