it's that easy, eh?
thank you for your input, it's much appreciated.
i will have a look at nlmer as well.
markus
On Thu, Feb 16, 2012 at 23:17, David Winsemius [via R] <
ml-node+s789695n4395634...@n4.nabble.com> wrote:
>
> On Feb 16, 2012, at 11:43 AM, protoplast wrote:
>
> > hello mailing list!
> > i still consider myself an R beginner, so please bear with me if my
> > questions seems strange.
> >
> > i'm in the field of biology, and have done consecutive hydraulic
> > conductivity measurements in three parallels ("Sample"), resulting
> > in three
> > sets of conductivity values ("PLC" for percent loss of conductivity,
> > relative to 100%) at multiple pressures ("MPa").
> >
> > ---
> > Sample MPa PLC
> > 1 -0.3498324 0.000000
> > 1 -1.2414770 15.207821
> > 1 -1.7993249 23.819995
> > 1 -3.0162866 33.598570
> > 1 -3.5184321 46.376933
> > 1 -3.9899791 67.532226
> > 1 -4.2731145 89.735541
> > 1 -4.7597244 99.836239
> > 2 -0.2754036 0.000000
> > 2 -1.2912619 12.476132
> > 2 -1.5128974 13.543273
> > ...
> > ---
> >
> > since each sample is a statistical unit, i have fitted each sample-
> > subset to
> > a sigmoid curve:
> >
> > ---
> > plot(
> > NA,
> > NA,
> > main="",
> > xlim=c(-20,0),
> > ylim=c(0,100),
> > xlab = "water potential [MPa]",
> > ylab = "percent loss of conductivity [%]",
> > xaxp = c(0,-20,4),
> > yaxp = c(0,100,5),
> > tck = 0.02,
> > mgp = c(1.5,0.1,0),
> > )
> >
> > for(i in 1:3){
> > x <- subset(curvedata,Group == i)$MPa
> > y <- subset(curvedata,Group == i)$PLC
> > name <- subset(curvedata,Group == i)$Sample
> > points(x,y)
> >
> > vlc <- nls(y ~ 100/(1+exp(a*(x-b))), start=c(a=1, b=-3),
> > data=list(x,y))
> >
> > curve(100/(1+exp(coef(vlc)[1]*(x-coef(vlc)[2]))), col=1, add = TRUE)
> >
> > Rsquared <- 1 - var(residuals(vlc))/var(y)
> >
> > summarizeall[i ,"Run"] <- i
> > summarizeall[i ,"Sample"] <- name[1]
> > summarizeall[i ,"a"] <- coef(vlc)[1]
> > summarizeall[i ,"b"] <- coef(vlc)[2]
> > summarizeall[i ,"R2"] <- Rsquared
> >
> > listnow <- data.frame(list(Run = c(i),Sample = c(name[1]), a =
> > c(coef(vlc)[1]), b = c(coef(vlc)[2]), R2 = c(Rsquared)))
> > print(listnow)
> >
> > i <- i+1
> > }
> > ---
> >
> > ...and get three slightly different curves with three different
> > estimatinos
> > of fit (r², Rsquared).
> >
> > ---
> >> summarizeall
> > Sample a b R2
> > 1 1 1.388352 -3.277755 0.9379886
> > 2 2 1.800007 -3.363075 0.9327164
> > 3 3 1.736857 -2.743972 0.9882998
> >
> >> average
> > Var n a b R2
> > 1 Mean 3 1.6417389 -3.1282673 NA
> > 2 SE . 0.1279981 0.1937197 NA
> > ---
> >
> > by averaging parameters a and b of the curve, i create a "mean
> > curve" that
> > is added to the plot (red curve in the attached image).
> >
> > http://r.789695.n4.nabble.com/file/n4394609/conductivity-curve.gif
> >
> > ---
> > meana <- average[1,"a"]
> > meanb <- average[1,"b"]
> > curve(), col=2, lwd=2, add = TRUE)
> > ---
> >
> > and now here's my problem:
> > i'd like to calculate R squared for all points on that mean curve.
> > since i have to average the curve parameters, i loose the curve's
> > residuals
> > that are stored in my variable vlc (the result of the nls function)
> > for
> > every sample.
> > just fitting one curve to all the data points is not good enough.
>
> So just calculate them?
> # pseudo-code: residual= actual - predicted
>
> gresid <- curvedata$PLC - 100/(1+exp(meana*(curvedata$MPa-meanb))
>
>
> If you are convinced that your formula for R^2 makes sense and this
> practice is generally accepted in your domain, then you can apply it
> across the whole dataset. I would have thought that a single
> regression model built with nlmer might have been more statistically
> sound. (But this is a bit outside my domain of comfort for giving
> advice.)
>
> >
> > an extensive google search over several days hasn't gotten me
> > anywhere, but
> > maybe someone here can help me?
> >
> > is there an efficient way to calculate r squared for a predefined
> > function
> > with "unrelated" data points?
> > (unrelated as in "not used directly for fitting")
> >
> >
> > thanks in advance
> > markus
>
> David Winsemius, MD
> West Hartford, CT
>
> ______________________________________________
> [hidden email] <http://user/SendEmail.jtp?type=node&node=4395634&i=0>mailing
> list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>
> ------------------------------
> If you reply to this email, your message will be added to the discussion
> below:
>
> http://r.789695.n4.nabble.com/how-to-get-r-squared-for-a-predefined-curve-or-function-with-other-data-points-tp4394609p4395634.html
> To unsubscribe from how to get r-squared for a predefined curve or
> function with "other" data points, click
> here<http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_code&node=4394609&code=bm9sZi5tYXJrdXNAZ21haWwuY29tfDQzOTQ2MDl8LTE0MDUwOTcyOTk=>
> .
> NAML<http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=macro_viewer&id=instant_html%21nabble%3Aemail.naml&base=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespace&breadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml>
>
--
View this message in context:
http://r.789695.n4.nabble.com/how-to-get-r-squared-for-a-predefined-curve-or-function-with-other-data-points-tp4394609p4396720.html
Sent from the R help mailing list archive at Nabble.com.
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
