R help - Jun 2008 - Maximum Likelihood Estimation

Using R, I would like to calculate algorithms to estimate coefficients á and â
within the gamma function: f(costij)=((costij)^á)*exp(â*costij).  I have its
logarithmic diminishing line data (Logarithmic Diminishing Line Data Table) and
have installed R¢s Maximum Likelihood Estimation package; however, I am unsure
which method to apply in order to calculate the algorisms (i.e., Newton-Raphson
Maximization, Nelder-Mead Maximization, etc.)  Any guidance you all could
provide would be appreciated.
 
Logarithmic Diminishing Line Data Table
1              0.983385666
2              0.578408021
3              0.421101402
4              0.334555838
5              0.278826347
6              0.239521701
7              0.210100828
8              0.187133124
9              0.168633944
10           0.15336976
11           0.140530624
12           0.129561113
13           0.120066673
14           0.111758725
15           0.104420943
16           0.097887721
17           0.092029955
18           0.086745393
19           0.081951911
20           0.077582726
21           0.073582928
22           0.069906916
23           0.066516481
24           0.06337934
25           0.060468017
26           0.057758964
27           0.055231875
28           0.052869138
29           0.050655392
30           0.048577177
31           0.046622642
32           0.044781306
33           0.043043865
34           0.041402028
35           0.03984838
36           0.038376265
37           0.036979694
38           0.03565326
39           0.034392065
40           0.033191667
41           0.03204802
42           0.030957435
43           0.029916538
44           0.02892224
45           0.027971702
46           0.027062313
47           0.026191669
48           0.025357547
49           0.024557893
50           0.023790804
51           0.023054514
52           0.022347383
53           0.021667883
54           0.021014591
55           0.02038618
56           0.01978141
57           0.019199121
58           0.018638226
59           0.018097707
60           0.017576607
61           0.017074029
62           0.016589127
63           0.016121106
64           0.015669216
65           0.015232751
66           0.014811043
67           0.014403462
68           0.014009411
69           0.013628326
70           0.013259673
71           0.012902944
72           0.01255766
73           0.012223364
74           0.011899623
75           0.011586025
76           0.011282179
77           0.010987712
78           0.01070227
79           0.010425513
80           0.01015712
81           0.009896783
82           0.009644208
83           0.009399116
84           0.009161239
85           0.008930321
86           0.008706116
87           0.008488392
88           0.008276923
89           0.008071496
90           0.007871904
91           0.00767795
92           0.007489446
93           0.007306209
94           0.007128067
95           0.006954853
96           0.006786404
97           0.006622569
98           0.006463198
99           0.00630815
100         0.006157287
 
Thanks,
Todd


      
	[[alternative HTML version deleted]]

Todd Brauer <toddbrauer <at> yahoo.com> writes:
> 
> Using R, I would like to calculate algorithms to estimate coefficients ?
and ?
within the gamma function:
> f(costij)=((costij)^?)*exp(?*costij).? I have its logarithmic diminishing
line
data
> (Logarithmic Diminishing Line Data Table) and have installed R?s Maximum
Likelihood Estimation
> package; however, I am unsure which method to apply in order to calculate
the
algorisms (i.e.,
> Newton-Raphson Maximization, Nelder-Mead Maximization, etc.)? Any guidance
you
all could provide
> would be appreciated.
  For a simple, well-behaved problem (like this one) it doesn't
really matter that much what algorithm you pick. optim()'s default
is Nelder-Mead (typically slightly more robust and less efficient).
Here's some sample code to do a least-squares fit:

x <- read.table("gammafit.txt")
names(x) <- c("x","y")
plot(x)
sum(x)
gfun <- function(p) {
    a <- p[1]
    b <- p[2]
    exp.y <- x$x^b*exp(-a*x$x)
    sum((exp.y-x$y)^2)
}
opt1 <- optim(gfun,par=c(a=0.1,b=1))
v <- opt1$par
curve(x^v["b"]*exp(-v["a"]*x),add=TRUE,col=2)

  However, if you really want to do MLE you may need to specify
your problem a little more carefully -- sum of squares will
find a reasonable answer but likelihood-based inference won't
be right.  Is this a survival problem (in which case you
might want to use the survival package)?  If you construct
your own likelihood function, you might need
to use pgamma(...,lower.tail=FALSE).

  Ben Bolker