| |
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]
Re: st: Iteratively re-weighted least squares (IRLS)
From |
Joseph Coveney <[email protected]> |
To |
Statalist <[email protected]> |
Subject |
Re: st: Iteratively re-weighted least squares (IRLS) |
Date |
Sat, 07 Apr 2007 16:04:06 +0900 |
Just to add to the earlier reply, the -lnlsq(0)- option actually fits what
pharmacokineticists call the "exponential" variance model, which is (or at
least it used to be) often used by them to approximate a "proportional"
variance model such as yours.
Perhaps a better approach would be to explicitly and directly fit a
"proportional" error model, such with as the Nelder-Mead simplex algorithm.
The simplex algorithm is implemented in Stata as a user-written
command, -amoeba-. I've illustrated its use in your problem in the do-file
below, and would recommend it over -nl , lnlsq()-, although the
approximation with the latter doesn't look too shabby in the illustrative
case and wouldn't need resampling techniques for confidence limits or
hypothesis testing.
Joseph Coveney
clear
set more off
set seed `=date("2007-04-07", "ymd")'
set obs 500
// C = (a + b*X^g)*e
local a 2
local b 3
local g 0.75
generate double X = uniform()
generate double e = 0.15 * invnormal(uniform())
generate double C = (`a' + `b' * X^`g') * (1 + e)
*
capture program drop proerrmod
program define proerrmod
version 9.2
tempvar intermediate residual
quietly {
generate double `intermediate' = ///
`1'[1,1] + `1'[1,2] * $xvar^(`1'[1,3])
generate double `residual' = ($yvar / `intermediate') - 1
replace `residual' = `residual' * `residual'
}
summarize `residual', meanonly
scalar define `2' = -r(sum)
end
*
matrix define init = (1, 1, 1)
scalar define obj = .
matrix define B = (., ., .)
matrix colnames init = a:_cons b:_cons g:_cons
global yvar C
global xvar X
amoeba proerrmod init obj B
exit
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/