[R] stats::power.t.test error
Martin Maechler
m@ech|er @end|ng |rom @t@t@m@th@ethz@ch
Tue Oct 15 17:59:15 CEST 2019
>>>>> Bert Gunter
>>>>> on Tue, 15 Oct 2019 07:41:35 -0700 writes:
> "...plausible sample sizes i.e. integers."
> ??
> f(...) = function that returns a real.
> ceiling(f(...)) = function that returns an integer.
> The problem is the "plausible" part.
Actually, power.t.test() does not return an integer for 'n'
typically in any case.
I found that it's actually quite easy to power.t.test() do the
diddling for us and return a number between 1 and 2.
What you do with that number is your decision, but formally it
solves the root finding problem :
> (ptt1 <- power.t.test(delta = 0.6, sd=0.00001, power=0.9 , sig.level=0.05))
Two-sample t test power calculation
n = 1.004283
delta = 0.6
sd = 1e-05
sig.level = 0.05
power = 0.9
alternative = two.sided
NOTE: n is number in *each* group
--------------------
As the change is small, and I see that Witold has good reasons
to prefer this to wrapping everything in try(.) or (better) tryCatch(.),
I propose to commit the change after a bit more testing.
Martin
> Cheers,
> Bert
> Bert Gunter
> "The trouble with having an open mind is that people keep coming along and
> sticking things into it."
> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
> On Tue, Oct 15, 2019 at 3:11 AM Witold E Wolski <wewolski using gmail.com> wrote:
>> Dear Peter,
>>
>> Yes, It is a technical issue and a matter of diddling around. And I
>> agree with your comment regarding the 2 observations.
>> I have several thousands variance estimates for which I need to
>> compute the sample sizes automatically. Using try statements is
>> typically the last thing I would like to resort too.
>> Is there an alternative implementation of power.t.test on CRAN which
>> could the diddling for me and return plausible sample sizes i.e.
>> integers.
>>
>> best regards
>> Witek
>>
>> On Fri, 4 Oct 2019 at 16:28, peter dalgaard <pdalgd using gmail.com> wrote:
>> >
>> > This is mainly a technical issue with uniroot trying to go outside of
>> its interval: (2, 1e7)
>> >
>> > It is fairly easy to find an approximate solution by diddling a little
>> by hand:
>> >
>> > > power.t.test(delta = 0.5849625, sd=0.01, n=1.04, sig.level=0.05)$power
>> > [1] 0.8023375
>> >
>> > Notice, however, that 1.04 observations in each group makes no sense at
>> all. In order to actually do a t-test you need at least 2 observations per
>> group (since we assume equal group sizes) or you have no variance estimate.
>> Already at sd=0.1, you are crossing the n=2 border, so for any smaller sd,
>> you will just get higher power with n=2. (Also, anything with single-digit
>> degrees of freedom for variance is probably expecting rather much regarding
>> to Gaussian distribution of your data.)
>> >
>> > -pd
>> >
>> > > On 4 Oct 2019, at 14:30 , Witold E Wolski <wewolski using gmail.com> wrote:
>> > >
>> > > Hi,
>> > >
>> > > power.t.test works for some range of input parameters but fails
>> otherwise.
>> > >
>> > >> power.t.test(delta = 0.5849625, sd=0.1, power=0.8, sig.level=0.05)$n
>> > > [1] 1.971668
>> > >> power.t.test(delta = 0.5849625, sd=0.05, power=0.8, sig.level=0.05)$n
>> > > [1] 1.620328
>> > >> power.t.test(delta = 0.5849625, sd=0.01, power=0.8, sig.level=0.05)$n
>> > > Error in uniroot(function(n) eval(p.body) - power, c(2, 1e+07), tol =
>> tol, :
>> > > did not succeed extending the interval endpoints for f(lower) *
>> f(upper) <= 0
>> > > In addition: Warning message:
>> > > In qt(sig.level/tside, nu, lower.tail = FALSE) : NaNs produced
>> > >
>> > > I guessing that sd is very small compared with delta, hence the
>> > > problem. But what are allowed values (ratios) of delta and sd?
>> > >
>> > > Best
>> > > Witek
>> > >
>> > > --
>> > > Witold Eryk Wolski
>> > >
>> >
>> > --
>> > Peter Dalgaard, Professor,
>> > Center for Statistics, Copenhagen Business School
>> > Solbjerg Plads 3, 2000 Frederiksberg, Denmark
>> > Phone: (+45)38153501
>> > Office: A 4.23
>> > Email: pd.mes using cbs.dk Priv: PDalgd using gmail.com
>>
>> --
>> Witold Eryk Wolski
More information about the R-help
mailing list