Failure in predicting parameters
Locked
 
|
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 ______________________________________________ [hidden email] 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. |
Re: Failure in predicting parameters
|
|
Do the negative values in your data make any sense? Note that if Fb must be
>0, Fc must be also. But I have *not* examined your code/equations in detail, so feel free to ignore if this is irrelevant. Cheers, Bert Bert Gunter "The trouble with having an open mind is that people keep coming along and sticking things into it." -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip ) On Sun, Mar 14, 2021 at 9:46 AM Luigi Marongiu <[hidden email]> 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 > > ______________________________________________ > [hidden email] 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. > [[alternative HTML version deleted]] ______________________________________________ [hidden email] 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. |
Re: Failure in predicting parameters
|
|
In reply to this post by Luigi
> 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 <[hidden email]> 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 > > ______________________________________________ > [hidden email] 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. ______________________________________________ [hidden email] 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. |
|
|
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 <[hidden email]> 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 <[hidden email]> 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 > > > > ______________________________________________ > > [hidden email] 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 ______________________________________________ [hidden email] 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. |
|
|
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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 > > > > > > ______________________________________________ > > > [hidden email] 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 ______________________________________________ [hidden email] 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. |
|
|
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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 > > > > > > > > ______________________________________________ > > > > [hidden email] 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 ______________________________________________ [hidden email] 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. |
Re: Failure in predicting parameters
|
|
Hello,
Maybe a bit late but there is a contributed package [1] for quantitative PCR fitting non-linear models with the Levenberg-Marquardt algorithm. estim and vector R below are your model and your fitted values vector. The RMSE of this fit is smaller than your model's. Isn't this simpler? library(qpcR) df1 <- data.frame(Cycles = seq_along(high), high) fit <- pcrfit( data = df1, cyc = 1, fluo = 2 ) summary(fit) coef(estim) coef(fit) sqrt(sum(resid(estim)^2)) #[1] 1724.768 sqrt(sum(resid(fit)^2)) #[1] 1178.318 highpred <- predict(fit, newdata = df1) plot(1:45, high, type = "l", col = "red") points(1:45, R, col = "blue") points(1:45, highpred$Prediction, col = "cyan", pch = 3) [1] https://CRAN.R-project.org/package=qpcR Hope this helps, Rui Barradas Às 06:51 de 18/03/21, Luigi Marongiu escreveu: > 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 >>>>> >>>>> ______________________________________________ >>>>> [hidden email] 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 > > > ______________________________________________ [hidden email] 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. |
|
|
Thank you, I'll try it!
On Thu, Mar 18, 2021 at 9:46 PM Rui Barradas <[hidden email]> wrote: > > Hello, > > Maybe a bit late but there is a contributed package [1] for quantitative > PCR fitting non-linear models with the Levenberg-Marquardt algorithm. > > estim and vector R below are your model and your fitted values vector. > The RMSE of this fit is smaller than your model's. > > > Isn't this simpler? > > > library(qpcR) > > df1 <- data.frame(Cycles = seq_along(high), high) > > fit <- pcrfit( > data = df1, > cyc = 1, > fluo = 2 > ) > summary(fit) > > coef(estim) > coef(fit) > > > sqrt(sum(resid(estim)^2)) > #[1] 1724.768 > sqrt(sum(resid(fit)^2)) > #[1] 1178.318 > > > highpred <- predict(fit, newdata = df1) > > plot(1:45, high, type = "l", col = "red") > points(1:45, R, col = "blue") > points(1:45, highpred$Prediction, col = "cyan", pch = 3) > > > [1] https://CRAN.R-project.org/package=qpcR > > Hope this helps, > > Rui Barradas > > Às 06:51 de 18/03/21, Luigi Marongiu escreveu: > > 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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 > >>>>> > >>>>> ______________________________________________ > >>>>> [hidden email] 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 ______________________________________________ [hidden email] 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. |
«
Return to R help
|
1 view|%1 views
| Free forum by Nabble | Edit this page |