R Under development (unstable) (2024-07-10 r86888 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(fitdistrplus)
Loading required package: MASS
Loading required package: survival
> 
> nsample <- 10
> 
> # (1) Fit of a Weibull distribution to serving size data by maximum 
> # goodness-of-fit estimation using all the distances available
> # 
> 
> data(groundbeef)
> serving <- groundbeef$serving
> mgedist(serving, "weibull", gof="CvM")
$estimate
    shape     scale 
 2.093204 82.660014 

$convergence
[1] 0

$value
[1] 0.6556672

$hessian
           shape      scale
shape 4.05295367 0.09244476
scale 0.09244476 0.02418777

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      65       NA 

$optim.message
NULL

$loglik
[1] -1255.623

$gof
[1] "CvM"

> mgedist(serving, "weibull", gof="CvM", silent = FALSE)
$estimate
    shape     scale 
 2.093204 82.660014 

$convergence
[1] 0

$value
[1] 0.6556672

$hessian
           shape      scale
shape 4.05295367 0.09244476
scale 0.09244476 0.02418777

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      65       NA 

$optim.message
NULL

$loglik
[1] -1255.623

$gof
[1] "CvM"

Warning messages:
1: In pweibull(c(10, 11.5, 17, 20, 20, 20, 20, 20, 20, 20, 20, 24,  :
  NaNs produced
2: In pweibull(c(10, 11.5, 17, 20, 20, 20, 20, 20, 20, 20, 20, 24,  :
  NaNs produced
> mgedist(serving, "weibull", gof="KS")
$estimate
    shape     scale 
 2.065634 81.450487 

$convergence
[1] 0

$value
[1] 0.112861

$hessian
           shape    scale
shape 122.668263 6.509057
scale   6.509057 7.599584

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
     127       NA 

$optim.message
NULL

$loglik
[1] -1255.975

$gof
[1] "KS"

> mgedist(serving, "weibull", gof="AD")
$estimate
    shape     scale 
 2.125473 82.890260 

$convergence
[1] 0

$value
[1] 3.501035

$hessian
           shape     scale
shape 29.4165108 0.1823375
scale  0.1823375 0.1354409

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      63       NA 

$optim.message
NULL

$loglik
[1] -1255.392

$gof
[1] "AD"

> mgedist(serving, "weibull", gof="ADR")
$estimate
    shape     scale 
 2.072087 82.761868 

$convergence
[1] 0

$value
[1] 1.610479

$hessian
           shape       scale
shape 13.5240921 -0.33242262
scale -0.3324226  0.07977375

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      65       NA 

$optim.message
NULL

$loglik
[1] -1255.836

$gof
[1] "ADR"

> mgedist(serving, "weibull", gof="ADL")
$estimate
    shape     scale 
 2.197498 82.016005 

$convergence
[1] 0

$value
[1] 1.845939

$hessian
           shape      scale
shape 15.3272022 0.54407117
scale  0.5440712 0.05549883

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      65       NA 

$optim.message
NULL

$loglik
[1] -1255.415

$gof
[1] "ADL"

> mgedist(serving, "weibull", gof="AD2R")
$estimate
   shape    scale 
 1.90328 81.33464 

$convergence
[1] 0

$value
[1] 11.56415

$hessian
          shape       scale
shape 334.61081 -10.4227495
scale -10.42275   0.5223167

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      67       NA 

$optim.message
NULL

$loglik
[1] -1259.112

$gof
[1] "AD2R"

> mgedist(serving, "weibull", gof="AD2L")
$estimate
    shape     scale 
 2.483836 78.252113 

$convergence
[1] 0

$value
[1] 9.786977

$hessian
           shape     scale
shape 113.511932 4.1108355
scale   4.110836 0.2341312

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      69       NA 

$optim.message
NULL

$loglik
[1] -1265.933

$gof
[1] "AD2L"

> mgedist(serving, "weibull", gof="AD2")
$estimate
    shape     scale 
 2.081168 85.281194 

$convergence
[1] 0

$value
[1] 26.95166

$hessian
         shape       scale
shape 534.9606 -10.5940982
scale -10.5941   0.7606462

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      67       NA 

$optim.message
NULL

$loglik
[1] -1256.313

$gof
[1] "AD2"

> 
> # (2) Fit of a uniform distribution using Cramer-von Mises or
> # Kolmogorov-Smirnov distance
> # 
> 
> set.seed(1234)
> u <- runif(nsample, min=5, max=10)
> mgedist(u, "unif", gof="CvM")
$estimate
     min      max 
5.151281 9.808278 

$convergence
[1] 0

$value
[1] 0.1217292

$hessian
          min       max
min 0.2131728 0.1603053
max 0.1603053 0.2961833

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      53       NA 

$optim.message
NULL

$loglik
[1] -Inf

$gof
[1] "CvM"

> mgedist(u, "unif", gof="KS")
$estimate
     min      max 
5.840091 9.452868 

$convergence
[1] 0

$value
[1] 0.2106888

$hessian
          min       max
min 179.94486  55.35421
max  55.35421 179.91068

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
     165       NA 

$optim.message
NULL

$loglik
[1] -Inf

$gof
[1] "KS"

> 
> # (3) Fit of a triangular distribution using Cramer-von Mises or
> # Kolmogorov-Smirnov distance
> # 
> require(mc2d)
Loading required package: mc2d
Loading required package: mvtnorm

Attaching package: 'mc2d'

The following objects are masked from 'package:base':

    pmax, pmin

> set.seed(1234)
> t <- rtriang(nsample, min=5, mode=6, max=10)
> mgedist(t, "triang", start = list(min=4, mode=6,max=9), gof="CvM")
$estimate
     min     mode      max 
4.602048 7.367066 8.682229 

$convergence
[1] 0

$value
[1] 0.08579645

$hessian
           min      mode       max
min  0.1914605 0.2739648 0.1892873
mode 0.2739648 0.7323464 0.4060730
max  0.1892873 0.4060730 0.2910040

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
     180       NA 

$optim.message
NULL

$loglik
[1] -11.84833

$gof
[1] "CvM"

Warning message:
In checkparamlist(arg_startfix$start.arg, arg_startfix$fix.arg,  :
  Some parameter names have no starting/fixed value but have a default value: mean.
> mgedist(t, "triang", start = list(min=4, mode=6,max=9), gof="KS")
$estimate
     min     mode      max 
5.191300 7.414283 8.306100 

$convergence
[1] 0

$value
[1] 0.1853719

$hessian
           min     mode       max
min  148.15115 110.9594  93.74732
mode 110.95941 292.0984 208.67195
max   93.74732 208.6720 208.66467

$optim.function
[1] "optim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
     298       NA 

$optim.message
NULL

$loglik
[1] -Inf

$gof
[1] "KS"

Warning message:
In checkparamlist(arg_startfix$start.arg, arg_startfix$fix.arg,  :
  Some parameter names have no starting/fixed value but have a default value: mean.
> 
> 
> # (4) scaling problem
> # the simulated dataset (below) has particularly small values, hence without scaling (10^0),
> # the optimization raises an error. The for loop shows how scaling by 10^i
> # for i=1,...,6 makes the fitting procedure work correctly.
> if(FALSE)
+ {
+   set.seed(1234)
+   x2 <- rnorm(nsample, 1e-4, 2e-4)
+   for(i in 6:0)
+     cat(i, try(mgedist(x2*10^i,"cauchy")$estimate, silent=TRUE), "\n")
+ }
> 
> # (5) scaling problem
> #
> if(FALSE)
+ {
+   x <- c(-0.00707717, -0.000947418, -0.00189753, 
+          -0.000474947, -0.00190205, -0.000476077, 0.00237812, 0.000949668, 
+          0.000474496, 0.00284226, -0.000473149, -0.000473373, 0, 0, 0.00283688, 
+          -0.0037843, -0.0047506, -0.00238379, -0.00286807, 0.000478583, 
+          0.000478354, -0.00143575, 0.00143575, 0.00238835, 0.0042847, 
+          0.00237248, -0.00142281, -0.00142484, 0, 0.00142484, 0.000948767, 
+          0.00378609, -0.000472478, 0.000472478, -0.0014181, 0, -0.000946522, 
+          -0.00284495, 0, 0.00331832, 0.00283554, 0.00141476, -0.00141476, 
+          -0.00188947, 0.00141743, -0.00236351, 0.00236351, 0.00235794, 
+          0.00235239, -0.000940292, -0.0014121, -0.00283019, 0.000472255, 
+          0.000472032, 0.000471809, -0.0014161, 0.0014161, -0.000943842, 
+          0.000472032, -0.000944287, -0.00094518, -0.00189304, -0.000473821, 
+          -0.000474046, 0.00331361, -0.000472701, -0.000946074, 0.00141878, 
+          -0.000945627, -0.00189394, -0.00189753, -0.0057143, -0.00143369, 
+          -0.00383326, 0.00143919, 0.000479272, -0.00191847, -0.000480192, 
+          0.000960154, 0.000479731, 0, 0.000479501, 0.000958313, -0.00383878, 
+          -0.00240674, 0.000963391, 0.000962464, -0.00192586, 0.000481812, 
+          -0.00241138, -0.00144963)
+   
+   
+   #only i == 0, no scaling, should not converge.
+   for(i in 6:0)
+     cat(i, try(mgedist(x*10^i,"cauchy")$estimate, silent=TRUE), "\n")
+ }
> 
> 
> # (6) test error messages
> #
> 
> dnorm2 <- pnorm2 <- function(x, a)
+   "NA"
> x <- rexp(10)
> 
> #should get a one-line error 
> res <- mgedist(x, "norm2", start=list(a=1))
<simpleError in theop - (2 * 1:n - 1)/(2 * n): non-numeric argument to binary operator>
> #as in 
> attr(try(log("a"), silent=TRUE), "condition")
<simpleError in log("a"): non-numeric argument to mathematical function>
> 
> 
> try(mgedist(c(serving, "a"), "gamma"))
Error in mgedist(c(serving, "a"), "gamma") : 
  data must be a numeric vector of length greater than 1 for non censored data
            or a dataframe with two columns named left and right and more than one line for censored data
> try(mgedist(c(serving, NA), "gamma"))
Error in checkUncensoredNAInfNan(data) : data contain NA values.
> try(mgedist(c(serving, Inf), "gamma"))
Error in checkUncensoredNAInfNan(data) : 
  data contain Inf (infinite) values.
> try(mgedist(c(serving, -Inf), "gamma"))
Error in checkUncensoredNAInfNan(data) : 
  data contain Inf (infinite) values.
> try(mgedist(c(serving, NaN), "gamma"))
Error in checkUncensoredNAInfNan(data) : 
  data contain NaN (not a numeric) values.
> 
> 
> # (7) test the component optim.message
> x <- rnorm(1000)
> #change parameter to obtain unsuccessful convergence
> mgedist(x, "norm", control=list(maxit=2), start=list(mean=1e5, sd=1), optim.method="L-BFGS-B", lower=0)
$estimate
 mean    sd 
1e+05 1e+00 

$convergence
[1] 0

$value
[1] 333.3333

$hessian
     mean sd
mean    0  0
sd      0  0

$optim.function
[1] "optim"

$optim.method
[1] "L-BFGS-B"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
       1        1 

$optim.message
[1] "CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL"

$loglik
[1] -Inf

$gof
[1] "CvM"

> 
> # (8) test bounds
> x <- rnorm(1000)
> mgedist(x, "norm", optim.method="L-BFGS-B", lower=c(-Inf, 0)) #optim and L-BFGS-B
$estimate
      mean         sd 
0.01134416 0.97945793 

$convergence
[1] 0

$value
[1] 0.02003669

$hessian
             mean          sd
mean 191.57400735 -0.05237777
sd    -0.05237777 62.19267021

$optim.function
[1] "optim"

$optim.method
[1] "L-BFGS-B"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
       8        8 

$optim.message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

$loglik
[1] -1399.894

$gof
[1] "CvM"

> mgedist(x, "norm", optim.method="Nelder", lower=c(-Inf, 0))
$estimate
      mean         sd 
0.01134315 0.97945779 

$convergence
[1] 0

$value
[1] 0.02003669

$hessian
NULL

$optim.function
[1] "constrOptim"

$optim.method
[1] "Nelder-Mead"

$fix.arg
NULL

$fix.arg.fun
NULL

$weights
NULL

$counts
function gradient 
      59       NA 

$optim.message
NULL

$loglik
[1] -1399.894

$gof
[1] "CvM"

> 
> 
> 
> proc.time()
   user  system elapsed 
   2.64    0.45    3.12