[R] Intersection for two curves
Muhammad Rahiz
muhammad.rahiz at ouce.ox.ac.uk
Sat Apr 24 02:06:46 CEST 2010
Thanks David & Peter,
The locator() works but not practical as I have to repeat the process
many times.
Does the code works on linear regression only? When i tried to find the
intersection at a non-linear curve, i get the following error
Error in optimize(f = function(x) abs(xyf(ds) - n), c(0, 13)) :
invalid function value in 'optimize'
I'd like to know what the error means and how I can correct it.
I have my sample data and code as follows;
* s1 cm
1 0 0.57223196
2 10 0.33110049
3 20 0.11163181
4 30 0.10242237
5 40 0.09254315
6 50 0.02739370
7 60 0.02567137
8 70 0.02492397
9 80 0.03206637
10 90 0.02487381
11 100 0.01747584
12 110 0.15977533
13 120 0.13317708
rm(list=ls())
x<- read.table("test.data.txt",header=TRUE)
ds <- x[,2] # distance
cr <- x[,3] # correlation values
plot(ds,cr)
n <- 1/2.71
abline(h=n)
fc <- function(x,a,b){a*exp(-b*x)} # where a & b are constants
fm <- nls(cr ~ fc(ds,a,b),start=c(a=1,b=0))
co <- coef(fm)
curve(fc(x,a=co[1],b=co[2]),add=TRUE,col="red",lwd=1.5)
int <- function(x) coef(fm)[1] + x*coef(fm)[2]
in1 <- optimize(f=function(x) c(0,120), abs(int(ds)-n))
abline(v=in1$minimize)
thanks,
Muhammad
On 04/23/2010 07:32 PM, Peter Ehlers wrote:
> On 2010-04-23 11:46, David Winsemius wrote:
>
>> On Apr 23, 2010, at 1:06 PM, Muhammad Rahiz wrote:
>>
>>
>>> Does anyone know of a method that I can get the intersection where the
>>> red and blue curves meet i.e. the value on the x-axis?
>>>
>>> x<- 1:10
>>> y<- 10:1
>>> plot(x,y)
>>> abline(lm(y~x),col="blue")
>>> abline(h=2.5,col="red")
>>>
>> Two ways :
>>
>> > xy<- lm(y~x)
>> > xyf<- function(x) coef(xy)[1] +x*coef(xy)[2]
>>
>> # absolute difference
>> > optimise(f=function(x) abs(xyf(x)-2.5), c(1,10) )
>> $minimum
>> [1] 8.49998
>>
>> $objective
>> (Intercept)
>> 1.932015e-05
>>
>> #N minimize squared difference
>> > optimise(f=function(x) (xyf(x)-2.5)^2, c(1,10) )
>> $minimum
>> [1] 8.5
>>
>> $objective
>> (Intercept)
>> 3.155444e-30
>>
>>
> Another (crude) way is to use locator(). I usually maximize
> the plot window for this.
>
>
More information about the R-help
mailing list