1
# This R package is free software; you can redistribute it and/or
2
# modify it under the terms of the GNU Library General Public
3
# License as published by the Free Software Foundation; either
4
# version 2 of the License, or (at your option) any later version.
6
# This R package is distributed in the hope that it will be useful,
7
# but WITHOUT ANY WARRANTY; without even the implied warranty of
8
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9
# GNU Library General Public License for more details.
11
# You should have received a copy of the GNU Library General
12
# Public License along with this R package; if not, write to the
13
# Free Foundation, Inc., 59 Temple Place, Suite 330, Boston,
17
################################################################################
18
# FUNCTIONS: DESCRIPTION:
19
# stableMode Computes the mode of the stable DF
20
################################################################################
22
##' Computes the mode of the alpha stable distribution
23
##' @title Mode of the stable distribution
26
##' @param beta.max for numerical purposes, values of beta too close to 1,
27
##' are set to beta.max
28
##' @param tol numerical tolerance used in optimize()
29
##' @return a number, the stable mode
30
##' @author Diethelm Wuertz and Martin Maechler
31
stableMode <- function(alpha, beta, beta.max = 1 - 1e-11,
32
tol = .Machine$double.eps^0.25)
34
stopifnot(0 < alpha, alpha <= 2, length(alpha) == 1,
35
-1 <= beta, beta <= 1, length(beta) == 1,
36
length(beta.max) == 1)
38
# # Test for values close to beta = 1
39
# alpha <- seq(0, 2, by = 0.1)
40
# ans <- matrix(NA, nrow=length(alpha), ncol = 4)
41
# for (i in 1:seq_along(alpha)) {
43
# stableMode(alpha = alpha[i], beta = 0.99 ),
44
# stableMode(alpha = alpha[i], beta = 0.99999 ),
45
# stableMode(alpha = alpha[i], beta = 0.99999999 ),
46
# stableMode(alpha = alpha[i], beta = 0.99999999999))
50
# alpha 0.99 0.99999 0.99999999 0.99999999999
51
# 0.0 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
52
# 0.2 -3.214142e-01 -3.246759e-01 -3.246787e-01 -3.246788e-01
53
# 0.4 -6.105318e-01 -6.158562e-01 -6.158616e-01 -6.158616e-01
54
# 0.6 -6.550106e-01 -6.594746e-01 -6.594790e-01 -6.594790e-01
55
# 0.8 -5.558811e-01 -5.590032e-01 -5.590063e-01 -5.590063e-01
56
# 1.0 -4.271033e-01 -4.293078e-01 -4.293099e-01 -4.293099e-01
57
# 1.2 -3.074015e-01 -3.090820e-01 -3.090804e-01 -3.090804e-01
58
# 1.4 -2.050956e-01 -2.063979e-01 -2.063951e-01 -2.063951e-01
59
# 1.6 -1.199623e-01 -1.208875e-01 -1.208853e-01 -1.208853e-01
60
# 1.8 -5.098617e-02 -5.145758e-02 -5.145639e-02 -5.145639e-02
61
# 2.0 -7.487432e-05 -7.487432e-05 -7.487432e-05 -7.487432e-05
66
if(beta > beta.max) beta <- beta.max
68
optimize(dstable, interval = c(-0.7, 0)*sign(beta),
69
alpha = alpha, beta = beta, pm = 0,
70
maximum = TRUE, tol = tol)$maximum
added
removed
R/dist-stableMode.R
