st: Re: generalized nonlinear models

["Joseph Coveney" <jcoveney@bigplanet.com>]

David C Airey wrote:

Thanks again for the additional suggestions.

What I was after was "generalized nonlinear models". -regress- is to -glm- as
-nl- is to generalized nonlinear models. SAS offers the poisson distribution in
the ML based PROC NLMIXED, which also offer mixed models.

I found an example to emulate: slide 36 to 40 in
<www.uiweb.uidaho.edu/ag/statprog/BillWSWS05.ppt>. Although my data are animal
(C. elegans) rather than plant, it looks just like my dose-response data, more
of less.

R gnm offers generalized linear models without random effects.

Finally, I'm told offline that AD Model Builder (http://admb-project.org) can be
used.

--------------------------------------------------------------------------------

Apologies for jumping-in late to this thread, but I don't understand what the 
hitch is.  Setting up these models in Stata is not any more complicated than 
programming them in SAS's PROC NLMIXED.  The do-file below illustrates this 
using an example taken after that in the PowerPoint slide show that Dave 
references (Slides 36 through 40).  Specifying and fitting the model takes only 
about four lines of actual unique code, the rest being boilerplate.

Joseph Coveney

version 11.1

clear *
set more off
set seed `=date("2010-11-14", "YMD")'

set obs 200
generate double dose = 0.125 * runiform()
generate double mu = 10 + (70 - 10) / ///
    (1 + exp(2 * ln(dose / 0.01)))
replace dose = 0 if runiform() < 0.05
replace mu = 70 if dose == 0
genpoisson response, mu(mu)

*
* BEGIN HERE
*

program define gnml4
    version 11.1
    args lnfj C D B I

    tempvar mu
    quietly {
		generate double `mu' = cond(dose == 0, `D', ///
            `C' + (`D' - `C') / (1 + ////
            exp(`B' * ln(dose / `I'))))
        replace `lnfj' = -`mu' + $ML_y1 * ln(`mu') - ///
            lngamma($ML_y1 + 1) // you can skip this term
    }
end

ml model lf gnml4 (C: response = ) (D:) (B:) (I:)
ml init /C = 70 /D = 10 /B = 0.8254 /I = 0.01
ml maximize

exit


*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/