[R] Dirichlet surface

[Kehl Dániel]

It helped a lot indeed, thank you very much!
Now I understand why it was a problem for persp!

Daniel

2011-03-30 10:31 keltezéssel, Gavin Simpson írta:
> On Wed, 2011-03-30 at 09:55 -0700, Kehl Dániel wrote:
>> Dear David,
>>
>> I think that is a small bug too, maybe because the function is constant?
>> is there a nice way to put the c(0,2.1) argument optionally, only if all
>> the parameters are 1?
>> Should I post the problem somewhere else (developers maybe?)
> I don't think this is a bug really; the code is having to compute limits
> of the z axis and you supplied it one bit of information: a 2. If your
> data are so degenerate then it is not unreasonable to expect some user
> intervention. Admittedly, persp doesn't seem to work like other R
> plotting functions.
>
> You could do something like:
>
> persp(x1, x2, z,
>        zlim = if(length(na.omit(unique(as.vector(z))))<  2){ c(0,2.1) }
> else { NULL})
>
> in your call to persp so it only uses user-defined limits if the number
> of numeric values in z is less than 2.
>
> HTH
>
> G
>
>> thanks:
>> Daniel
>>
>> 2011-03-30 04:42 keltezéssel, David Winsemius írta:
>>> On Mar 29, 2011, at 4:45 PM, Kehl Dániel wrote:
>>>
>>>> Dear list members,
>>>>
>>>> I want to draw surfaces of Dirichlet distributions with different
>>>> parameter settings.
>>>> My code is the following:
>>>> #<begin code>
>>>> a1<- a2<- a3<- 2
>>>> #a2<- .5
>>>> #a3<- .5
>>>> x1<- x2<- seq(0.01, .99, by=.01)
>>>>
>>>> f<- function(x1, x2){
>>>>       term1<- gamma(a1+a2+a3)/(gamma(a1)*gamma(a2)*gamma(a3))
>>>>       term2<- x1^(a1-1)*x2^(a2-1)*(1-x1-x2)^(a3-1)
>>>>       term3<- (x1 + x2<  1)
>>>>       term1*term2*term3
>>>>       }
>>>>
>>>> z<- outer(x1, x2, f)
>>>> z[z<=0]<- NA
>>>>
>>>> persp(x1, x2, z,
>>>>       main = "Dirichlet Distribution",
>>>>       col = "lightblue",
>>>>       theta = 50,
>>>>       phi = 20,
>>>>       r = 50,
>>>>       d = 0.1,
>>>>       expand = 0.5,
>>>>       ltheta = 90,
>>>>       lphi = 180,
>>>>       shade = 0.75,
>>>>       ticktype = "detailed",
>>>>       nticks = 5)
>>>> #<end code>
>>>>
>>>> It works fine (I guess), except for a1=a2=a3=1. In that case I get
>>>> the error: Error in persp.default...  invalid 'z' limits.
>>>> The z matrix has only elements 2 and NA.
>>> Might be a buglet in persp. If you set the zlim argument to c(0,2.1),
>>> you get the expected constant plane at z=2 for those arguments.
>>>> Any ideas are appreciated.
>>>>
>>>> Thank you:
>>>> Daniel
>>>> University of Pécs
>>>>
>>>> ______________________________________________
>>>> R-help at r-project.org 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, MD
>>> Heritage Laboratories
>>> West Hartford, CT
>>>
>>>
>>>
>> ______________________________________________
>> R-help at r-project.org 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.