Hi Phil,
thanks for reviewing the multiple linear regression implementations and
setting up the R/NIST data tests. I finally got around to installing R
and can now run them too.
Phil Steitz wrote:
While clear and elegant from a matrix algebra standpoint, the "nailve"
implementation in OLSMultipleLinearRegression has bad numerical
qualities. It is well known that solving the normal equations directly
does not give good numerics. I just added some tests to actually verify
parameter values, using the classic "Longly" dataset, for which NIST
provides certified statistics. This is a "hard" design matrix. R was
able to get to within 1E-8 of the certified parameter values.
OLSMultipleLinearRegression can only get 1E-1.
The OLS implementation has been added as a simple by-product of the GLS
case - which is the main one I have needed for hypothesis testing - as
it came "for free" with unitary covariance.
True - the emphasis was on clarity and formulaic simplicity. And also
following the old Donald Knuth maxim "optimization is the root of all
evil". But it seems like there is a need for refinement of the
implementation - the devil raised his head :-)
We have talked in the past about providing an implementation based on QR
decomposition. Anyone up for using the QR decomposition that we now
have to do this? I really think we need to do it (or something else to
improve numerics) before releasing this class. I will get to it
eventually, but am a little pegged at the moment. I will review and
apply patches if someone is willing to do the implementation. I can
also explain here or offline how the R tests and NIST datasets work, as
these are useful in validating code.
I'd be happy to improve the impl. I'm getting my head around R and
NIST, but perhaps a chat offline would not hurt!
Another thing that we should think about before releasing any of this
stuff is the completeness of the API. Many standard regression
statistics are missing. If we are going to stick with the Interface /
Implementation setup, we need to get the right stuff into the
interface. It is also awkward to have to insert "1"'s in the design
matrix to get an intercept term computed. This is convenient for
implementation, but awkward for users. A more natural setup (IMHO)
would be to expose a "noIntercept" or "hasIntercept" property for the
model.
No problem with adding other statistics - let's just decide on what is
the stardard regression API.
And finally, how do you see the no/hasIntercept model working?
Cheers
---------------------------------------------------------------------
To unsubscribe, e-mail: [EMAIL PROTECTED]
For additional commands, e-mail: [EMAIL PROTECTED]