Re: [OctDev] regress and rcoplot (partially written, please help)

[Paul K. <pki...@us...>]

On Mar 7, 2006, at 8:51 PM, William Poetra Yoga Hadisoeseno wrote:

> On 2/21/06, William Poetra Yoga Hadisoeseno <wil...@gm...> 
> wrote:
>> On 2/20/06, Paul Kienzle <pki...@us...> wrote:
>>>
>>> Okay, I added a test of the Longley data set from the NIST 
>>> Statistical
>>> Reference Database.
>>>
>>> The performance is pretty bad (norm of answer-certified value is 
>>> about
>>> 10^-7), but that's a bit better than the competition.
>>>
>>> Someone with more R experience than me will have to see if it 
>>> performs
>>> any better.
>>>
>>> Also, if you can convert bint to standard error values that would be
>>> convenient.
>>>
>>>         - Paul
>>>
>>> ----
>>>
>>> % Longley data from the NIST Statistical Reference Dataset
>>> Z = [  60323    83.0   234289   2356     1590    107608  1947
>>>         61122    88.5   259426   2325     1456    108632  1948
>>>         60171    88.2   258054   3682     1616    109773  1949
>>>         61187    89.5   284599   3351     1650    110929  1950
>>>         63221    96.2   328975   2099     3099    112075  1951
>>>         63639    98.1   346999   1932     3594    113270  1952
>>>         64989    99.0   365385   1870     3547    115094  1953
>>>         63761   100.0   363112   3578     3350    116219  1954
>>>         66019   101.2   397469   2904     3048    117388  1955
>>>         67857   104.6   419180   2822     2857    118734  1956
>>>         68169   108.4   442769   2936     2798    120445  1957
>>>         66513   110.8   444546   4681     2637    121950  1958
>>>         68655   112.6   482704   3813     2552    123366  1959
>>>         69564   114.2   502601   3931     2514    125368  1960
>>>         69331   115.7   518173   4806     2572    127852  1961
>>>         70551   116.9   554894   4007     2827    130081  1962 ];
>>> % Results certified by NIST using 100 digit arithmetic
>>> % b and standard error in b
>>> V = [  -3482258.63459582         890420.383607373
>>>          15.0618722713733         84.9149257747669
>>>         -0.358191792925910E-01    0.334910077722432E-01
>>>         -2.02022980381683         0.488399681651699
>>>         -1.03322686717359         0.214274163161675
>>>         -0.511041056535807E-01    0.226073200069370
>>>          1829.15146461355         455.478499142212 ];
>>> Rsq = 0.995479004577296;
>>> F = 330.285339234588;
>>> y = Z(:,1); X = [ones(rows(Z),1), Z(:,2:end)];
>>> [b, bint, r, rint, stats] = regress (y, X, 0.05);
>>> assert(b, V(:,1), 6e-7)
>>> assert(stats(1), Rsq, 7e-14)
>>> assert(stats(2), F, 5e-9)
>>>
>>>
>>
>> Yes, you can convert bint to the standard error values with (0.05 is
>> alpha, 9 is the degree of freedom):
>> ((bint(:, 1) - bint(:, 2)) / 2) / tinv (0.05 / 2, 9)
>> but the accuracy will be affected. (It's not so accurate anyway...)
>>
>
> I still don't know how to improve the accuracy...
> I've added an assertion for the standard error values, but I've got to
> set the tolerance to 7.6e-3 (that's not good). Also, I think one
> decimal point for the tolerance would be good, what do you think?

I generally use only one digit for the tolerance.

> By the way, in the octave-forge version (in the repository), the
> functions tinv and fcdf are replaced by t_inv and f_cdf. Why?

Because I tinv and fcdf were not available in my version of Octave.
These can change once we abandon support for 2.1.x


- Paul