Help using Base SAS procedures

Comb Function Limitations

Accepted Solution Solved
Reply
New Contributor
Posts: 2
Accepted Solution

Comb Function Limitations

I'm trying to use an equation that involves combination calculations (the COMB function in SAS). The equation results in very large numbers, but in the end the value is between 0 and 1. SAS can't handle large values within the COMB function so I keep receiving errors. For example, COMB(1045,600) is acceptable, but SAS doesn't like COMB(1046,600) since I assume it is too large. Anyone know a way to get around this? The equation I'm using is:

COMB(k,i)*COMB(N-k,n-i)/COMB(N,n)

Where k = 1044, i = 600, N = 1046, n = 600

Thanks


Accepted Solutions
Solution
‎11-06-2014 09:05 AM
SAS Super FREQ
Posts: 3,475

Re: Comb Function Limitations

First, SAS is a case-insensitive language, so N and n represent the same variable in your expression. In the code below I replace n by nn.

Second, one of the basic tenets of numerical analysis is that you should avoid computing factorials and combinatorial expressions in the naive way, because these expressions are prone to overflow.

Third, expressions like this occur often in probability models. Your expression looked familiar to me, and I fortunately recalled that this expresion is the probability density function for the hypergeometric distribution.  Thus your problem is solved in SAS by calling the built-in PDF function, like this:

N=1046;

k=1044;

nn = 600;

i = 600;

h = pdf("Hypergeom", i, N, k, nn);  /* = COMB(k,i)*COMB(N-k,nn-i)/COMB(N,nn) */

I think an important lesson to learn is that it is helpful for the OP to describe what he is doing and why he is doing it. The context is often useful in guiding those who are trying to answer the question.

View solution in original post


All Replies
Super User
Posts: 10,476

Re: Comb Function Limitations

I get Comb(N,n) as 9.176E307. when N=1045. Which is right at the limit of the system. Increasing N>1045 and n=600 exceeds my system capability. (1.798E308)

You can see what your limit is with

data _null_;

x= constant('BIG');

put x;

run;

Respected Advisor
Posts: 3,777

Re: Comb Function Limitations

PROC Star
Posts: 7,357

Re: Comb Function Limitations

As always, data_null came up with what I think is an excellent suggestion for addressing your problem. However, in the example you provided, you have an instance where N-n equals zero. That won't make sense for a comb function, but one also can't obtain the log of 0.

Additionally, are you on a system that can differentiate between upper and lower case variable names? You may want to rename n to be n2.

Trusted Advisor
Posts: 1,300

Re: Comb Function Limitations

This is because the answer goes outside the bounds of the double data type.

Trusted Advisor
Posts: 1,300

Re: Comb Function Limitations

proc groovy;

submit;

def binomial(N, K) {

    ret = BigInteger.ONE;

    for (k = 0; k < K; k++) {

        ret = ret.multiply(BigInteger.valueOf(N-k))

                 .divide(BigInteger.valueOf(k+1));

    }

    return ret;

}

System.out.println(binomial(1046, 600));

endsubmit;

quit;

215214135164432318969664923845300683858831944672653837461277170046182567034643018347761262897990464760845987368908944618207241161102715206382228014600474147009141057204830063147454042567661669186964749079231644777086378723140490665422937054079891989213282500050827623980250811932790782038439267537485639619200

Respected Advisor
Posts: 4,641

Re: Comb Function Limitations

Data_Null_ provided the most useful suggestion, but also:

depending on the relative size of your variables, you may get away with the equivalent expression

comb(N-n, k-i) / comb(N, k) * comb(n, i)

Notice that it is also a good idea to divide before you multiply.

PG

PG
Respected Advisor
Posts: 3,777

Re: Comb Function Limitations

How would the OPs expression or your equivalent need to be changed to use LCOMB?

Respected Advisor
Posts: 4,641

Re: Comb Function Limitations

For the sake of completeness, that would be:

EXP( LCOMB(k,i) + LCOMB(N-k,n-i) - LCOMB(N,n) )


PG

PG
Solution
‎11-06-2014 09:05 AM
SAS Super FREQ
Posts: 3,475

Re: Comb Function Limitations

First, SAS is a case-insensitive language, so N and n represent the same variable in your expression. In the code below I replace n by nn.

Second, one of the basic tenets of numerical analysis is that you should avoid computing factorials and combinatorial expressions in the naive way, because these expressions are prone to overflow.

Third, expressions like this occur often in probability models. Your expression looked familiar to me, and I fortunately recalled that this expresion is the probability density function for the hypergeometric distribution.  Thus your problem is solved in SAS by calling the built-in PDF function, like this:

N=1046;

k=1044;

nn = 600;

i = 600;

h = pdf("Hypergeom", i, N, k, nn);  /* = COMB(k,i)*COMB(N-k,nn-i)/COMB(N,nn) */

I think an important lesson to learn is that it is helpful for the OP to describe what he is doing and why he is doing it. The context is often useful in guiding those who are trying to answer the question.

Trusted Advisor
Posts: 1,300

Re: Comb Function Limitations

Great post Rick.

PROC Star
Posts: 7,357

Re: Comb Function Limitations

I was impressed by all 4 alternative solutions!

Trusted Advisor
Posts: 1,300

Re: Comb Function Limitations

I would say my post was more for amusement than any real practicality.

SAS Super FREQ
Posts: 3,475

Re: Comb Function Limitations

In the same vein, here is a blog post from 2011 that uses SAS to define and manipulate "BigIntegers" for computing the factorials of  numbers when the factorial value exceeds the value that can be stored in a double (1.798E308).  I could write a similar program for combinations.  Notice that I use the LFACT function, which is the analog to the LCOMB function that Pierre and Data_Null use here.  Like Fried Egg's post, this way of solving the problem is more for amusement than for practical scientidic computing, so I think the post uses all of the features that Art likes.

New Contributor
Posts: 2

Re: Comb Function Limitations

This solved the issue I was having. Thank you!

☑ This topic is SOLVED.

Need further help from the community? Please ask a new question.

Discussion stats
  • 14 replies
  • 551 views
  • 7 likes
  • 7 in conversation