Take the 2-minute tour ×
Stack Overflow is a question and answer site for professional and enthusiast programmers. It's 100% free, no registration required.

A measurement shows a signal that is formed like a square root function with offset and a factor. How can I find the coefficients and plot the raw data and the fitted curve in one plot?

require(ggplot2)
require(nlmrt)  # may be this will help later..?

# generate simulated measurement data

time <- seq(-10,20,0.2)

signal <- sqrt(time + 2) # generate sqrt signal; includes some NA
signal[is.na(signal)] <- 0 # set all NA to zero
signal <- signal + rnorm(length(time)) * 0.1  # add noise

df <- data.frame(x=time, y=signal)

# find coefficiants for y ~ b * sqrt(x - a)
# no idea how...

# plot raw data and fitted curve in one ggplot diagram

ggplot()+
    geom_point(data=df, aes(x=x, y=y))

enter image description here

share|improve this question
1  
Just a hint, you probably want to fit y ~ b * sqrt(x - a), as opposed to a + b*sqrt(c*x). –  Hong Ooi Jun 14 '13 at 17:33
1  
If you can't specify the cutoff point a priori (i.e. if it's a fitted parameter), you're probably going to have to fit this separately (it wouldn't fit into a linear model framework, and there are particularly tricky aspects to fitting cutoffs since the goodness-of-fit is flat over the intervals between data points). If you can specify the cutoff then I think you can construct a dummy variable and do this ... –  Ben Bolker Jun 14 '13 at 18:21
    
As first step I can define a constant for cutoff manually. Later I would like to find this automatic. But this might be topic of another question then. –  Jonas Stein Jun 14 '13 at 18:28
add comment

1 Answer 1

up vote 4 down vote accepted

Provided you know where the cutpoint is and that the value before the cutpoint is zero:

sfun <- function(x,brk,a=1) {
    ifelse(x<brk,0,suppressWarnings(a*sqrt(x-brk)))
}

(suppressWarnings() is there because ifelse evaluates both the if and the else cases for all values of x, and we don't want warnings about taking the square root of negative numbers)

Test (not shown):

curve(sfun(x,1,1),from=0,to=10) ## test (not shown)

Simulate some data:

x <- seq(0,10,length=101)
set.seed(1)
y <- rnorm(length(x),sfun(x,1,1),sd=0.25)
DF <- data.frame(x,y)

Since all we need to figure out is how the square root function is scaled, we can do this with a regression through the origin (take out the -1 if you want to allow the value below the cutpoint to be non-zero):

library("ggplot2")
theme_set(theme_bw())
ggplot(DF,aes(x,y))+geom_point()+
    geom_smooth(method="lm",
                formula=y~sfun(x,brk=1)-1)
ggsave("truncsqrt.png")

enter image description here

share|improve this answer
    
oops, I didn't even notice that you'd provided your own simulation code. Anyone who wants can feel free to edit in the OP's example. –  Ben Bolker Jun 14 '13 at 21:35
    
Perfect! How is the width of the grey band defined? And can I even get the coefficients as numbers? –  Jonas Stein Jun 14 '13 at 22:52
1  
the grey band is the 95% confidence interval. If you want the coefficients, just do lm(y~sfun(x,brk=1)-1,data=DF) –  Ben Bolker Jun 14 '13 at 23:28
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.