Dear all, I am trying to use the rstable(n, alpha, beta, gamma = 1, delta = 0, pm = c(0, 1, 2))) function to generate positive stable random numbers. For positive stable distribution, beta==1 and alpha is in (0,1), which defines random variables with support (0, infinity). So, I used rstable(100, 0.5, 1) for an example. I found that this gives me some negative numbers. For example,> rstable(10, 0.5, 1)[1] 6.3016252 399.3659030 11.2735789 1.9550625 -0.6762333 1.6810761 [7] 0.9091360 1.9100991 -0.7593737 24.2788471 Does anybody know why this should happen? Thanks a lot, Julie
Peter Dalgaard
2005-Aug-30 19:55 UTC
[R] problem in generating positive stable random numbers
"X. Cong" <xcong at stat.rice.edu> writes:> Dear all, > > I am trying to use the > rstable(n, alpha, beta, gamma = 1, delta = 0, pm = c(0, 1, 2))) > > function to generate positive stable random numbers. For positive stable > distribution, beta==1 and alpha is in (0,1), which defines random variables > with support (0, infinity). So, I used rstable(100, 0.5, 1) for an example. > I found that this gives me some negative numbers. For example, > > > rstable(10, 0.5, 1) > [1] 6.3016252 399.3659030 11.2735789 1.9550625 -0.6762333 1.6810761 > [7] 0.9091360 1.9100991 -0.7593737 24.2788471 > > Does anybody know why this should happen?Doesn't sound like it should...>From where did you get rstable()? It isn't part of R itself, andpoking around shows at least three different sources for it (CircStats, fBasics, stable (Lambert/Lindsey, not on CRAN)). I seem to recall that the third one also has a distribution function which is almost, but not quite monotone. The usual advice of taking package problems to the package maintainer applies. -- O__ ---- Peter Dalgaard ??ster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907
Diethelm Wuertz
2005-Aug-30 22:32 UTC
[R] problem in generating positive stable random numbers
X. Cong wrote:>Dear all, > >I am trying to use the >rstable(n, alpha, beta, gamma = 1, delta = 0, pm = c(0, 1, 2))) > > function to generate positive stable random numbers. For positive stable >distribution, beta==1 and alpha is in (0,1), which defines random variables >with support (0, infinity). >Just a remark - rstable() is from R-package fBasics ... 1) I think the support for beta=1 and alpha = 1/2 is (-1, infinity), isn't it? - Then everything is fine. 2) beta ==1 is a difficult value for numerical computations, try also beta = 1-1e-8! 3) Use the program stable.exe from http://academic2.american.edu/~jpnolan/stable/stable.html for comparison! 4) Take care in which parametrization you work ... Best regards Diethelm Remark for 3) STABLE 3.14.02 (2005/02/28) Serial number 131 Copyright 1997-2003 John P. Nolan (jpnolan at american.edu) Output file: stable.out Current tolerance settings: -1 debug (F) 1 relative error for pdf ( 0.1200000000E-13) 2 relative error for cdf ( 0.1200000000E-13) 3 relative error for quantiles ( 0.1200000000E-13) 4 alpha and beta rounding ( 0.1000000000E-01) 5 x tolerance near zeta ( 0.5000000000E-02) 6 exponential cutoff ( 200.0000000 ) 7 peak/strim location tolerance ( 0.1000000000E-13) 8 strim tolerance ( 0.1000000000E-50) 9 minimum alpha ( 0.1000000000 ) 10 minimum xtol ( 0.1000000000E-12) 11 threshold for quantile search ( 0.1000000000E-09) 8/30/2005 23:59:34.72 Simulation of stable random variables n= 10 iseed= -1 iparam= 0 alpha beta gamma delta 0.50000 1.00000 1.00000 0.00000 39.6516320412757 0.247648247004956 12.3195079629881 -0.146866807337654 0.640815989111718 0.219559949961403 3.05027005351445 -0.667188061257615 4.02760809108209 -0.588547074136516>So, I used rstable(100, 0.5, 1) for an example. >I found that this gives me some negative numbers. For example, > > > >>rstable(10, 0.5, 1) >> >> > [1] 6.3016252 399.3659030 11.2735789 1.9550625 -0.6762333 1.6810761 > [7] 0.9091360 1.9100991 -0.7593737 24.2788471 > >Does anybody know why this should happen? > >Thanks a lot, >Julie > >______________________________________________ >R-help at stat.math.ethz.ch mailing list >https://stat.ethz.ch/mailman/listinfo/r-help >PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html > > >