R › R help

calculate quantiles of a custom function

Classic

List

Threaded

13 messages Options

Gerhard

calculate quantiles of a custom function

Hi,

I guess that my problem has an obvious answer, but I have not been able to
find it.

Suppose I create a custom function, consisting of two beta-distributions:

myfunction <- function(x) {
dbeta(x,2,6) + dbeta(x,6,2)
}

How can I calculate the quantiles of myfunction?

I have not seen any continous function treated in the docs, and applying the
"quantile function" gives me an error (since it seems only to be defined on
lists and atoms).

Thank you in advance,

Gerhard

______________________________________________
[hidden email] mailing list
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.

Bert Gunter

Re: calculate quantiles of a custom function

Gerhard:

Strictly speaking, it's quantiles of a custom "distribution", not function.

There may be some way to handle your example easily, but, in general,
you would need to solve the resulting integral equation. This is hard
-- closed form solutions rarely exist; good approximations require
work.

So a standard approach is: simulate. Indeed, many simulation tricks
(under the rubric of "variance reduction") have been developed exactly
for such monte carlo integration. Consult a good reference or
knowledgeable person for details.

-- Bert

On Tue, Jan 3, 2012 at 4:24 AM, Gerhard <[hidden email]> wrote:

> Hi,
>
> I guess that my problem has an obvious answer, but I have not been able to
> find it.
>
> Suppose I create a custom function, consisting of two beta-distributions:
>
> myfunction <- function(x) {
> dbeta(x,2,6) + dbeta(x,6,2)
> }
>
> How can I calculate the quantiles of myfunction?
>
> I have not seen any continous function treated in the docs, and applying the
> "quantile function" gives me an error (since it seems only to be defined on
> lists and atoms).
>
> Thank you in advance,
>
> Gerhard
>
> ______________________________________________
> [hidden email] mailing list
> 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.

--

Bert Gunter
Genentech Nonclinical Biostatistics

Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm

______________________________________________
[hidden email] mailing list
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.

David Winsemius

Re: calculate quantiles of a custom function

In reply to this post by Gerhard

On Jan 3, 2012, at 7:24 AM, Gerhard wrote:

> Hi,
>
> I guess that my problem has an obvious answer, but I have not been
> able to
> find it.
>
> Suppose I create a custom function, consisting of two beta-
> distributions:
>
> myfunction <- function(x) {
> dbeta(x,2,6) + dbeta(x,6,2)
> }
>

Given the symmetry of the beta function, I suspect that dividing by
two would make that a real distribution.

> How can I calculate the quantiles of myfunction?

You could build a CDF using 'integrate'.

>
> I have not seen any continous function treated in the docs, and
> applying the
> "quantile function" gives me an error (since it seems only to be
> defined on
> lists and atoms).

There is also a family of packages, all of whose names begin with
"distr". They provide a variety of facilities for developing new
distributions, both discrete and continuous.

--

David Winsemius, MD
West Hartford, CT

______________________________________________
[hidden email] mailing list
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.

Gerhard

Re: calculate quantiles of a custom function

Am Dienstag, 3. Januar 2012, 11:05:11 schrieben Sie:

>
> The quick way is to look at the structure with 'str':
>
> str(integrate(myfunction,0,.9))
> List of 5
> $ value : num 1.85
> $ abs.error : num 2.05e-14
> $ subdivisions: int 1
> $ message : chr "OK"
> $ call : language integrate(f = myfunction, lower = 0, upper
> = 0.9)
> - attr(*, "class")= chr "integrate"
>
> The longer way would be to look at the help)integrate) page where the
> "Value" section would have told you that a list with named components
>
> was returned:
> > integrate(myfunction,0,.9)$value
>
> [1] 1.850299
>
> In this instance you could also probably use 'cumsum' (which is just
> doing the numerical integration by hand
>
> plot(cumsum( myfunction(seq(0,1,by=0.001)))/
> sum(myfunction(seq(0,1,by=0.001)) ) )

Thanks a lot; now I now where I need to look for help before bothering you
again.

Best,

Gerhard

>
> > I need to hand on the value to another function, and then I get an
> > error
> > message, because the output of integrate is not a numerical argument.
> >
> >>> I have not seen any continous function treated in the docs, and
> >>> applying the
> >>> "quantile function" gives me an error (since it seems only to be
> >>> defined on
> >>> lists and atoms).
> >>
> >> There is also a family of packages, all of whose names begin with
> >> "distr". They provide a variety of facilities for developing new
> >> distributions, both discrete and continuous.
> >
> > That looks quite frightening. But if it turns out to be necessary.
>
> But a good way to improve you knowledge of statistical concepts.
>
> > Thank you,
> >
> > Gerhard
>
> David Winsemius, MD
> West Hartford, CT

______________________________________________
[hidden email] mailing list
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.

VictorDelgado

Re: calculate quantiles of a custom function

In reply to this post by Gerhard

Gerhard wrote

Suppose I create a custom function, consisting of two beta-distributions:

myfunction <- function(x) {
dbeta(x,2,6) + dbeta(x,6,2)
}

How can I calculate the quantiles of myfunction?

Thank you in advance,

Gerhard

Gehard, if do you want to know the quantiles of the new distribution created by "myfunction". Maybe you can also do:

x <- seq(0,1,.01) # insert your 'x'
q <- myfunction(x)
# And:
quantile(x)

0% 25% 50% 75% 100%
0.000000 1.476177 2.045389 2.581226 2.817425

# This gives the sample quantiles. You can also look foward to simulations (like Bert Gunter had suggested) to know better the properties of distributions quantiles obtained after 'myfunction'.

Victor Delgado cedeplar.ufmg.br P.H.D. student UFOP assistant professor

VictorDelgado

Re: calculate quantiles of a custom function

VictorDelgado wrote

quantile(x)

Correcting to

quantile(q)

Victor Delgado cedeplar.ufmg.br P.H.D. student UFOP assistant professor

jones-2

Re: calculate quantiles of a custom function

In reply to this post by VictorDelgado

What do quantiles mean here? If you have a mixture density, say

myf <- function(x,p0) p0*dbeta(x,2,6) + (1-p0)*dbeta(x,6,2)

then I know what quantiles mean. To find the Pth quantile use uniroot
to solve for the x such that myf(x,p0) - P =0.

albyn

Quoting VictorDelgado <[hidden email]>:

>
> Gerhard wrote
>>
>>
>> Suppose I create a custom function, consisting of two beta-distributions:
>>
>> myfunction <- function(x) {
>> dbeta(x,2,6) + dbeta(x,6,2)
>> }
>>
>> How can I calculate the quantiles of myfunction?
>>
>> Thank you in advance,
>>
>> Gerhard
>>
>>
>
> Gehard, if do you want to know the quantiles of the new distribution created
> by "myfunction". Maybe you can also do:
>
> x <- seq(0,1,.01) # insert your 'x'
> q <- myfunction(x)
> # And:
> quantile(x)
>
> 0% 25% 50% 75% 100%
> 0.000000 1.476177 2.045389 2.581226 2.817425
>
> # This gives the sample quantiles. You can also look foward to simulations
> (like Bert Gunter had suggested) to know better the properties of
> distributions quantiles obtained after 'myfunction'.
>
>
>
> -----
> Victor Delgado
> cedeplar.ufmg.br P.H.D. student
> www.fjp.mg.gov.br reseacher
> --
> View this message in context:
> http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp4256887p4257551.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> [hidden email] mailing list
> 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.
>
>

Gerhard Schaden

Re: calculate quantiles of a custom function

In reply to this post by VictorDelgado

Am Dienstag, 3. Januar 2012, 08:50:44 schrieb VictorDelgado:
> VictorDelgado wrote
>
> > quantile(x)
>
> Correcting to
>
> quantile(q)
>
> -----

Dear Victor,

thank you for your answer.

Best,

Gerhard

> Victor Delgado
> cedeplar.ufmg.br P.H.D. student
> www.fjp.mg.gov.br reseacher
> --
> View this message in context:
> http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp42
> 56887p4257575.html Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> [hidden email] mailing list
> 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.

Duncan Murdoch-2

Re: calculate quantiles of a custom function

In reply to this post by jones-2

On 03/01/2012 1:33 PM, Albyn Jones wrote:
> What do quantiles mean here? If you have a mixture density, say
>
> myf<- function(x,p0) p0*dbeta(x,2,6) + (1-p0)*dbeta(x,6,2)
>
> then I know what quantiles mean. To find the Pth quantile use uniroot
> to solve for the x such that myf(x,p0) - P =0.

You would want to integrate first...

Duncan Murdoch

> albyn
>
> Quoting VictorDelgado<[hidden email]>:
>
> >
> > Gerhard wrote
> >>
> >>
> >> Suppose I create a custom function, consisting of two beta-distributions:
> >>
> >> myfunction<- function(x) {
> >> dbeta(x,2,6) + dbeta(x,6,2)
> >> }
> >>
> >> How can I calculate the quantiles of myfunction?
> >>
> >> Thank you in advance,
> >>
> >> Gerhard
> >>
> >>
> >
> > Gehard, if do you want to know the quantiles of the new distribution created
> > by "myfunction". Maybe you can also do:
> >
> > x<- seq(0,1,.01) # insert your 'x'
> > q<- myfunction(x)
> > # And:
> > quantile(x)
> >
> > 0% 25% 50% 75% 100%
> > 0.000000 1.476177 2.045389 2.581226 2.817425
> >
> > # This gives the sample quantiles. You can also look foward to simulations
> > (like Bert Gunter had suggested) to know better the properties of
> > distributions quantiles obtained after 'myfunction'.
> >
> >
> >
> > -----
> > Victor Delgado
> > cedeplar.ufmg.br P.H.D. student
> > www.fjp.mg.gov.br reseacher
> > --
> > View this message in context:
> > http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp4256887p4257551.html
> > Sent from the R help mailing list archive at Nabble.com.
> >
> > ______________________________________________
> > [hidden email] mailing list
> > 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.
> >
> >
>
> ______________________________________________
> [hidden email] mailing list
> 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.

Prof. Dr. Matthias Kohl

Re: calculate quantiles of a custom function

In reply to this post by jones-2

Dear Gerhard,

you could also use package "distr"; e.g.

library(distr)

## use generating function "AbscontDistribution"
D <- AbscontDistribution(d = function(x) dbeta(x, 2, 6) + dbeta(x,6,2),
low = 0, up = 1, withStand = TRUE)

## quantiles
q(D)(seq(0,1,0.1))

Best
Matthias

On 03.01.2012 19:33, Albyn Jones wrote:

> What do quantiles mean here? If you have a mixture density, say
>
> myf <- function(x,p0) p0*dbeta(x,2,6) + (1-p0)*dbeta(x,6,2)
>
> then I know what quantiles mean. To find the Pth quantile use uniroot to
> solve for the x such that myf(x,p0) - P =0.
>
> albyn
>
> Quoting VictorDelgado <[hidden email]>:
>
>>
>> Gerhard wrote
>>>
>>>
>>> Suppose I create a custom function, consisting of two
>>> beta-distributions:
>>>
>>> myfunction <- function(x) {
>>> dbeta(x,2,6) + dbeta(x,6,2)
>>> }
>>>
>>> How can I calculate the quantiles of myfunction?
>>>
>>> Thank you in advance,
>>>
>>> Gerhard
>>>
>>>
>>
>> Gehard, if do you want to know the quantiles of the new distribution
>> created
>> by "myfunction". Maybe you can also do:
>>
>> x <- seq(0,1,.01) # insert your 'x'
>> q <- myfunction(x)
>> # And:
>> quantile(x)
>>
>> 0% 25% 50% 75% 100%
>> 0.000000 1.476177 2.045389 2.581226 2.817425
>>
>> # This gives the sample quantiles. You can also look foward to
>> simulations
>> (like Bert Gunter had suggested) to know better the properties of
>> distributions quantiles obtained after 'myfunction'.
>>
>>
>>
>> -----
>> Victor Delgado
>> cedeplar.ufmg.br P.H.D. student
>> www.fjp.mg.gov.br reseacher
>> --
>> View this message in context:
>> http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp4256887p4257551.html
>>
>> Sent from the R help mailing list archive at Nabble.com.
>>
>> ______________________________________________
>> [hidden email] mailing list
>> 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.
>>
>>
>
> ______________________________________________
> [hidden email] mailing list
> 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.

--
Prof. Dr. Matthias Kohl
www.stamats.de

______________________________________________
[hidden email] mailing list
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.

jones-2

Re: calculate quantiles of a custom function

In reply to this post by Duncan Murdoch-2

right. replace dbetas with pbetas.

albyn

Quoting Duncan Murdoch <[hidden email]>:

> On 03/01/2012 1:33 PM, Albyn Jones wrote:
>> What do quantiles mean here? If you have a mixture density, say
>>
>> myf<- function(x,p0) p0*dbeta(x,2,6) + (1-p0)*dbeta(x,6,2)
>>
>> then I know what quantiles mean. To find the Pth quantile use uniroot
>> to solve for the x such that myf(x,p0) - P =0.
>
> You would want to integrate first...
>
> Duncan Murdoch
>
>> albyn
>>
>> Quoting VictorDelgado<[hidden email]>:
>>
>>>
>>> Gerhard wrote
>>>>
>>>>
>>>> Suppose I create a custom function, consisting of two beta-distributions:
>>>>
>>>> myfunction<- function(x) {
>>>> dbeta(x,2,6) + dbeta(x,6,2)
>>>> }
>>>>
>>>> How can I calculate the quantiles of myfunction?
>>>>
>>>> Thank you in advance,
>>>>
>>>> Gerhard
>>>>
>>>>
>>>
>>> Gehard, if do you want to know the quantiles of the new
>>> distribution created
>>> by "myfunction". Maybe you can also do:
>>>
>>> x<- seq(0,1,.01) # insert your 'x'
>>> q<- myfunction(x)
>>> # And:
>>> quantile(x)
>>>
>>> 0% 25% 50% 75% 100%
>>> 0.000000 1.476177 2.045389 2.581226 2.817425
>>>
>>> # This gives the sample quantiles. You can also look foward to simulations
>>> (like Bert Gunter had suggested) to know better the properties of
>>> distributions quantiles obtained after 'myfunction'.
>>>
>>>
>>>
>>> -----
>>> Victor Delgado
>>> cedeplar.ufmg.br P.H.D. student
>>> www.fjp.mg.gov.br reseacher
>>> --
>>> View this message in context:
>>>
>>> http://r.789695.n4.nabble.com/calculate-quantiles-of-a-custom-function-tp4256887p4257551.html
>>> Sent from the R help mailing list archive at Nabble.com.
>>>
>>> ______________________________________________
>>> [hidden email] mailing list
>>> 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.
>>>
>>>
>>
>> ______________________________________________
>> [hidden email] mailing list
>> 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.
>
>
>

Gerhard

Re: calculate quantiles of a custom function

In reply to this post by Prof. Dr. Matthias Kohl

Am Dienstag, 3. Januar 2012, 19:51:36 schrieb Prof. Dr. Matthias Kohl:
> D <- AbscontDistribution(d = function(x) dbeta(x, 2, 6) + dbeta(x,6,2),
> low = 0, up = 1, withStand = TRUE)

Dear all,

thank you all again for your help.

So, summing up, (in case this might be useful to other beginners - like me)
this is how it can be done:

############################
library(distr)

dcustom <- function(x) {
(dbeta(x,2,6) + dbeta(x,6,2))/2 # I need to divide by 2 to get 1 as
# result of integration;
}

pcustom <- function(x) {
integrate(dmyspeaker,0,x)$value
}

D <- AbscontDistribution(d = dcustom, low = 0, up = 1, withStand = TRUE)

qcustom <- function(x){
q(D)(x)
}
############################

Best,

Gerhard

______________________________________________
[hidden email] mailing list
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.

jones-2

Re: calculate quantiles of a custom function

FWIW, the integral of a mixture density is the same mixture of the
CDFs, so you can use the pbeta functions:

pcustom <- function(x) (pbeta(x,2,6) + pbeta(x,6,2))/2

albyn

Quoting Gerhard <[hidden email]>:

> Am Dienstag, 3. Januar 2012, 19:51:36 schrieb Prof. Dr. Matthias Kohl:
>> D <- AbscontDistribution(d = function(x) dbeta(x, 2, 6) + dbeta(x,6,2),
>> low = 0, up = 1, withStand = TRUE)
>
> Dear all,
>
> thank you all again for your help.
>
> So, summing up, (in case this might be useful to other beginners - like me)
> this is how it can be done:
>
> ############################
> library(distr)
>
> dcustom <- function(x) {
> (dbeta(x,2,6) + dbeta(x,6,2))/2 # I need to divide by 2 to get 1 as
> # result of integration;
> }
>
> pcustom <- function(x) {
> integrate(dmyspeaker,0,x)$value
> }
>
> D <- AbscontDistribution(d = dcustom, low = 0, up = 1, withStand = TRUE)
>
> qcustom <- function(x){
> q(D)(x)
> }
> ############################
>
> Best,
>
> Gerhard
>
> ______________________________________________
> [hidden email] mailing list
> 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.
>
>