Um texto embutido e sem conjunto de caracteres especificado associado... Nome: n?o dispon?vel Url: https://stat.ethz.ch/pipermail/r-help/attachments/20071202/9577be0a/attachment.pl
Res: fitting "power model" in nls()
4 messages · Milton Cezar Ribeiro, Gabor Grothendieck
OK. Since the model is linear except for A lets use brute force to repeatedly evaluate the sum of squares for values of A between -2 and 2 proceeding in steps of .01 solving the other parameters using lm. That will give us better starting values and we should be able to use nls on that.
x <- seq(-2, 2, .01) ss <- sapply(x, function(A) sum(resid(lm(richness ~ I(area^A)))^2)) plot(ss ~ x) x[which.min(ss)]
[1] -0.45
model.lm <- lm(richness ~ I(area^-0.45)) # use starting values based on lm and A = -0.45 st <- c(Const = coef(model.lm)[[1]], B = coef(model.lm)[[2]], A = x[which.min(ss)]) nls(richness ~ Const+B*(area^A), st = st)
Nonlinear regression model model: richness ~ Const + B * (area^A) data: parent.frame() Const B A 33.9289 -33.4595 -0.4464 residual sum-of-squares: 8751 Number of iterations to convergence: 2 Achieved convergence tolerance: 3.368e-06 Note that our A value is suspiciously close to A = -0.5 and sqrt(area) is length so I wonder if there is an argument based on units of measurement that might support a model of the form: richness = Const + B / sqrt(area)
On Dec 2, 2007 3:39 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear Gabor, Thank you for your reply. In fact I am ajusting several models at same time, like linear, log-linear, log-log, piecewise etc. One of the models are the power model. I really need to fit a power model because it one of the hypothesis which have been suggested on literature. In addition, there are other variables which are beeing tested as explanatory. Kind regards, miltinho ----- Mensagem original ---- De: Gabor Grothendieck <ggrothendieck at gmail.com> Para: Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> Cc: R-help <r-help at stat.math.ethz.ch> Enviadas: Domingo, 2 de Dezembro de 2007 17:28:23 Assunto: Re: [R] fitting "power model" in nls() Is that really the model we want? When we have problems sometimes its just a sign that the model is not very good in the first place. plot(richness ~ area) shows most of the points crowded the left and just a few points out to the right. This does not seem like a very good pattern for model fitting. plot(richness ~ log(area)) plot(log(richness) ~ log(area)) both look nicer. On Dec 2, 2007 2:08 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear all, I am still fighting against my "power model". I tryed several times to use nls() but I can?t run it. I am sending my variables and also the model which I would like to fit. As you can see, this "power model" is not the best model to be fit, but I
really need also to fit it.
The model which I would like to fit is Richness = B*(Area^A)
richness<-c(44,36,31,39,38,26,37,33,34,48,25,22,44,5,9,13,17,15,21,10,16,22,13,20,9,15,14,21,23,23,32,29,20,
26,31,4,20,25,24,32,23,33,34,23,28,30,10,29,40,10,8,12,13,14,56,47,44,37,27,17,32,31,26,23,31,34,
37,32,26,37,28,38,35,27,34,35,32,27,22,23,13,28,13,22,45,33,46,37,21,28,38,21,18,21,18,24,18,23,22,
38,40,52,31,38,15,21)
area<-c(26.22,20.45,128.68,117.24,19.61,295.21,31.83,30.36,13.57,60.47,205.30,40.21,
7.99,1.18,5.40,13.37,4.51,36.61,7.56,10.30,7.29,9.54,6.93,12.60, 2.43,18.89,15.03,14.49,28.46,36.03,38.52,45.16,58.27,67.13,92.33,1.17,
29.52,84.38,87.57,109.08,72.28,66.15,142.27,76.41,105.76,73.47,1.71,305.75,
325.78,3.71,6.48,19.26,3.69,6.27,1689.67,95.23,13.47,8.60,96.00,436.97,
472.78,441.01,467.24,1169.11,1309.10,1905.16,135.92,438.25,526.68,88.88,31.43,21.22,
640.88,14.09,28.91,103.38,178.99,120.76,161.15,137.38,158.31,179.36,214.36,187.05,
140.92,258.42,85.86,47.70,44.09,18.04,127.84,1694.32,34.27,75.19,54.39,79.88,
63.84,82.24,88.23,202.66,148.93,641.76,20.45,145.31,27.52,30.70)
plot(richness~area)
I tryed to fit the following model:
m1<-nls(richness ~ Const+B*(area^A))
Thanks a lot,
miltinho
Brazil.
para armazenamento!
[[alternative HTML version deleted]]
______________________________________________ R-help at 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. ________________________________ Abra sua conta no Yahoo! Mail, o ?nico sem limite de espa?o para armazenamento!
Also the fitted values satisfy Const = -B = 33 (approximately) so we could try:
plot(richness ~ area) nls(richness ~ C * (1 - 1/sqrt(area)), start = c(C = 33))
Nonlinear regression model
model: richness ~ C * (1 - 1/sqrt(area))
data: parent.frame()
C
32.85
residual sum-of-squares: 8764
Number of iterations to convergence: 1
Achieved convergence tolerance: 5.595e-10
simple.nls <- .Last.value points(fitted(simple.nls) ~ area, pch = "+", col = "red")
On Dec 2, 2007 4:06 PM, Gabor Grothendieck <ggrothendieck at gmail.com> wrote:
OK. Since the model is linear except for A lets use brute force to repeatedly evaluate the sum of squares for values of A between -2 and 2 proceeding in steps of .01 solving the other parameters using lm. That will give us better starting values and we should be able to use nls on that.
x <- seq(-2, 2, .01) ss <- sapply(x, function(A) sum(resid(lm(richness ~ I(area^A)))^2)) plot(ss ~ x) x[which.min(ss)]
[1] -0.45
model.lm <- lm(richness ~ I(area^-0.45)) # use starting values based on lm and A = -0.45 st <- c(Const = coef(model.lm)[[1]], B = coef(model.lm)[[2]], A = x[which.min(ss)]) nls(richness ~ Const+B*(area^A), st = st)
Nonlinear regression model model: richness ~ Const + B * (area^A) data: parent.frame() Const B A 33.9289 -33.4595 -0.4464 residual sum-of-squares: 8751 Number of iterations to convergence: 2 Achieved convergence tolerance: 3.368e-06 Note that our A value is suspiciously close to A = -0.5 and sqrt(area) is length so I wonder if there is an argument based on units of measurement that might support a model of the form: richness = Const + B / sqrt(area) On Dec 2, 2007 3:39 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear Gabor, Thank you for your reply. In fact I am ajusting several models at same time, like linear, log-linear, log-log, piecewise etc. One of the models are the power model. I really need to fit a power model because it one of the hypothesis which have been suggested on literature. In addition, there are other variables which are beeing tested as explanatory. Kind regards, miltinho ----- Mensagem original ---- De: Gabor Grothendieck <ggrothendieck at gmail.com> Para: Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> Cc: R-help <r-help at stat.math.ethz.ch> Enviadas: Domingo, 2 de Dezembro de 2007 17:28:23 Assunto: Re: [R] fitting "power model" in nls() Is that really the model we want? When we have problems sometimes its just a sign that the model is not very good in the first place. plot(richness ~ area) shows most of the points crowded the left and just a few points out to the right. This does not seem like a very good pattern for model fitting. plot(richness ~ log(area)) plot(log(richness) ~ log(area)) both look nicer. On Dec 2, 2007 2:08 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear all, I am still fighting against my "power model". I tryed several times to use nls() but I can?t run it. I am sending my variables and also the model which I would like to fit. As you can see, this "power model" is not the best model to be fit, but I
really need also to fit it.
The model which I would like to fit is Richness = B*(Area^A)
richness<-c(44,36,31,39,38,26,37,33,34,48,25,22,44,5,9,13,17,15,21,10,16,22,13,20,9,15,14,21,23,23,32,29,20,
26,31,4,20,25,24,32,23,33,34,23,28,30,10,29,40,10,8,12,13,14,56,47,44,37,27,17,32,31,26,23,31,34,
37,32,26,37,28,38,35,27,34,35,32,27,22,23,13,28,13,22,45,33,46,37,21,28,38,21,18,21,18,24,18,23,22,
38,40,52,31,38,15,21)
area<-c(26.22,20.45,128.68,117.24,19.61,295.21,31.83,30.36,13.57,60.47,205.30,40.21,
7.99,1.18,5.40,13.37,4.51,36.61,7.56,10.30,7.29,9.54,6.93,12.60, 2.43,18.89,15.03,14.49,28.46,36.03,38.52,45.16,58.27,67.13,92.33,1.17,
29.52,84.38,87.57,109.08,72.28,66.15,142.27,76.41,105.76,73.47,1.71,305.75,
325.78,3.71,6.48,19.26,3.69,6.27,1689.67,95.23,13.47,8.60,96.00,436.97,
472.78,441.01,467.24,1169.11,1309.10,1905.16,135.92,438.25,526.68,88.88,31.43,21.22,
640.88,14.09,28.91,103.38,178.99,120.76,161.15,137.38,158.31,179.36,214.36,187.05,
140.92,258.42,85.86,47.70,44.09,18.04,127.84,1694.32,34.27,75.19,54.39,79.88,
63.84,82.24,88.23,202.66,148.93,641.76,20.45,145.31,27.52,30.70)
plot(richness~area)
I tryed to fit the following model:
m1<-nls(richness ~ Const+B*(area^A))
Thanks a lot,
miltinho
Brazil.
para armazenamento!
[[alternative HTML version deleted]]
______________________________________________ R-help at 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. ________________________________ Abra sua conta no Yahoo! Mail, o ?nico sem limite de espa?o para armazenamento!
I played around with this a bit more and noticed that the "plinear" algorithm of nls converged using nearly every starting value I tried. In fact A = 0 was the only starting value that I could find that did not converge. Note that with "plinear" you only specify the starting values for non-linear parameters, in this case A, while the unnamed linear parameters are implied as coefficients of the columns of the matrix defined in the rhs.
nls(richness ~ cbind(1, area^A), start = c(A = 1), algorithm = "plinear")
Nonlinear regression model
model: richness ~ cbind(1, area^A)
data: parent.frame()
A .lin1 .lin2
-0.4464 33.9290 -33.4595
residual sum-of-squares: 8751
Number of iterations to convergence: 6
Achieved convergence tolerance: 4.968e-07
On Dec 2, 2007 4:06 PM, Gabor Grothendieck <ggrothendieck at gmail.com> wrote:
OK. Since the model is linear except for A lets use brute force to repeatedly evaluate the sum of squares for values of A between -2 and 2 proceeding in steps of .01 solving the other parameters using lm. That will give us better starting values and we should be able to use nls on that.
x <- seq(-2, 2, .01) ss <- sapply(x, function(A) sum(resid(lm(richness ~ I(area^A)))^2)) plot(ss ~ x) x[which.min(ss)]
[1] -0.45
model.lm <- lm(richness ~ I(area^-0.45)) # use starting values based on lm and A = -0.45 st <- c(Const = coef(model.lm)[[1]], B = coef(model.lm)[[2]], A = x[which.min(ss)]) nls(richness ~ Const+B*(area^A), st = st)
Nonlinear regression model model: richness ~ Const + B * (area^A) data: parent.frame() Const B A 33.9289 -33.4595 -0.4464 residual sum-of-squares: 8751 Number of iterations to convergence: 2 Achieved convergence tolerance: 3.368e-06 Note that our A value is suspiciously close to A = -0.5 and sqrt(area) is length so I wonder if there is an argument based on units of measurement that might support a model of the form: richness = Const + B / sqrt(area) On Dec 2, 2007 3:39 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear Gabor, Thank you for your reply. In fact I am ajusting several models at same time, like linear, log-linear, log-log, piecewise etc. One of the models are the power model. I really need to fit a power model because it one of the hypothesis which have been suggested on literature. In addition, there are other variables which are beeing tested as explanatory. Kind regards, miltinho ----- Mensagem original ---- De: Gabor Grothendieck <ggrothendieck at gmail.com> Para: Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> Cc: R-help <r-help at stat.math.ethz.ch> Enviadas: Domingo, 2 de Dezembro de 2007 17:28:23 Assunto: Re: [R] fitting "power model" in nls() Is that really the model we want? When we have problems sometimes its just a sign that the model is not very good in the first place. plot(richness ~ area) shows most of the points crowded the left and just a few points out to the right. This does not seem like a very good pattern for model fitting. plot(richness ~ log(area)) plot(log(richness) ~ log(area)) both look nicer. On Dec 2, 2007 2:08 PM, Milton Cezar Ribeiro <milton_ruser at yahoo.com.br> wrote:
Dear all, I am still fighting against my "power model". I tryed several times to use nls() but I can?t run it. I am sending my variables and also the model which I would like to fit. As you can see, this "power model" is not the best model to be fit, but I
really need also to fit it.
The model which I would like to fit is Richness = B*(Area^A)
richness<-c(44,36,31,39,38,26,37,33,34,48,25,22,44,5,9,13,17,15,21,10,16,22,13,20,9,15,14,21,23,23,32,29,20,
26,31,4,20,25,24,32,23,33,34,23,28,30,10,29,40,10,8,12,13,14,56,47,44,37,27,17,32,31,26,23,31,34,
37,32,26,37,28,38,35,27,34,35,32,27,22,23,13,28,13,22,45,33,46,37,21,28,38,21,18,21,18,24,18,23,22,
38,40,52,31,38,15,21)
area<-c(26.22,20.45,128.68,117.24,19.61,295.21,31.83,30.36,13.57,60.47,205.30,40.21,
7.99,1.18,5.40,13.37,4.51,36.61,7.56,10.30,7.29,9.54,6.93,12.60, 2.43,18.89,15.03,14.49,28.46,36.03,38.52,45.16,58.27,67.13,92.33,1.17,
29.52,84.38,87.57,109.08,72.28,66.15,142.27,76.41,105.76,73.47,1.71,305.75,
325.78,3.71,6.48,19.26,3.69,6.27,1689.67,95.23,13.47,8.60,96.00,436.97,
472.78,441.01,467.24,1169.11,1309.10,1905.16,135.92,438.25,526.68,88.88,31.43,21.22,
640.88,14.09,28.91,103.38,178.99,120.76,161.15,137.38,158.31,179.36,214.36,187.05,
140.92,258.42,85.86,47.70,44.09,18.04,127.84,1694.32,34.27,75.19,54.39,79.88,
63.84,82.24,88.23,202.66,148.93,641.76,20.45,145.31,27.52,30.70)
plot(richness~area)
I tryed to fit the following model:
m1<-nls(richness ~ Const+B*(area^A))
Thanks a lot,
miltinho
Brazil.
para armazenamento!
[[alternative HTML version deleted]]
______________________________________________ R-help at 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. ________________________________ Abra sua conta no Yahoo! Mail, o ?nico sem limite de espa?o para armazenamento!