[R] Dirichlet surface

Kehl Dániel kehld at ktk.pte.hu
Wed Mar 30 20:12:20 CEST 2011


Actually, it works for the a=1 case, not for the others. It still gives 
the invalid 'zlim argument' error.
I'll try to work it out maybe instead of NULL giving a c which is 
dependent on the max(z).

Daniel

2011-03-30 10:42 keltezéssel, Kehl Dániel írta:
> 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.
>
> ______________________________________________
> 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.



More information about the R-help mailing list