turn on suggestions

Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

Find a Community

Topic Options

- RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page

- Mark as New
- Bookmark
- Subscribe
- RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

12-03-2017 12:35 PM

Hello,

I'm a student and for a school project, i need to compute many loglikelihood of a logit model based on Spector and Mazzeo's article, and this in order to produce a 3D plot.

So for this i created a sequence of parameters and i try to compute log likelihood for each couple of parameters and used a loop for this.

But programe refuses to execute at 28th iteration and i don't understand why.

Here is my code :

```
data spector;
input grade const tuce gpa ;
cards;
0 1 20 2.66
0 1 22 2.89
0 1 24 3.28
0 1 12 2.92
0 1 17 2.86
0 1 17 2.76
0 1 21 2.87
0 1 25 3.03
0 1 20 2.63
0 1 23 3.32
0 1 23 3.57
0 1 26 3.53
0 1 19 2.74
0 1 25 2.75
0 1 19 2.83
0 1 23 3.12
0 1 22 2.06
0 1 14 2.89
0 1 26 3.51
0 1 24 2.67
0 1 21 3.10
1 1 21 4.00
1 1 29 3.92
1 1 25 3.26
1 1 25 3.16
1 1 28 3.62
1 1 24 3.54
1 1 27 2.83
1 1 17 3.39
1 1 21 3.65
1 1 23 4.00
1 1 19 2.39
;
run;
PROC IML ;
use spector;
read all var{const tuce gpa} into X;
read all var{grade} into Y ;
close spector;
print X Y ;
N = nrow(X);
*step 1 : beta sequence ;
X2=do(-1.5,1.5,0.5)`;
X3=do(1,4,0.5)`;
nX2=nrow(X2);
nX3=nrow(X3);
print nX2, nX3;
*step 2 : beta occurence;
rX2=repeat(X2,nX3,1);
rX3=repeat(X3,nX2,1);
rcst=repeat(-10.656002,nrow(rX2),1);
print rX2, rX3;
*step 3 sort;
call sort(rX3);
print rX3;
*step 4 concat;
theta=rcst||rX2||rX3;
ntheta = nrow(theta);
print theta ntheta;
*step 5 compute loglikelihood for each sequence of beta;
*Compute log likelihood of the sample;
do j=1 to 49;
XB = X*(theta`)[,j];
expXB = exp(XB);
P = expXB/(1+expXB);
logli1 = (Y#log(P))[+] ;
b=1-P;
a = 1-Y;
logli2 = (log(b##a))[+];
logli = logli1 + logli2 ;
logln = logln//logli ;
end;
print logln;
```

I know that " b##a " is not a good way but it's the only one i've found.

Someone can help me please ?

Thank you in advance for helping.

- Mark as New
- Bookmark
- Subscribe
- RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

Posted in reply to willi_m

12-03-2017 04:45 PM - edited 12-03-2017 06:47 PM

log(a * b) = log(a) + log(b).

log(a ## b) = b * log(a)

Avoid the multiplication or exponentiation followed by log to keep the numbers smaller and more manageable.

Log(a ## b) might overflow whereas b * log(a) will be better behaved.

Thanks, Rick. I corrected my typos.

- Mark as New
- Bookmark
- Subscribe
- RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

Posted in reply to willi_m

12-03-2017 06:30 PM - edited 12-04-2017 05:52 AM

Warren's advice is good, although he meant to write that

log(a # b) = log(a) + log(b)

log(a ## b) = b # log(a)

Here are some articles for you to read:

- "Two simple ways to construct a log-likelihood function in SAS"
- "Two ways to compute maximum likelihood estimates in SAS"
- To evaluate the log-likelihood at on a uniform grid of points, use the EXPANDGRID function as shown in "How to find an initial guess for an optimization"
- To graph a surface plot, see "Create a surface plot in SAS"

- Mark as New
- Bookmark
- Subscribe
- RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

Posted in reply to Rick_SAS

12-04-2017 01:57 AM - edited 12-04-2017 01:58 AM

Yes i read this somewhere in the forum.

But in writting my log as " b*log(a) " there is an error of compute.

```
332 logli2 = (a*log(b))[+];
ERROR: (execution) Matrices do not conform to the operation.
operation : * at line 332 column 16
operands : a, _TEM1001
a 32 rows 1 col (numeric)
_TEM1001 32 rows 1 col (numeric)
statement : ASSIGN at line 332 column 5
```

There are no error when i write " a#log(b) " and out of the loop

```
XB = X*(theta`)[,1];
expXB = exp(XB);
P = expXB/(1+expXB);
logli1 = (Y#log(P))[+] ;
b=1-P;
a = 1-Y;
logli2 = (a#log(b))[+]; /* HERE */
logli = logli1 + logli2 ;
logln = logln//logli ;
```

But when i want to implemant this in a loop there is a problem

```
345 do j=1 to 49;
346 XB = X*(theta`)[,j];
346!
347 expXB = exp(XB);
347!
348 P = expXB/(1+expXB);
348!
349 logli1 = (Y#log(P))[+] ;
350 b=1-P;
351 a = 1-Y;
352 logli2 = (a#log(b))[+];
353
354 logli = logli1 + logli2 ;
355 logln = logln//logli ;
356 end;
ERROR: (execution) Invalid argument to function.
count : number of occurrences is 2
operation : LOG at line 352 column 20
operands : b
b 32 rows 1 col (numeric)
statement : ASSIGN at line 352 column 5
```

Moreover, my loop is ok for a less number of iteration (13)

```
429 do j=1 to 13; /*13 iterations OK*/
430 XB = X*(theta`)[,j];
430!
431 expXB = exp(XB);
431!
432 P = expXB/(1+expXB);
432!
433 logli1 = (Y#log(P))[+] ;
434 b=1-P;
435 a = 1-Y;
436 logli2 = (a#log(b))[+];
437
438 logli = logli1 + logli2 ;
439 logln = logln//logli ;
440 end;
```

And with 14 iterations : problems

```
441 do j=1 to 14; /* 14 iterations PROBLEMS */
442 XB = X*(theta`)[,j];
442!
443 expXB = exp(XB);
443!
444 P = expXB/(1+expXB);
444!
445 logli1 = (Y#log(P))[+] ;
446 b=1-P;
447 a = 1-Y;
448 logli2 = (a#log(b))[+];
449
450 logli = logli1 + logli2 ;
451 logln = logln//logli ;
452 end;
ERROR: (execution) Invalid argument to function.
count : number of occurrences is 2
operation : LOG at line 448 column 20
operands : b
b 32 rows 1 col (numeric)
statement : ASSIGN at line 448 column 5
```

So i don't understand where is the problem, and this is problematic because after i would like to compute more of 160 000 loglikelihood for my 3D plot in order to obtain a smoother graph.

- Mark as New
- Bookmark
- Subscribe
- RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

Posted in reply to willi_m

12-06-2017 09:31 AM

The error seems to indicate that you are trying to evaluate the LOG function for an invalid parameter. Note that if P=0 then LOG(P) is undefined and if P=1 then LOG(1-P) is undefined. In your case, I think you are running into the P=1 problem for certain values of j.

I don't completely understand what you are doing, but it looks like you are evaluating the logistic log-likelihood function for linear combinations of the variables for a grid of (X2, X3) coffficients. For some values of j, the value of XB is so large that logistic(XB) is numerically equal to 1. Try using a more restricted range of coefficients, such as

X2=do(-1,1,0.5)`;