[R] Problem with NLSstClosestX; and suggested fix
Keith Jewell
k.jewell at campden.co.uk
Wed Oct 14 18:11:58 CEST 2009
Problem is demonstrated with this code, intended to find the approximate 'x'
at which the 'y' is midway between the left and right asymptotes. This
particular data set returns NA, which is a bit silly!
--------------
sXY <- structure(list(x = c(0, 24, 27, 48, 51, 72, 75, 96, 99), y =
c(4.98227,
6.38021, 6.90309, 7.77815, 7.64345, 7.23045, 7.27875, 7.11394,
6.95424)), .Names = c("x", "y"), row.names = c(NA, 9L), class =
c("sortedXyData",
"data.frame"))
a <- NLSstLfAsymptote(sXY)
d <- NLSstRtAsymptote(sXY)
NLSstClosestX(sXY, (a+d)/2)
------------------------
I think the problem arises when the target y value is exactly equal to one
of the y values in sXY and can be fixed by trapping that situation thus:
--------------------
NLSstClosestX.sortedXyData <- function (xy, yval)
{
deviations <- xy$y - yval
if (any(deviations==0)) xy$x[match(0, deviations)] else { # new line
inserted
if (any(deviations <= 0)) {
dev1 <- max(deviations[deviations <= 0])
lim1 <- xy$x[match(dev1, deviations)]
if (all(deviations <= 0)) {
return(lim1)
}
}
if (any(deviations >= 0)) {
dev2 <- min(deviations[deviations >= 0])
lim2 <- xy$x[match(dev2, deviations)]
if (all(deviations >= 0)) {
return(lim2)
}
}
dev1 <- abs(dev1)
dev2 <- abs(dev2)
lim1 + (lim2 - lim1) * dev1/(dev1 + dev2)
} # new line inserted
}
-------------------
Comments/corrections welcome.
Keith Jewell
=====================================================
> sessionInfo()
R version 2.9.2 (2009-08-24)
i386-pc-mingw32
locale:
LC_COLLATE=English_United Kingdom.1252;LC_CTYPE=English_United
Kingdom.1252;LC_MONETARY=English_United
Kingdom.1252;LC_NUMERIC=C;LC_TIME=English_United Kingdom.1252
attached base packages:
[1] stats graphics grDevices datasets tcltk utils methods
base
other attached packages:
[1] nlme_3.1-93 xlsReadWrite_1.3.3 svSocket_0.9-43 svMisc_0.9-48
TinnR_1.0.3 R2HTML_1.59-1
[7] Hmisc_3.6-1
loaded via a namespace (and not attached):
[1] cluster_1.12.0 grid_2.9.2 lattice_0.17-25 stats4_2.9.2
tools_2.9.2 VGAM_0.7-9
More information about the R-help
mailing list