stats::power.t.test error

classic Classic list List threaded Threaded
7 messages Options
Reply | Threaded
Open this post in threaded view
|

stats::power.t.test error

Witold E Wolski
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

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

R help mailing list-2
Think about this. What is the null hypothesis? What is the alternative?
What are their distributions? What is the probability that you get a value
from the alternative when the null hypothesis holds and vice versa? Then
think again about the relevance of your alternative hypothesis. You'll get
a better understanding of power calculation by doing such exercise.

ir. Thierry Onkelinx
Statisticus / Statistician

Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE AND
FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
[hidden email]
Havenlaan 88 bus 73, 1000 Brussel
www.inbo.be

///////////////////////////////////////////////////////////////////////////////////////////
To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey
///////////////////////////////////////////////////////////////////////////////////////////

<https://www.inbo.be>


Op vr 4 okt. 2019 om 14:30 schreef Witold E Wolski <[hidden email]>:

> 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
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

Peter Dalgaard-2
In reply to this post by Witold E Wolski
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 <[hidden email]> 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
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

--
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Office: A 4.23
Email: [hidden email]  Priv: [hidden email]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

Witold E Wolski
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 <[hidden email]> 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 <[hidden email]> 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
> >
> > ______________________________________________
> > [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
> --
> Peter Dalgaard, Professor,
> Center for Statistics, Copenhagen Business School
> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> Phone: (+45)38153501
> Office: A 4.23
> Email: [hidden email]  Priv: [hidden email]
>
>
>
>
>
>
>
>
>


--
Witold Eryk Wolski

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

Bert Gunter-2
"...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.

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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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
> > >
> > > ______________________________________________
> > > [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> > > https://stat.ethz.ch/mailman/listinfo/r-help
> > > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > > and provide commented, minimal, self-contained, reproducible code.
> >
> > --
> > Peter Dalgaard, Professor,
> > Center for Statistics, Copenhagen Business School
> > Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> > Phone: (+45)38153501
> > Office: A 4.23
> > Email: [hidden email]  Priv: [hidden email]
> >
> >
> >
> >
> >
> >
> >
> >
> >
>
>
> --
> Witold Eryk Wolski
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

Peter Dalgaard-2
In reply to this post by Witold E Wolski
You don't really want the diddling, since it gives meaningless values anyway...

For a pragmatic strategy, how about this?:

(a) calculate the power at n=2, if bigger than target power, done, else
(b) calculate n to reach target power, now guaranteed to have n > 2. Round upwards.

Peter D.

> On 15 Oct 2019, at 12:07 , Witold E Wolski <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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
>>>
>>> ______________________________________________
>>> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>
>> --
>> Peter Dalgaard, Professor,
>> Center for Statistics, Copenhagen Business School
>> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
>> Phone: (+45)38153501
>> Office: A 4.23
>> Email: [hidden email]  Priv: [hidden email]
>>
>>
>>
>>
>>
>>
>>
>>
>>
>
>
> --
> 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: [hidden email]  Priv: [hidden email]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: stats::power.t.test error

Martin Maechler
In reply to this post by Bert Gunter-2
>>>>> 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 <[hidden email]> 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 <[hidden email]> 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 <[hidden email]> 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: [hidden email]  Priv: [hidden email]
    >>
    >> --
    >> Witold Eryk Wolski

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.