Hello Jim,

Thanks very much for providing the article. I wrote to Wolfram, and they replied back. They thanked me and said they will work on the issue (contrary to some comments I have read on the web that they ignore bug reports). Of course, they have not fixed it yet.

I will read the article in detail and will probably use it. I found another article (which I have in my work computer, but which I can send if you are interested). That article provided two good approximations which I could then use as a starting point in FindRoot to get fairly good results. We need these values at work and what is published in tables is not always enough or is too inconvenient.

I looked into the issue myself and it is fairly technical. Some features of this distributions are known exactly, for others only bounds are known, but those could be used as a first guess or bracketing step in a solution method. I am not surprised that someone devoted a fundamental article just to the handling of this problem. I am a bit surprised that Wolfram let it slip. They could have been rushing through that development. Hopefully they will flush it out; they got smart people there.

Thanks a lot,

OL.

Excellent! That is a concise example that shows the issue. (One can format Mathematica code using the "Code Sample" icon. Just hover over the editing icons just below the "Reply to this discussion" text to see approximately what each means.) Here is your code doing so:

a = InverseCDF[NormalDistribution[0, 1], 0.95];
f[n_, cl_] := InverseCDF[NoncentralStudentTDistribution[n - 1, Sqrt[n]*a], cl]
fInv[n_, k_] := CDF[NoncentralStudentTDistribution[n - 1, Sqrt[n]*a], k]
Column[(AbsoluteTiming[{#, fInv[#, f[#, 0.95]]}] &) /@Range[415, 420]]

(Note that there were missing underscores in the function definitions that I've added in.)

I get the same result as you with Mathematica 10.1 (Windows 7):

{218.647, {415, 0.95}},
{326.461, {416, 0.95}},
{222.949, {417, 0.95}},
{339.589, {418, 0.95}},
{80.7641, {419, 0.50088}},
{116.356, {420, 0.50338}}

The good news is that it doesn't take as long to get a wrong answer. (Sorry, couldn't help myself.)

So...Is this a bug or just pushing the function too far? If the latter, should there not be a warning?

Thanks again. You are right that it takes less than a tenth of a second. It is also interesting that the message prints when calculating the CDF rather than the inverse (first time I see it that way).

You are also right that my original posting was long.

I still think that if I can improve this with a simple bisection code, then there should be room for further improvement.

Thanks,

OL.

It takes less than a tenth of a second. I wasn't trying to simplify things. While it was clear from your write-up as to what the general problem was, if I might be so bold, your message rambled a bit and I was just trying to get something more specific to shoot at. You'll just be more likely to get help with concise set of statements (using the "Code Sample" formatting).

Hello Jim,

Thanks for replying to my post.

Yes, it captures the essence of the issue. But your result is accurate (even if it printed a warning).

You simplified it slightly by making the non-centrality parameter equal to zero, which basically means you are using the central student T distribution (even if you invoke the non-central one). In my case I need the non-centrality parameter to be given by the formula in my original post, which makes it proportional to the Sqrt of the sample size. The central distribution is symmetric about the peak and peaks at zero. The non-central one is not only shifted away from zero, but also becomes asymmetric about the peak.

I notice is that version 10 seems to be behaving similarly (you got a warning and probably a somewhat long computation time). However, you got an accurate result. In my case, I have noticed that the result was inaccurate for sample sizes starting around 419.

Thanks,

OL.