s-news
[Top] [All Lists]

Problems porting R selfStart SSweibull to S

To: "'s-news@lists.biostat.wustl.edu'" <s-news@lists.biostat.wustl.edu>
Subject: Problems porting R selfStart SSweibull to S
From: "Yang, Richard" <dyang@NRCan.gc.ca>
Date: Tue, 20 Nov 2001 18:18:06 -0500
S++ers;

I experienced difficulties porting the SSweibull function in R to S (Splus 6
or 2000). The code looks like the following:

        >wb <- deriv( ~ Asym - Drop*exp(-exp(lrc)*x^pwr), c("Asym", "Drop",
"lrc", "pwr"),
    function(x, Asym, Drop, lrc, pwr) {})

        >wbInit <-             ## from Dr. Doug Bates
        function (mCall, data, LHS) 
        {
    xy <- sortedXyData(mCall[["x"]], LHS, data)
    if (nrow(xy) < 5) {
        stop("Too few distinct input values to fit the Weibull growth
model")
    }
    if (any(xy[["x"]] < 0)) {
        error("All x values must be non-negative to fit the Weibull growth
model")
    }
    Rasym <- NLSstRtAsymptote(xy)
    Lasym <- NLSstLfAsymptote(xy)
    pars <- coef(lm(log(-log((Rasym - y)/(Rasym - Lasym))) ~ 
        log(x), data = xy, subset = x > 0))
    val <- coef(nls(y ~ cbind(1, -exp(-exp(lrc) * x^pwr)), data = xy, 
        alg = "plinear", start = c(lrc = pars[[1]], pwr = pars[[2]])))[c(3, 
        4, 1, 2)]
    names(val) <- mCall[c("Asym", "Drop", "lrc", "pwr")]
    val
        }

        >SSwb <- selfStart(wb, initial = wbInit)

Following the creation of SSwb selfStart function, a test with the Orange
data resulted in an error:

        > getInitial(circumference ~ SSwb(age, Asym, Drop, lrc, pwr), data =

        Orange)
        Error: Object "Rasym" not found

It sounded like the NLSstRtsymptote() was not working. However,

        >Orange.sortAvg <- sortedXyData("age", "circumference", Orange)

        >NLSstRtAsymptote(Orange.sortAvg)
        >[1] 193.9
        > NLSstLfAsymptote(Orange.sortAvg)
        >[1] 12.9

The SSweibull function works fine with the Orange dataset in R, but the same
error occurred in both Splus 6 and Splus 2000. 

        What did I miss? Please advise. 

        TIA

Richard


<Prev in Thread] Current Thread [Next in Thread>