Hi Folks,
I note that, while the "chisq" functions
dchisq(x, df, ncp=0, log = FALSE)
pchisq(q, df, ncp=0, lower.tail = TRUE, log.p = FALSE)
qchisq(p, df, ncp=0, lower.tail = TRUE, log.p = FALSE)
rchisq(n, df, ncp=0)
all have a slot for the non-centrality parameter "ncp", of
the functions for the t and F distributions:
dt(x, df, log = FALSE)
pt(q, df, ncp=0, lower.tail = TRUE, log.p = FALSE)
qt(p, df, lower.tail = TRUE, log.p = FALSE)
rt(n, df)
df(x, df1, df2, log = FALSE)
pf(q, df1, df2, ncp=0, lower.tail = TRUE, log.p = FALSE)
qf(p, df1, df2, lower.tail = TRUE, log.p = FALSE)
rf(n, df1, df2)
only the CDF functions 'pt' and 'pf' allow this parameter to
be set. (If you try in the others, you get the message
"unused argument(s) (ncp ...)").
Why is this? Being able to set it would be just as useful ...
Thanks, and best wishes to all,
Ted.
[and apologies if I should have been able to read the answer
for myself, somewhere ... ]
--------------------------------------------------------------------
E-Mail: (Ted Harding) <Ted.Harding at nessie.mcc.ac.uk>
Fax-to-email: +44 (0)870 167 1972
Date: 17-Oct-02 Time: 16:24:31
------------------------------ XFMail ------------------------------
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
(Ted Harding) <Ted.Harding at nessie.mcc.ac.uk> writes:> only the CDF functions 'pt' and 'pf' allow this parameter to > be set. (If you try in the others, you get the message > "unused argument(s) (ncp ...)"). > > Why is this? Being able to set it would be just as useful ...We don't have any references on how to calculate them! (Except for the brute-force approaches of numeric differentiation, root finding, and transformation of uniform distributions.) -- O__ ---- Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907 -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Ted Harding says:> -----Original Message----- > From: Ted.Harding at nessie.mcc.ac.uk > [mailto:Ted.Harding at nessie.mcc.ac.uk] > Sent: Friday, October 18, 2002 2:14 AM > To: Peter Dalgaard BSA > Cc: r-help at stat.math.ethz.ch > Subject: Re: [R] Non-central distributions > > Thanks, Peter! (Must try to give this some thought ...).[WNV] The density functions for the t and F distributions are in fact quite easy and only require hypergeometric functions in addition to standard things. These could be useful anyway for all sorts of things. As far as I know, percentage points of the non-central distributions are not much used, but what would be very useful would be to have the percentage points (with respect to the non-centrality parameter) of the distribution function G(delta) = 1-P(X^2, n, delta), (i.e. you take the upper tail area as defining a distribution function in delta. Such a distributon has a finite probability at the origin, of course. These are the quantities you need, for example, for things like sample size determination and power calculations. Random numbers from the non-central distributions are easy enough to generate, of course, using the central ones. Again, I'm not sure just how much slick versions of them would be useful, though.> Anyway, in this respect R is still ahead of S-Plus, which > doesn't seem to carry ANY non-centrality as standard! > (Except possibly obscurely tucked away in some add-on library).[WNV] Tsk tsk, Ted. They are there for pf and pchisq, at least. Bill Venables.> Ted. > > On 17-Oct-02 Peter Dalgaard BSA wrote: > > (Ted Harding) <Ted.Harding at nessie.mcc.ac.uk> writes: > > > >> only the CDF functions 'pt' and 'pf' allow this parameter to > >> be set. (If you try in the others, you get the message > >> "unused argument(s) (ncp ...)"). > >> > >> Why is this? Being able to set it would be just as useful ... > > > > We don't have any references on how to calculate them! (Except for the > > brute-force approaches of numeric differentiation, root finding, and > > transformation of uniform distributions.) > > -------------------------------------------------------------------- > E-Mail: (Ted Harding) <Ted.Harding at nessie.mcc.ac.uk> > Fax-to-email: +44 (0)870 167 1972 > Date: 17-Oct-02 Time: 17:13:30 > ------------------------------ XFMail ------------------------------ > -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-. > -.-.- > r-help mailing list -- Read > http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html > Send "info", "help", or "[un]subscribe" > (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch > _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._. > _._._-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
> all have a slot for the non-centrality parameter "ncp", of > the functions for the t and F distributions: > > dt(x, df, log = FALSE) > pt(q, df, ncp=0, lower.tail = TRUE, log.p = FALSE) > qt(p, df, lower.tail = TRUE, log.p = FALSE) > rt(n, df)Here is a function for calculating the density of a non-central t. dtnc <- function(x, df, ncp=0) { y <- -ncp*x/sqrt(df+x^2) a <- (-y + sqrt(y^2 + 4*df)) / 2 Hhmy <- 1/gamma(df+1) * a^df * exp(-0.5*(a+y)^2) * sqrt(2*pi*a^2/(df+a^2)) * (1 - 3*df/(4*(df+a^2)^2) + 5*df^2/(6*(df+a^2)^3)) gamma(df + 1) / (2^((df-1)/2) * gamma(df/2) * sqrt(pi*df)) * exp(-0.5*df*ncp^2/(df+x^2)) * (df/(df+x^2))^((df+1)/2) * Hhmy } This is an approximation based on Resnikoff & Lieberman (1957). But it's quite accurate. The exact pdf can either be expressed as an infinite sum or requires the evaluation of an integral. I tried to implement the latter using integrate(), but the results were very erratic. I checked the accuracy of the approximation in various ways: Comparing dtnc to dt in case ncp = 0 results in a maximum error of .0046 at the point x = 0 and df = 1. For df = 10, that error is already down to .000021. I also compared pt with integrate(dtnc, lower=-Inf, upper=x, df=df, ncp=ncp) for a wide range of df's and ncp's and found that the results closely matched to at least three decimal places even for small df. For 12 df, differences only showed up in the 5th decimal. Obviously, use of integrate() might contribute a little bit to the discrepancies. Also, I checked whether dtnc integrates to 1 with integrate(dtnc, lower=-Inf, upper=Inf, df=df, ncp=ncp) for a wide range of df's and ncp's. Differences from 1 only showed up in the second decimal places for df = 1, third decimal place for df = 2, and the fourth decimal place for df = 4, and so on ... So, it looks like the approximation provides pretty accurate results. -- Wolfgang Viechtbauer -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._