From svetosch@gmx.net Sun Dec  6 18:22:52 2009
From: Sven Schreiber <svetosch@gmx.net>
To: gretl-devel@gretlml.univpm.it
Subject:
 Re: [Gretl-devel] Goodness of Fit test to normal distribution and QQ plot
Date: Mon, 07 Dec 2009 00:22:27 +0100
Message-ID: <4B1C3CB3.6090402@gmx.net>
In-Reply-To: Pine.A41.4.58.0911302143220.2568566@f1n11.sp2net.wfu.edu
MIME-Version: 1.0
Content-Type: multipart/mixed; boundary="===============0465242305800350481=="

--===============0465242305800350481==
Content-Type: text/plain; charset="utf-8"
Content-Transfer-Encoding: 7bit

Allin Cottrell schrieb:
> On Mon, 30 Nov 2009, Sven Schreiber wrote:
> 
>> Costia schrieb:
>>> 2. I haven't found a Q-Q plot test in GRETL, (like in SAS e.g.). I think
>>> it would be a good thing to add it.
>> I agree. It's not difficult to use a quick-and-dirty script to get it
>> done, but properly handling missing values etc. needs more care.
> 
> I haven't messed much with Q-Q plots.  Can anyone confirm if the
> following script is correct?  If so, then I guess we could do it
> as a built-in at some point.

sorry for not answering earlier...

> 
> <script>
> nulldata 500
> scalar n = $nobs
> # artificial series to test
> genr y = sort(normal())
> series Qemp, Qnorm
> scalar p, N, nl, nh
> 
> loop i=1..n -q
>   p = i/(n+1)

>   N = (n + 1) * p - 1;

given your def. of p, wouldn't that be equiv. to N = i-1 ?

>   nl = floor(N);

accordingly redundant? (nl = i-1)

>   nh = ceil(N);

ditto?

>   if (nh == 0 || nh == n)
>     Qemp[i] = NA

don't understand the reason for this: here it seems the first item in
the empirical quantiles will always be set to NA?

>   else
>     Qemp[i] = quantile(y, p)
>   endif
>   Qnorm[i] = critical(N, 1-p)

or rather Qnorm[i] = 1-critical(N,p)? (could make a difference for
discrete distributions, not for continuous as here)

> endloop
> 

and the 45-degree line is needed for comparison

altogether I think something like this also would do the trick (bugs
probably included, even if it's only two lines...):

matrix mp = transp(seq(1,n))/n	
# (col vector of cumulative relative frequencies for n obs)
matrix mQtheory = invcdf(N,mp)	
# BTW, the help about invcdf() says P(X<x), but shouldn't it be P(X<=x)?

and your y (the sorted input series) can already be plotted against
mQtheory to give the QQ plot (if I'm not missing something right now
which may well be the case since I'm a bit tired). You might argue that
this doesn't really deal with ties, but from the various possible
definitions I've seen of quantiles per se and how to deal with ties I
think it's actually legitimate.

Why not just copy the convention and maybe the code from R?

good night,
sven

--===============0465242305800350481==--