# [R] Minimizing a Function with three Parameters

Uwe Ligges ligges at statistik.uni-dortmund.de
Fri Dec 2 09:04:19 CET 2005

voodooochild at gmx.de wrote:

> voodooochild at gmx.de wrote:
>
>
>>Hi,
>>
>>I'm trying to get maximum likelihood estimates of \alpha, \beta_0 and
>>\beta_1, this can be achieved by solving the following three equations:
>>
>>n / \alpha + \sum\limits_{i=1}^{n} ln(\psihat(i)) -
>>\sum\limits_{i=1}^{n} ( ln(x_i + \psihat(i)) ) = 0
>>
>>\alpha \sum\limits_{i=1}^{n} 1/(psihat(i)) - (\alpha+1)
>>\sum\limits_{i=1}^{n} ( 1 / (x_i + \psihat(i)) ) = 0
>>
>>\alpha \sum\limits_{i=1}^{n} ( i / \psihat(i) ) - (\alpha + 1)
>>\sum\limits_{i=1}^{n} ( i / (x_i + \psihat(i)) ) = 0
>>
>>where \psihat=\beta_0 + \beta_1 * i. Now i want to get iterated values
>>for \alpha, \beta_0 and \beta_1, so i used the following implementation
>>
>># first equation
>>l1 <- function(beta0,beta1,alpha,x) {
>> n<-length(x)
>> s2<-length(x)
>>   for(i in 1:n) {
>>   s2[i]<-log(beta0+beta1*i)-log(x[i]+beta0+beta1*i)
>>   }
>> s2<-sum(s2)
>> return((n/alpha)+s2)
>>}
>>
>># second equation
>>l2 <- function(beta0,beta1,alpha,x) {
>> n<-length(x)
>> s1<-length(x)
>> s2<-length(x)
>>   for(i in 1:n) {
>>   s1[i]<-1/(beta0+beta1*i)
>>   s2[i]<-1/(beta0+beta1*i+x[i])
>>   }
>> s1<-sum(s1)
>> s2<-sum(s2)
>> return(alpha*s1-(alpha+1)*s2)
>>}
>>
>>#third equation
>>l3 <- function(beta0,beta1,alpha,x) {
>> n<-length(x)
>> s1<-length(x)
>> s2<-length(x)
>>   for(i in 1:n) {
>>   s1[i]<-i/(beta0+beta1*i)
>>   s2[i]<-i/(x[i]+beta0+beta1*i)
>>   }
>> s1<-sum(s1)
>> s2<-sum(s2)
>> return(alpha*s1-(alpha+1)*s2)
>>}
>>
>># all equations in one
>>gl <- function(beta0,beta1,alpha,x) {
>> l1(beta0,beta1,alpha,x)^2 + l2(beta0,beta1,alpha,x)^2 +
>>l3(beta0,beta1,alpha,x)^2
>>}
>>
>>#iteration with optim
>>optim(c(1,1,1),gl,x)
>>
>>i get always an error massage. Is optim anyway the 'right' method to get
>>all three parameters iterated at the same time?
>>
>>best regards
>>Andreas
>>
>>______________________________________________
>>R-help at stat.math.ethz.ch mailing list
>>https://stat.ethz.ch/mailman/listinfo/r-help
>>
>>
>>
>>
>
> hi sundar,
>
>
> now i have another model where instead of i    i2 is used, but i don't
> now way i got so large estimates?
>
> x<-c(10,8,14,17,15,22,19,27,35,40)
>
> # 1.Gleichung
>
> l1 <- function(beta0,beta1,alpha,x) {
>  n<-length(x)
>  s2<-length(x)
>    for(i in 1:n) {
>    s2[i]<-log(beta0+beta1*i2)-log(x[i]+beta0+beta1*i2)
>    }
>  s2<-sum(s2)
>  return((n/alpha)+s2)
> }
>
>
> # 2.Gleichung
>
> l2 <- function(beta0,beta1,alpha,x) {
>  n<-length(x)
>  s1<-length(x)
>  s2<-length(x)
>    for(i in 1:n) {
>    s1[i]<-1/(beta0+beta1*i2)
>    s2[i]<-1/(beta0+beta1*i2+x[i])
>    }
>  s1<-sum(s1)
>  s2<-sum(s2)
>  return(alpha*s1-(alpha+1)*s2)
> }
>
> # 3.Gleichung
>
> l3 <- function(beta0,beta1,alpha,x) {
>  n<-length(x)
>  s1<-length(x)
>  s2<-length(x)
>    for(i in 1:n) {
>    s1[i]<-(i2)/(beta0+beta1*i2)
>    s2[i]<-(i2)/(x[i]+beta0+beta1*i2)
>    }
>  s1<-sum(s1)
>  s2<-sum(s2)
>  return(alpha*s1-(alpha+1)*s2)
> }
>
> # Zusammenfügen aller Teile
>
> gl <- function(beta,x) {
>  beta0<-beta[1]
>  beta1<-beta[2]
>  alpha<-beta[3]
>  v1<-l1(beta0,beta1,alpha,x)2
>  v2<-l2(beta0,beta1,alpha,x)2
>  v3<-l3(beta0,beta1,alpha,x)2
>  v1+v2+v3
> }
>
>
>
> optim(c(20000,6000,20000),gl,x=x,control=list(reltol=1e-12))
>
> the values should be alpha=20485, beta0=19209 and beta1=6011
>
> and another point is, what is a good method to find good starting values
> for 'optim'. it seems, that i only get the desired values when the
> starting values are in the same region. I used
> control=list(reltol=1e-12), but it seems, that then it is also important
> to have the starting values in the same region as the the desired values.

Depends on the response surface. If the latter behaves well, you can be
lucky, if it does behave really ill, you will get grey hair during your
further research on this topic ...

Uwe Ligges

> regards
> andreas
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help