...Additionally...your set of parameters match none of the curves in figure 4.

I think the ordering of the parameters as listed on the graphs is different 
than in the text of the article.

The 'v' parameter controls the location of the 'elbow' and should be near 
log(x) in each graph, while the 't' parameter is the sharpness of the elbow.  
So eyeballing the elbows:

   log(c(80,150,300))= 4.382027 5.010635 5.703782

These appear to match positional parameter #4 in the legends, which you apply 
as parameter t in your function.  


# editing your function a bit:

ddpln <- function(x, a=2.5, b=0.01, v=log(100), t=v/10){
  # Density of Double Pareto LogNormal distribution
  # from "b. alzahrani" <cs_2...@hotmail.com> email of 2013-11-13
  # per formula 8 from http://cs.stanford.edu/people/jure/pubs/dpln-kdd08.pdf

   c <- (a * b /(a+b))

   norm1<-pnorm((log(x)-v-(a*t^2))/t)
   norm2<-pnorm((log(x)-v+(b*t^2))/t)
   expo1<-  a*v+(a^2*t^2/2)
   expo2<- -b*v+(b^2*t^2/2)  # edited from the paper's eqn 8  with: s/t/v/

   z<- (x^(-a-1)) * exp(expo1)*(norm1)
   y<- (x^(b-1))  * exp(expo2)*(1-norm2)  # 1-norm is the complementary CDF of 
N(0,1)

   c*(z+y)
}

x<-10^seq(0,5,by=0.1) # 0-10k

plot(x,ddpln(x,a=2.8,b=.01,v=3.8,t=0.35),log='xy',type='l')  # Similar to fig 4 
left.

plot(x,ddpln(x,a=2.5,b=.01,v=log(2)),log='xy',type='l')
plot(x,ddpln(x,a=2.5,b=.01,v=log(10)),log='xy',type='l')
plot(x,ddpln(x,a=2.5,b=.01,v=log(100)),log='xy',type='l')
plot(x,ddpln(x,a=2.5,b=.01,v=log(1000)),log='xy',type='l')
plot(x,ddpln(x,a=2.5,b=.01,v=log(10000)),log='xy',type='l')

Dave

On Nov 13, 2013, at 11:43 AM, "b. alzahrani" <cs_2...@hotmail.com>
 wrote:

> 
> Hi
> 
> I found this paper http://cs.stanford.edu/people/jure/pubs/dpln-kdd08.pdf 
> that models the DPLN distribution as in equation 8. I implemented this in R 
> but cannot get the same curve as in Figure 4. can you please check if my code 
> below is correct: e.g. is the use of pnorm() correct here?
> 
> ddlpn <- function(x){
>  a=2.8
>  b=0.01
>  v=0.45
>  t=6.5
>  j <- (a * b /(a+b))
> 
>  norm1<-pnorm((log(x)-v-(a*t^2))/t)
>  expo1<- a*v+(a^2*t^2/2)
> 
>  z<-exp(expo1)*(x^(-a-1))*(norm1)
> 
>  norm2<-pnorm((log(x)-v+(b*t^2))/t)
> expo2<- -b*t+(b^2*t^2/2)
> 
>  y<- x^(b-1)*exp(expo2)*(1-norm2)  # 1-norm is the complementary CDF of N(0,1)
>  j*(z+y)
> }
> ******************************************************************
> Bander Alzahrani, Teacher Assistant
> Information Systems Department
> Faculty of Computing & Information Technology
> King Abdulaziz University
> 
> *************************************
> 
> 
> 
>> From: d...@vims.edu
>> To: cs_2...@hotmail.com
>> CC: r-h...@stat.math.ethz.ch
>> Subject: Re: [R] Double Pareto Log Normal Distribution
>> Date: Tue, 12 Nov 2013 16:51:22 +0000
>> 
>> 
>> http://www.math.uvic.ca/faculty/reed/dPlN.3.pdf is the original ref and has 
>> the equations.
>> 
>> library(VGAM) for *pareto() and library(stats) for *lnorm() should get you 
>> most of the way there.
>> 
>> On Nov 12, 2013, at 10:47 AM, "b. alzahrani" <cs_2...@hotmail.com>
>> wrote:
>> 
>>> Hi guys
>>> I would like to generate random number Double Pareto Log Normal 
>>> Distribution (DPLN). does anyone know how to do this in R or if there is 
>>> any built-in function.
>>> 
>>> Thanks
>>> 
>>> ******************************************************************
>>> Bander 
>>> *************************************
>>> 
>>> 
>>>                                       
>>>     [[alternative HTML version deleted]]
>>> 
>>> ______________________________________________
>>> R-help@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.
>> 
>> --
>> Dr. David Forrest
>> d...@vims.edu
>> 
>> 
>> 
>> 
>> 
>                                         
>       [[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help@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.

--
Dr. David Forrest
d...@vims.edu




Attachment: signature.asc
Description: Message signed with OpenPGP using GPGMail

______________________________________________
R-help@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.

Reply via email to