EZ Study

Actuarial Biology Chemistry Economics Calculators Confucius Engineer
Physics
C.S.

Newton-Raphson Algorithm in SAS
Find the Bayesian posterior MLE estimate in SAS

Introduction: In numerical analysis, Newton's method(also known as the Newton
-Raphson method
), named after Isaac Newton and Joseph Raphson, is a method for finding successively better approximations to the roots (or zeroes) of a real-valued function. How can we apply that in SAS? Where or when we can apply? Find the Bayesian posterior estimate in SAS might be a good place.
Example: In this tutorial, we want to apply the Newton Raphson algorithm to find the Bayesian posterior MLE estimate.

Based on previous tutorial, on Bayesian analysis, the conjugate prior for binormial
distribution is beta distribution
. Assume the batting probability(p)) of a group of players follows beta distribution with parameters (a,b), to find the MLE estimate of those parameters, we can apply the Newton-Raphson algorithm.

Given the batting results of the group of players, say, player i, hiting K_i out of total N_i battings; the corresponding batting probability for player i is: K_i/N_i. We can simply apply univariate procedure to get the beta estimators, just as in the previous tutorial. In that case, 2/20 will have the same effect as 20/200 since they are both p=0.2. To make good use of those K_i and N_i, we have to write down the likelihood function, find the MLE of (a,b).
``````

%macro Newton_beta(datain,n,s);
**  get starting values for (a) and (b) using the mean and variance *;

data Datain;
set &datain(where=(&s<=&n and &n>0));
rate = &s/&n;
run;
proc means data=Datain(where=(&n>0)) noprint;
var rate;
output out=Dataout mean=mean stddev=stddev;
run;

data _null_;
set Dataout;
if stddev**2>mean*(1-mean) then stddev=sqrt(0.9*mean*(1-mean));
a = mean*((mean/(stddev**2))*(1-mean)-1);
b = a*(1-mean)/mean;
call symput("a",a);
call symput("b",b);
run;
``````
Acknowledgement: The original SAS code is from the work by Dr. Speights.