karsten_krug at gmx.net
2007-Jan-04 10:24 UTC
[Rd] problem with function 'optimise' (PR#9438)
Full_Name: Karsten Krug Version: 2.4.0 OS: Open Suse 10.0, Windows XP Submission from: (NULL) (88.134.13.50) I found a problem in the 'optimise' function for one dimensional optimisation. Example 1: Try to find a maximum of the function below with the use of 'optimise' in the interval [0,0.5]. The function follows a parabola and has two local maxima located at the margins of the interval [0,0.5]. Please try the following code which performs this optimisation, produces a plot of the function within the given interval and adds the computed maximum returned by 'optimise()' ex1 <- function(x) log(4 * x^2 - 2 * x + 0.5) opt1 <- optimise(ex1, lower=0, upper=0.5, maximum=T) x <- seq(0,0.5,0.001) plot(x, ex1(x), type="l") abline(v=opt1$maximum, col="red") The returned value of the maximum is 0.309017 which is definitely wrong. Example 2: The next function has also two equivalent maxima this time not located at the margins of the interval [0, 9/41]. Thus, the maxima are real analytic extrema. The function 'optimise()' again returns a wrong value. ex2 <- function(x) log((18/41) * x - 2 * x^2) + 16 * log(4 * x^2 - (36/41) * x + (9/41)) + 24 * log((23/82) + (18/41) * x - 2 * x^2) opt2 <- optimise(ex2, lower=0, upper=(9/41), maximum=T) x <- seq(0,9/41,0.001) plot(x,ex2(x),type="l") abline(v=opt2$maximum, col="red") The two functions are not only of theoretical interest but have to be optimised as likelihood functions when dealing with SNP data. I observed this problem on Windows XP and Open Suse 10.0. Other platforms have to be tested. Some informations: Open Suse 10.0 platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major 2 minor 4.0 year 2006 month 10 day 03 svn rev 39566 language R version.string R version 2.4.0 (2006-10-03) Windows XP: platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 4.0 year 2006 month 10 day 03 svn rev 39566 language R version.string R version 2.4.0 (2006-10-03) With best regards, Karsten Krug Student at University of Leipzig
the issue here is a numerical one; for instance in your first example
check:
eps <- sqrt(.Machine$double.eps)
opt1 <- optimise(ex1, lower = 0 - eps, upper = 0.5, maximum = TRUE)
# or
opt1 <- optimise(ex1, lower = 0 + eps, upper = 0.5, maximum = TRUE)
which gives the correct results.
I hope it helps.
Best,
Dimitris
----
Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven
Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/(0)16/336899
Fax: +32/(0)16/337015
Web: http://med.kuleuven.be/biostat/
http://www.student.kuleuven.be/~m0390867/dimitris.htm
----- Original Message -----
From: <karsten_krug at gmx.net>
To: <r-devel at stat.math.ethz.ch>
Cc: <R-bugs at biostat.ku.dk>
Sent: Thursday, January 04, 2007 11:24 AM
Subject: [Rd] problem with function 'optimise' (PR#9438)
> Full_Name: Karsten Krug
> Version: 2.4.0
> OS: Open Suse 10.0, Windows XP
> Submission from: (NULL) (88.134.13.50)
>
>
>
> I found a problem in the 'optimise' function for one dimensional
> optimisation.
>
>
> Example 1:
> Try to find a maximum of the function below with the use of
> 'optimise' in the
> interval [0,0.5]. The function follows a parabola and has two local
> maxima
> located at the margins of the interval [0,0.5]. Please try the
> following code
> which performs this optimisation, produces a plot of the function
> within the
> given interval and adds the computed maximum returned by
> 'optimise()'
>
> ex1 <- function(x) log(4 * x^2 - 2 * x + 0.5)
> opt1 <- optimise(ex1, lower=0, upper=0.5, maximum=T)
> x <- seq(0,0.5,0.001)
> plot(x, ex1(x), type="l")
> abline(v=opt1$maximum, col="red")
>
> The returned value of the maximum is 0.309017 which is definitely
> wrong.
>
>
> Example 2:
> The next function has also two equivalent maxima this time not
> located at the
> margins of the interval [0, 9/41]. Thus, the maxima are real
> analytic extrema.
> The function 'optimise()' again returns a wrong value.
>
> ex2 <- function(x) log((18/41) * x - 2 * x^2) + 16 * log(4 * x^2 -
> (36/41) * x
> + (9/41)) + 24 * log((23/82) + (18/41) * x - 2 * x^2)
> opt2 <- optimise(ex2, lower=0, upper=(9/41), maximum=T)
> x <- seq(0,9/41,0.001)
> plot(x,ex2(x),type="l")
> abline(v=opt2$maximum, col="red")
>
>
> The two functions are not only of theoretical interest but have to
> be optimised
> as likelihood functions when dealing with SNP data.
>
> I observed this problem on Windows XP and Open Suse 10.0. Other
> platforms have
> to be tested.
>
>
> Some informations:
>
> Open Suse 10.0
>
> platform i686-pc-linux-gnu
> arch i686
> os linux-gnu
> system i686, linux-gnu
> status
> major 2
> minor 4.0
> year 2006
> month 10
> day 03
> svn rev 39566
> language R
> version.string R version 2.4.0 (2006-10-03)
>
> Windows XP:
>
> platform i386-pc-mingw32
> arch i386
> os mingw32
> system i386, mingw32
> status
> major 2
> minor 4.0
> year 2006
> month 10
> day 03
> svn rev 39566
> language R
> version.string R version 2.4.0 (2006-10-03)
>
>
>
> With best regards,
>
> Karsten Krug
> Student at University of Leipzig
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm