CODEWITHSUNDEEP
rcurve-fittingpredictionnon-linear-regression
Vint
Vint

Reputation: 499

Predicting data from a power curve manually

I have a series of data I have fit a power curve to, and I use the predict function in R to allow me predict y values based on additional x values.

set.seed(1485)
len <- 24
x <- runif(len)
y <- x^3 + rnorm(len, 0, 0.06)
ds <- data.frame(x = x, y = y)


mydata=data.frame(x,y)
z <- nls(y ~ a * x^b, data = mydata, start = list(a=1, b=1))
#z is same as M!

power <- round(summary(z)$coefficients[1], 3)
power.se <- round(summary(z)$coefficients[2], 3)
plot(y ~ x, main = "Fitted power model", sub = "Blue: fit; green: known")
s <- seq(0, 1, length = 100)
lines(s, s^3, lty = 2, col = "green")
lines(s, predict(z, list(x = s)), lty = 1, col = "blue")
text(0, 0.5, paste("y =x^ (", power, " +/- ", power.se,")", sep = ""), pos = 4)

Instead of using the predict function here, how could I manually calculate estimated y values based on additional x values based on this power function. If this were just a simple linear regression, I would calculate the slope and y intercept and calculate my y values by 

y= mx + b

Is there a similar equation I can use from the output of z that will allow me to estimate y values from additional x values? 

> z
Nonlinear regression model
  model: y ~ a * x^b
   data: mydata
    a     b 
1.026 3.201 
 residual sum-of-squares: 0.07525

Number of iterations to convergence: 5 
Achieved convergence tolerance: 5.162e-06

Upvotes: 0

Views: 479

Answers (2)

G. Grothendieck
G. Grothendieck

Reputation: 270160

Here are some ways.

1) with This evaluates the formula with respect to the coefficients: 

x <- 1:2 # input

with(as.list(coef(z)), a * x^b)
## [1] 1.026125 9.437504

2) attach We could also use attach although it is generally frowned upon:

attach(as.list(coef(z)))
a * x^b
## [1] 1.026125 9.437504

3) explicit Explicit definition:

a <- coef(z)[["a"]]; b <- coef(z)[["b"]]
a * x^b
## [1] 1.026125 9.437504

4) eval This one extracts the formula from z so that we don't have to specify it again. formula(z)[[3]] is the right hand side of the formula used to produce z. Use of eval is sometimes frowned upon but this does avoid the redundant specification of the formula.

eval(formula(z)[[3]], as.list(coef(z)))
## [1] 1.026125 9.437504

Upvotes: 0

Ian Wesley
Ian Wesley

Reputation: 3624

You would do it the same way except you use the power equation you modeled. You can access the parameters the model calculated using z$m$getPars()

Here is a simple example to illustrate:

predict(z, list(x = 1))

Results in: 1.026125

Which equals the results of 

z$m$getPars()["a"] * 1 ^ z$m$getPars()["b"]

Which is equivalet to y = a * x^b

Upvotes: 1

Related Questions