nls factor
A more compact way to code factors in nls is to use the syntax factor[]. Here's an example using a simplified version of Ravi's example: n <- 200 set.seed(123) x <- runif(n) a <- gl(n=2, k=n/2) # a two-level factor eps <- rnorm(n, sd=0.5) y <- as.numeric(a) * x^.5 + eps nls(y ~ a[]*x^b, start=list(a=c(1,1), b=0.25))
Wed, 2009-04-15 at 10:22 -0400, Ravi Varadhan wrote:
Hi,
Here is one-way to do it (the following code shows a simulation example):
n <- 200
set.seed(123)
x <- runif(n)
f <- gl(n=2, k=n/2) # a two-level factor
x1 <- x * (f == 1)
x2 <- x * (f == 2)
a <- c(rep(2, n/2), rep(5, n/2))
b <- 0.5
nsim <- 100
nls.coef <- matrix(NA, nsim, 3)
for (i in 1:nsim) {
eps <- rnorm(n, sd=0.5)
y <- a * x^b + eps
ans.nls <- try(nls(y ~ a1*x1^b + a2*x2^b, start=list(a1=1, a2=1, b=0.25)),
silent=TRUE)
if (class(ans.nls) != "try-error") nls.coef[i, ] <- coef(ans.nls)
}
apply(nls.coef, 2, summary)
Hope this helps,
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: rvaradhan at jhmi.edu
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
----------------------------------------------------------------------------
--------
-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
Behalf Of Manuel Gutierrez
Sent: Wednesday, April 15, 2009 8:27 AM
To: r-help at r-project.org
Subject: [R] nls factor
I want to fit the model y=a*x^b using nls; where "a" should be different for
each level of a factor.
What is the easiest way to fit it? Can i do it with nls?
I've looked the help pages and the MASS example in page 249 but the formula
is different and I don't know how to specify it for my model.
Thanks,
Manuel
[[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. ______________________________________________ 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.