From luke at stat.uiowa.edu  Fri Dec  9 17:05:26 2005
From: luke at stat.uiowa.edu (Luke Tierney)
Date: Fri, 9 Dec 2005 10:05:26 -0600 (CST)
Subject: [Rd] qt for df < 1
In-Reply-To: <x2k6efbgva.fsf@viggo.kubism.ku.dk>
References: <FD963C8B-71E1-4322-AE17-98303ED6B0D8@uiuc.edu>
	<x2k6efbgva.fsf@viggo.kubism.ku.dk>
Message-ID: <Pine.LNX.4.63.0512091003150.20896@nokomis.stat.uiowa.edu>

On Thu, 8 Dec 2005, Peter Dalgaard wrote:

> roger koenker <rkoenker at uiuc.edu> writes:
>
>> I was experimenting yesterday with a binomial make.link option
>> for estimating student t binary response models, tentatively
>> called gossit, and I noticed eventually that the R qt function doesn't
>> like df < 1.  Vaguely recalling that Splus didn't seem to mind such
>> weirdness,  I checked on our soon to be defunct Splus6.2 and
>> sure enough, it produced plausible answers instead of R's NA's.
>> Of course, I have no way of judging the quality of these answers,
>> but I'm curious about whether someone has already looked into
>> this can of worms.
>
> Well the help page has:
>
> For 'qt' only values of at least one are currently supported.
>
> and someone must have written that...
>
> R does have pt for df < 1, so a temporary fix using uniroot() seems
> doable.
>
>

Something like

     qqt<-function(p,df) sign(p-0.5)*sqrt(qf(1-2*pmin(p,1-p),1,df))

seems to do reasonably, at least in terms of consistency with pt, down
to 0.2 or mayby 0.1 df based on

     f<-function(d, n = 101) {
 	x<-seq(0, 1, len = n)
 	max(abs(pt(qqt(x, d), d) - x))
     }
     plot(function(d) log10(sapply(d, f)), .01,1)

luke

-- 
Luke Tierney
Chair, Statistics and Actuarial Science
Ralph E. Wareham Professor of Mathematical Sciences
University of Iowa                  Phone:             319-335-3386
Department of Statistics and        Fax:               319-335-3017
    Actuarial Science
241 Schaeffer Hall                  email:      luke at stat.uiowa.edu
Iowa City, IA 52242                 WWW:  http://www.stat.uiowa.edu



[Rd] qt for df < 1
Luke Tierney