[R] Symbolic Maximum Likelihood in R

[Gabor Grothendieck]

On 7/3/05, Doran, Harold <HDoran at air.org> wrote:
> Dear List:
> 
> Is any one aware of a package that would extend the D() function and allow for one to maximize a likelihood function symbolically? Something akin to Solve[x==0, parameter] function in Mathematica?
> 
> Clearly R has the capacity to _compute_ MLEs given a set of data. But, I'm looking for a package that would allow for me to define the likelihood function, find the 1st order partial derivative of this function (which can already be handled in D()), and then symbolically maximize this function?
> 

There are a number of free symbolic mathematics packages such as
Maxima and Yacas.  Here is a simple example in Yacas.  (Also
see http://finzi.psych.upenn.edu/R/Rhelp02a/archive/31418.html
for transferring functions from yacas to R).

C:\usr\yacas> yacas
[...various startup messages...]
In> f(x) := Exp(-(x-mu)^2)/Sqrt(2*Pi)
Out> True;
In> logf(x) := Ln(f(x));
Out> True;
In> loglik := logf(x1) + logf(x2) + logf(x3);
Out> Ln(Exp(-(x1-mu)^2)/Sqrt(2*Pi))+Ln(Exp(-(x2-mu)^2)/Sqrt(2*Pi))+Ln(Exp(-(x3-m
u)^2)/Sqrt(2*Pi));
In> score := D(mu) loglik;
Out> (-(-2)*(x2-mu)*Exp(-(x2-mu)^2)*2*Pi)/(2*Pi*Exp(-(x2-mu)^2))-((-2)*(x1-mu)*E
xp(-(x1-mu)^2)*2*Pi)/(2*Pi*Exp(-(x1-mu)^2))-((-2)*(x3-mu)*Exp(-(x3-mu)^2)*2*Pi)/
(2*Pi*Exp(-(x3-mu)^2));
In> score := Simplify(score);
Out> 2*(x2+(-3)*mu+x1+x3);
In> Solve(score == 0, mu);
Out> (x2+x1+x3)/3;