Are macro variables supported as input of SAS/IML modules and routines?

[Season]

I am conducting numerical integration (or more specifically, calculating the value of incomplete gamma function) via the QUAD routine of SAS/IML. Given that I have a batch of functions to evaluate, I composed a macro in which the integrand is defined repetitively via a module, like this:

%macro fun;
%do i=1 %to 2;
p&i.=p[&i.];/*p is a vector housing variables closely connected to the power of one of the terms of the integrand of the incomplete gamma function*/
p&i.a=p&i.-1;
start fun(u);
y=u**(eval(1*p&i.a))*exp(-u);/*incomplete gamma function*/
return (y);
finish;
a=j(1,1,0)||up[&i.];/*up is a vector of the upper limits of integration calculated previously*/
call quad(g&i.,'fun',a);
%end;
%mend;
proc iml;
%fun;
print g1 g2;
quit;

Only two loops were written in the macro %fun for I was troubleshooting. In practice, more loops are needed.

However, the SAS log reports that g1 and g2 had not been set to values. I read the log carefully and found a message:

p1a      0 row       0 col     (type ?, size 0)
p2a      0 row       0 col     (type ?, size 0)

So it seems that the scalars named p1a and p2a and served as the power of the independent variable of one of the terms of the integrand of the incomplete gamma function did not pass to either the fun module or the QUAD routine. Of course, it might also be the case that the value was passed to none of them.

Given that a do loop is not acceptable in a module of SAS/IML (I tried this method yesterday and it is because of this that I attempted to compose this macro), what should I do next?

Many thanks!

[Rick_SAS · Accepted solution]

Short answer: Yes they are. However, your code contains some errors:

1. There is not any function called EVAL in IML. You might be intending to use %EVAL, but that macro isn't necessary

2. If you want to pass parameters to the integrated used in the QUAD routine, you need to use the GLOBAL clause.

> Given that a do loop is not acceptable in a module of SAS/IML, what should I do next?

This statement is false. What you should do is use the IML language instead of macro loops. An IML program rarely requires a macro loop because the language supports loops and arrays.

I think you are asking how to write one integrated function that can be used in a loop to integrate a sequence of functions. For each call, you want the parameter and the upper bound of integration to vary according to values in some array.  Please study the following:

proc iml;
/* loop over parameters and upper bounds for integration */
start fun(u) global(g_alpha);
   y=u**(g_alpha-1)*exp(-u);/*incomplete gamma function*/
   return (y);
finish;

up = {3, 4, 5};
p = {1.1, 2.2, 0.8};
IncGamma = j(nrow(p), 1,.);
do i = 1 to nrow(p);
   g_alpha = p[I];   /* copy parameter into a global variable to the integrand */
   call quad(g, "fun", 0 || up[I]);  /* set the domain of integration */
   IncGamma[i] = g;  /* copy the answer into a vector */
end;

print p up IncGamma;

If I may make a suggestion, you don't need to use the QUAD function to solve this problem. Mathematically, the integral you are trying to compute is related to the CDF of the gamma distribution, which is built into SAS. Thus, you can solve your problem by using the following :

proc iml;
/* it's not necessary to call QUAD; use the CDF of the gamma distribution.
   See https://blogs.sas.com/content/iml/2021/01/25/incomplete-gamma-function-sas.html*/
start LowerIncGamma(x, alpha); /* lower incomplete gamma function */
   return Gamma(alpha) # CDF('Gamma', x, alpha);
finish;

up = {3, 4, 5};
p = {1.1, 2.2, 0.8};
IncGamma2 = LowerIncGamma(up, p); 
print p up IncGamma2;

Some references:

• GLOBAL clause: Evaluate an iterated integral in SAS
• How to compute the incomplete gamma function in SAS - The DO Loop
• Macros and loops in the SAS/IML language - The DO Loop

[Season]

Thank you, Rick, for your rapid, detailed and helpful reply! It takes me time to digest the information you provided in your thread. However, there is one issue that I can point out:

---

@Rick_SAS wrote:

> Given that a do loop is not acceptable in a module of SAS/IML, what should I do next?

This statement is false. What you should do is use the IML language instead of macro loops. An IML program rarely requires a macro loop because the language supports loops and arrays.

I made a mistake in my post. The correct message to convey is "building modules are not supported in a do loop". For instance, the following code makes SAS report an error in the log:

proc iml;
do i=1 to 10;
start fun(u);
y=u**2;
return(y);
finish;
end;
quit;

I will post my questions if there are any after reading your reply carefully. Thanks!

[Tom]

Poor example, why would you want to define the same function 10 times?

But dynamically generating the code that defines a function is not part of the language. Or really of much value since normally if you want a function to behave differently you give it different inputs.

However if you did need to build something were the definition of a function changes then just use some type of code generation.  Such as macro code or call execute() or writing the code to a file and using %INCLUDE.

[Season]

---

@Tom wrote:

Poor example, why would you want to define the same function 10 times?

---

The original code I composed revolved around the incomplete gamma function. As the integrand varied from iteration to iteration and I was unaware of the correct usage of the GLOBAL option of the START statement before I consulted @Rick_SAS, I tried to use the do loop to define the function differently in each loop, but I failed and therefore deleted the original code. I did not want to spend time on recalling and typing the real code I used before, so I generated the same function again and again. I think this is not a major concern to focus on. The real major concern is that building modules are not supported in the DO loop, an issue that seems to have been overlooked by IML textbooks, as least the one I referred to. Therefore, I wrote down my experience to alert to readers who may be willing to adopt a similar approach.

---

@Tom wrote:

However if you did need to build something were the definition of a function changes then just use some type of code generation.  Such as macro code or call execute() or writing the code to a file and using %INCLUDE.

---

@Rick_SAS gave an example of the GLOBAL option. I think that is good enough. I had just read a blog by @Rick_SAS (Macros and loops in the SAS/IML language - The DO Loop) in which he said that he did not support using macros in SAS/IML that much. Neither do I. I am not that familiar with macros.

[Rick_SAS]

> The real major concern is that building modules are not supported in the DO loop, an issue that seems to have been overlooked by IML textbooks, as least the one I referred to.

I would be interested to learn the names of the IML textbooks that you consulted.

Some of my writings about SAS/IML modules are compiled in the article, https://blogs.sas.com/content/iml/2015/06/17/everything-about-modules.html,
which states, 'All modules are global in scope. The SAS/IML language does not support "hidden modules" or "local modules" or "nested modules."'

[Season]

I have just finished reading the book A SAS/IML Companion for Linear Models | SpringerLink from cover to cover. That is the only IML textbook I have read till now. I had known that you had written an IML book. But by the time I was selecting which book to read, I was not in a dire need for the software. I learned it because I wished to have more "knowledge (or technology) reserve" in response to demanding requests that may be proposed by reviewers after submission of my research papers to journals. Given that your book was thicker, I chose the book with less pages.

[Season]

I read your blog and found that given there was an extra Γ(p) term to divide after I compute the incomplete gamma function, the whole thing I had been calculating was nothing but the cumulative distribution function (CDF) of Γ(up,p). I tried the numerical integration method by the QUAD routine first and found out that for small values of p, the QUAD routine failed potentially as a result of lack of specification of the SCALE= option. I therefore read your blog and found out the conclusion I gave at the beginning of this thread.

Still, I encountered a problem when I was actually computing the CDF's. I directly applied the CDF function element-wise and found an error. The code looked like this:

proc iml;
Gamma=j(nrow(p),1,1);
up = {3, 4, 5};
p = {1.1, 2.2, 0.8};
do i=1 to nrow(p);
Gamma[i]=cdf('gamma',up[i],p[i]);
end;
print p up Gamma;
quit;

The code failed. The log read:

ERROR: (execution) Invalid argument to function.

 operation : CDF at line 358 column 9
 operands  : *LIT1008, _TEM1001, _TEM1002

*LIT1008      1 row       1 col     (character, size 5)

 gamma

_TEM1001      1 row       1 col     (numeric)

 -2.326E19

_TEM1002      1 row       1 col     (numeric)

 -4.868E17

 statement : ASSIGN at line 358 column 1

Why is it the case?

[Rick_SAS]

> Why is it the case?

In your program, you have not yet defined the p vector when you attempt to allocate the Gamma vector. In the first line

Gamma=j(nrow(p),1,1);

the expression nrow(p) evaluates to 0. Therefore the vector has "zero rows."  Inside the DO loop, the expression Gamma[I] gives an error because the subscript Gamma[1] is invalid.

Solution: Move the line 

Gamma=j(nrow(p),1,1);

to the third line in the program, not the first.

[Season]

Ah, yes... That's a mistake that shouldn't had been make... But the problem lingers after I rearranged the sequence of codes. The error-reporting log was similar to the one I posted earlier. It still said there was an invalid argument to function. So, ..., why is the case now?

[Rick_SAS]

I cannot guess what your code looks like. The following works:

proc iml;
up = {3, 4, 5};
p = {1.1, 2.2, 0.8};
Gamma=j(nrow(p),1,1);
do i=1 to nrow(p);
   Gamma[i]=cdf('gamma',up[i],p[i]);
end;
print p up Gamma;

[Season]

That's interesting. My code exactly looks like the one you provided but it still kept failing. I will copy my log once again.

ERROR: (execution) Invalid argument to function.

 operation : CDF at line 180 column 9
 operands  : *LIT1009, _TEM1001, _TEM1002

*LIT1009      1 row       1 col     (character, size 5)

 GAMMA

_TEM1001      1 row       1 col     (numeric)

 -2.326E19

_TEM1002      1 row       1 col     (numeric)

 -4.868E17

 statement : ASSIGN at line 180 column 1

I looked at the entries of up and p and found minuscule numbers at the end of the two vectors (e.g.,1.4066E19). Is it because of the presence of the letter "E" that causes the CDF function to falsely identify it as a character variable and the consequent failure of the code?

[Season]

Now I see where the problem seems to be. The entire but the last few entires of the vector of results (i.e., the CDF's I wish to calculate) are in fact 1.

For instance, one of the calculations is as follows:

a=1817.9997;
b=38.043169;
gamma=cdf('GAMMA',a,b);

But given that I had pre-emptively filled the vector gamma with a batch of 1's with the J function, I was not able to discover that SAS had already conducted calculations on some of the entires. Setting the default values to 0 or missing helps me make this discovery.

But the last few entries of up and p contained negative values, which is why SAS kept reporting errors. I will go through my code and search for reasons of the appearance of negative values.