It worked. I re-written the equation as: ``` rutledge_param <- function(p, x, y) ( (p$M / ( 1 + exp(-(x-p$m)/p$s)) ) + p$B ) - y ``` and used Desmos to estimate the slope, so: ``` estim <- nls.lm(par = list(m = halfCycle, s = 2.77, M = MaxFluo, B = high[1]), fn = rutledge_param, x = 1:45, y = high) summary(estim) R <- rutledge(list(half_fluorescence = 27.1102, slope = 2.7680, max_fluorescence = 11839.7745, back_fluorescence = -138.8615) , 1:45) points(1:45, R, type="l", col="red") ```
Thanks On Tue, Mar 16, 2021 at 8:29 AM Luigi Marongiu <marongiu.lu...@gmail.com> wrote: > > Just an update: > I tried with desmos and the fitting looks good. Desmos calculated the > parameters as: > Fmax = 11839.8 > Chalf = 27.1102 (with matches with my estimate of 27 cycles) > k = 2.76798 > Fb = -138.864 > I forced R to accept the right parameters using a single named list > and re-written the formula (it was a bit unclear in the paper): > ``` > rutledge <- function(p, x) { > m = p$half_fluorescence > s = p$slope > M = p$max_fluorescence > B = p$back_fluorescence > y = (M / (1+exp( -((x-m)/s) )) ) + B > return(y) > } > ``` > but when I apply it I get a funny graph: > ``` > desmos <- rutledge(list(half_fluorescence = 27.1102, slope = 2.76798, > max_fluorescence = 11839.8, back_fluorescence > = -138.864) , high) > ``` > > On Mon, Mar 15, 2021 at 7:39 AM Luigi Marongiu <marongiu.lu...@gmail.com> > wrote: > > > > Hello, > > the negative data comes from the machine. Probably I should use raw > > data directly, although in the paper this requirement is not reported. > > The p$x was a typo. Now I corrected it and I got this error: > > ``` > > > > > rutledge_param <- function(p, x, y) ((p$M / (1 + exp(-1*(x-p$m)/p$s))) + > > > p$B) - y > > > estim <- nls.lm(par = list(m = halfFluo, s = slopes, M = MaxFluo, B = > > > high[1]), > > + fn = rutledge_param, x = 1:45, y = high) > > Error in dimnames(x) <- dn : > > length of 'dimnames' [2] not equal to array extent > > ``` > > Probably because 'slopes' is a vector instead of a scalar. Since the > > slope is changing, I don't think is right to use a scalar, but I tried > > and I got: > > ``` > > > estim <- nls.lm(par = list(m = halfFluo, s = 1, M = MaxFluo, B = high[1]), > > + fn = rutledge_param, x = 1:45, y = high) > > > estim > > Nonlinear regression via the Levenberg-Marquardt algorithm > > parameter estimates: 6010.94, 1, 12021.88, 4700.49288888889 > > residual sum-of-squares: 1.14e+09 > > reason terminated: Relative error in the sum of squares is at most `ftol'. > > ``` > > The values reported are the same I used at the beginning apart from > > the last (the background parameter) which is 4700 instead of zero. If > > I plug it, I get an L shaped plot that is worse than that at the > > beginning: > > ``` > > after = init = rutledge(halfFluo, 1, MaxFluo, 4700.49288888889, high) > > points(1:45, after, type="l", col="blue") > > ``` > > What did I get wrong here? > > Thanks > > > > On Sun, Mar 14, 2021 at 8:05 PM Bill Dunlap <williamwdun...@gmail.com> > > wrote: > > > > > > > rutledge_param <- function(p, x, y) ((p$M / (1 + > > > > exp(-1*(p$x-p$m)/p$s))) + p$B) - y > > > > > > Did you mean that p$x to be just x? As is, this returns numeric(0) > > > for the p that nls.lm gives it because p$x is NULL and NULL-aNumber is > > > numeric(). > > > > > > -Bill > > > > > > On Sun, Mar 14, 2021 at 9:46 AM Luigi Marongiu <marongiu.lu...@gmail.com> > > > wrote: > > > > > > > > Hello, > > > > I would like to use the Rutledge equation > > > > (https://pubmed.ncbi.nlm.nih.gov/15601990/) to model PCR data. The > > > > equation is: > > > > Fc = Fmax / (1+exp(-(C-Chalf)/k)) + Fb > > > > I defined the equation and another that subtracts the values from the > > > > expectations. I used minpack.lm to get the parameters, but I got an > > > > error: > > > > ``` > > > > > > > > > library("minpack.lm") > > > > > h <- c(120.64, 66.14, 34.87, 27.11, 8.87, -5.8, 4.52, -7.16, -17.39, > > > > + -14.29, -20.26, -14.99, -21.05, -20.64, -8.03, -21.56, -1.28, > > > > 15.01, > > > > + 75.26, 191.76, 455.09, 985.96, 1825.59, 2908.08, 3993.18, > > > > 5059.94, > > > > + 6071.93, 6986.32, 7796.01, 8502.25, 9111.46, 9638.01, 10077.19, > > > > + 10452.02, 10751.81, 11017.49, 11240.37, 11427.47, 11570.07, > > > > 11684.96, > > > > + 11781.77, 11863.35, 11927.44, 11980.81, 12021.88, 12058.35, > > > > 12100.63, > > > > + 12133.57, 12148.89, 12137.09) > > > > > high <- h[1:45] > > > > > MaxFluo <- max(high) > > > > > halfFluo <- MaxFluo/2 > > > > > halfCycle = 27 > > > > > find_slope <- function(X, Y) { > > > > + Slope <- c(0) > > > > + for (i in 2:length(X)) { > > > > + delta_x <- X[i] - X[i-1] > > > > + delta_y <- Y[i] - Y[i-1] > > > > + Slope[i] <- delta_y/delta_x > > > > + } > > > > + return(Slope) > > > > + } > > > > > slopes <- find_slope(1:45, high) > > > > > > > > > > rutledge <- function(m, s, M, B, x) { > > > > + divisor = 1 + exp(-1* ((x-m)/s) ) > > > > + y = (M/divisor) + B > > > > + return(y) > > > > + } > > > > > rutledge_param <- function(p, x, y) ((p$M / (1 + > > > > > exp(-1*(p$x-p$m)/p$s))) + p$B) - y > > > > > > > > > > > > > > > init = rutledge(halfFluo, slopes, MaxFluo, 0, high) > > > > > points(1:45, init, type="l", col="red") > > > > > estim <- nls.lm(par = list(m = halfFluo, s = slopes, M = MaxFluo, B = > > > > > high[1]), > > > > + fn = rutledge_param, x = 1:45, y = high) > > > > Error in nls.lm(par = list(m = halfFluo, s = slopes, M = MaxFluo, B = > > > > high[1]), : > > > > evaluation of fn function returns non-sensible value! > > > > ``` > > > > > > > > Where could the error be? > > > > > > > > > > > > -- > > > > Best regards, > > > > Luigi > > > > > > > > ______________________________________________ > > > > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > > > > 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. > > > > > > > > -- > > Best regards, > > Luigi > > > > -- > Best regards, > Luigi -- Best regards, Luigi ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.