# OW distribution

#### 2022-12-21

library(RelDists)
library(EstimationTools)
library(gamlss)
#> Loading required package: splines
#> Loading required package: gamlss.data
#>
#> Attaching package: 'gamlss.data'
#> The following object is masked from 'package:datasets':
#>
#>     sleep
#> Loading required package: gamlss.dist
#> Loading required package: MASS
#> Loading required package: nlme
#>
#> Attaching package: 'nlme'
#> The following object is masked from 'package:BBmisc':
#>
#>     collapse
#>  **********   GAMLSS Version 5.4-10  **********
#> For more on GAMLSS look at https://www.gamlss.com/
#> Type gamlssNews() to see new features/changes/bug fixes.

# Odd Weibull distribution

This distribution was proposed by Cooray (2006). The probability density function $$f(x)$$ and cumulative density function $$F(x)$$ are given by:

$f(x) = \left( \frac{\sigma\nu}{x} \right) (\mu x)^\sigma e^{(\mu x)^\sigma} \left(e^{(\mu x)^{\sigma}}-1\right)^{\nu-1} \left[ 1 + \left(e^{(\mu x)^{\sigma}}-1\right)^\nu \right]^{-2},$

and

$F(x) = 1 - \left[ 1 + \left( e^{(\mu x)^{\sigma}} - 1 \right)^\nu \right]^{-1}, \quad x>0,$

respectively, where $$\mu>0, \quad \sigma\nu>0$$. $$\mu$$ is the scale parameter, $$\sigma$$ and $$\nu$$ are the shape parameters. Next figure shows possible shapes of the $$f(x)$$ and $$F(x)$$ for several values of the parameters.

The hazard function (hf) of the OW distribution is:

$h(x) = \left( \frac{\sigma\nu}{x} \right) (\mu x)^\sigma e^{(\mu x)^\sigma} \left(e^{(\mu x)^{\sigma}}-1\right)^{\nu-1} \left[ 1 + \left(e^{(\mu x)^{\sigma}}-1\right)^\nu \right]^{-1}, x>0,$

where the hf can take the following shapes:

• Increasing if $$\sigma>1$$ and $$\sigma\nu>1$$.

• Decreasing if $$\sigma<1$$ and $$\sigma\nu<1$$.

• Unimodal shaped if either $$\sigma<0$$ and $$\nu<0$$ or $$\sigma<1$$ and $$\sigma\nu\geq 1$$.

• Bathtub shaped if $$\sigma>1$$ and $$\sigma\nu\geq 1$$.

The figure shows possible shapes of the hf mentioned above:

The following figure illustrate the regions corresponding to the different hazard shapes:

# Application examples

## Time to failure on electronic equipment

Cooray (2015) used the following data provided by Wang (2000) in order to fit an OW distribution through maximum likelihood estimation (MLE):

5, 11, 21, 31, 46, 75, 98, 122, 145, 165, 195, 224, 245, 293, 321, 330, 350, 420.

The data above is the time to failure of an electronic device in hours. As an alternative to classical MLE, We used the function to fit an only-intercept model in order to estimate parameters of OW distribution without covariates. Using our initValuesOW(), we can obtain an initial guess and the valid region.

data("equipment")
my_initial_guess <- initValuesOW(formula = equipment ~ 1)

summary(my_initial_guess)
#> --------------------------------------------------------------------
#> Initial Values
#> sigma = 5
#> nu = 0.1
#> --------------------------------------------------------------------
#> Search Regions
#> For sigma: all(sigma > 1)
#> For nu: all(nu < 1/sigma)
#> --------------------------------------------------------------------
#> Hazard shape: Bathtub

initValuesOW() function detected the Bathtub hazard shape, which corresponds to a convex-then-concave shape of total time on test (TTT) plot

old_par <- par(mfrow = c(1, 1)) # save previous graphical parameters

par(mar = c(3.7, 3.7, 1, 10), mgp = c(2.5, 1, 0))
plot(my_initial_guess, las = 1)
legend.HazardShape(x = 1.07, y = 1.04, xpd = TRUE)


par(old_par) # restore previous graphical parameters

Therefore, we define the search region according to initValuesOW() outputs

# Custom search region
myvalues <- list(sigma = "all(sigma > 1)",
nu = "all(nu < 1/sigma)")

and we perform the fit using gamlss()

# gamlss set up
con.out <-gamlss.control(n.cyc = 300, trace=TRUE)
myOW <- myOW_region(family = OW(sigma.link='identity'),
valid.values = myvalues, initVal = my_initial_guess)

param_ee <- gamlss(equipment ~ 1, sigma.fo = ~ 1, nu.fo = ~ 1,
sigma.start = 5, nu.start = 0.1,
control = con.out, family = myOW)
#> GAMLSS-RS iteration 1: Global Deviance = 220.3699
#> GAMLSS-RS iteration 2: Global Deviance = 218.2162
#> GAMLSS-RS iteration 3: Global Deviance = 217.5282
#> GAMLSS-RS iteration 4: Global Deviance = 217.3291
#> GAMLSS-RS iteration 5: Global Deviance = 217.2142
#> GAMLSS-RS iteration 6: Global Deviance = 217.1047
#> GAMLSS-RS iteration 7: Global Deviance = 216.9932
#> GAMLSS-RS iteration 8: Global Deviance = 216.8872
#> GAMLSS-RS iteration 9: Global Deviance = 216.7914
#> GAMLSS-RS iteration 10: Global Deviance = 216.707
#> GAMLSS-RS iteration 11: Global Deviance = 216.6337
#> GAMLSS-RS iteration 12: Global Deviance = 216.5706
#> GAMLSS-RS iteration 13: Global Deviance = 216.5168
#> GAMLSS-RS iteration 14: Global Deviance = 216.471
#> GAMLSS-RS iteration 15: Global Deviance = 216.4324
#> GAMLSS-RS iteration 16: Global Deviance = 216.3997
#> GAMLSS-RS iteration 17: Global Deviance = 216.3723
#> GAMLSS-RS iteration 18: Global Deviance = 216.3492
#> GAMLSS-RS iteration 19: Global Deviance = 216.3299
#> GAMLSS-RS iteration 20: Global Deviance = 216.3136
#> GAMLSS-RS iteration 21: Global Deviance = 216.3
#> GAMLSS-RS iteration 22: Global Deviance = 216.2885
#> GAMLSS-RS iteration 23: Global Deviance = 216.2789
#> GAMLSS-RS iteration 24: Global Deviance = 216.2708
#> GAMLSS-RS iteration 25: Global Deviance = 216.264
#> GAMLSS-RS iteration 26: Global Deviance = 216.2583
#> GAMLSS-RS iteration 27: Global Deviance = 216.2533
#> GAMLSS-RS iteration 28: Global Deviance = 216.2488
#> GAMLSS-RS iteration 29: Global Deviance = 216.245
#> GAMLSS-RS iteration 30: Global Deviance = 216.2418
#> GAMLSS-RS iteration 31: Global Deviance = 216.2391
#> GAMLSS-RS iteration 32: Global Deviance = 216.2369
#> GAMLSS-RS iteration 33: Global Deviance = 216.2351
#> GAMLSS-RS iteration 34: Global Deviance = 216.2336
#> GAMLSS-RS iteration 35: Global Deviance = 216.232
#> GAMLSS-RS iteration 36: Global Deviance = 216.2306
#> GAMLSS-RS iteration 37: Global Deviance = 216.2294
#> GAMLSS-RS iteration 38: Global Deviance = 216.2284
summary(param_ee)
#> ******************************************************************
#> Family:  c("OW", "Odd Weibull")
#>
#> Call:  gamlss(formula = equipment ~ 1, sigma.formula = ~1,
#>     nu.formula = ~1, family = myOW, sigma.start = 5,
#>     nu.start = 0.1, control = con.out)
#>
#> Fitting method: RS()
#>
#> ------------------------------------------------------------------
#> Mu link function:  log
#> Mu Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept)   -5.240      0.205  -25.56 8.79e-14 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> Sigma link function:  identity
#> Sigma Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept)    3.283      1.423   2.307   0.0357 *
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> Nu link function:  log
#> Nu Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept)  -1.2715     0.5019  -2.534   0.0229 *
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> No. of observations in the fit:  18
#> Degrees of Freedom for the fit:  3
#>       Residual Deg. of Freedom:  15
#>                       at cycle:  38
#>
#> Global Deviance:     216.2284
#>             AIC:     222.2284
#>             SBC:     224.8996
#> ******************************************************************

In the following table we summarize the results and compare them with those gotten by Cooray (2015)

Parameter MLE (Cooray 2015) GAMLSS
$$\mu$$ 5.35e-03 5.30e-03
$$\sigma$$ 3.22388 3.2828
$$\nu$$ 0.28424 0.2804

## Mortality in mice exposed to radiation

Cooray (2006) used a dataset with 208 data points provided by Kimball (1960) with the ages at death in weeks for male mice exposed to 240r of gamma radiation. Again, we implement a workflow for parameter estimation with myOW_region and gamlss functions.

# Do not forget to load 'RelDists' package
data("mice")
init_vals <- initValuesOW(formula = mice ~ 1)

summary(init_vals)
#> --------------------------------------------------------------------
#> Initial Values
#> sigma = 2
#> nu = 6
#> --------------------------------------------------------------------
#> Search Regions
#> For sigma: all(sigma > 1)
#> For nu: all(nu > 1/sigma)
#> --------------------------------------------------------------------
#> Hazard shape: Increasing

With initValuesOW() function we identified an increasing hazard shape, as well as was stated by Cooray (2006), because TTT plot is concave.

old_par <- par(mfrow = c(1, 1)) # save previous graphical parameters

par(mar = c(3.7, 3.7, 1, 10), mgp = c(2.5, 1, 0))
plot(init_vals, las = 1)
legend.HazardShape(x = 1.07, y = 1.04, xpd = TRUE)


par(old_par) # restore previous graphical parameters

Hence, we implement the estimation procedure

# gamlss set up
myOW <- myOW_region(initVal = init_vals)

# Custom search region
# Do not forget to load 'gamlss' library
param_mm <- gamlss(mice ~ 1, sigma.fo = ~ 1, nu.fo = ~ 1,
sigma.start = 2, nu.start = 6,
control = con.out,
family = myOW)
#> GAMLSS-RS iteration 1: Global Deviance = 2031.511
#> GAMLSS-RS iteration 2: Global Deviance = 2006.257
#> GAMLSS-RS iteration 3: Global Deviance = 2005.168
#> GAMLSS-RS iteration 4: Global Deviance = 2005.063
#> GAMLSS-RS iteration 5: Global Deviance = 2004.941
#> GAMLSS-RS iteration 6: Global Deviance = 2004.817
#> GAMLSS-RS iteration 7: Global Deviance = 2004.69
#> GAMLSS-RS iteration 8: Global Deviance = 2004.558
#> GAMLSS-RS iteration 9: Global Deviance = 2004.423
#> GAMLSS-RS iteration 10: Global Deviance = 2004.283
#> GAMLSS-RS iteration 11: Global Deviance = 2004.139
#> GAMLSS-RS iteration 12: Global Deviance = 2003.991
#> GAMLSS-RS iteration 13: Global Deviance = 2003.838
#> GAMLSS-RS iteration 14: Global Deviance = 2003.68
#> GAMLSS-RS iteration 15: Global Deviance = 2003.517
#> GAMLSS-RS iteration 16: Global Deviance = 2003.348
#> GAMLSS-RS iteration 17: Global Deviance = 2003.174
#> GAMLSS-RS iteration 18: Global Deviance = 2002.993
#> GAMLSS-RS iteration 19: Global Deviance = 2002.807
#> GAMLSS-RS iteration 20: Global Deviance = 2002.613
#> GAMLSS-RS iteration 21: Global Deviance = 2002.413
#> GAMLSS-RS iteration 22: Global Deviance = 2002.205
#> GAMLSS-RS iteration 23: Global Deviance = 2001.989
#> GAMLSS-RS iteration 24: Global Deviance = 2001.766
#> GAMLSS-RS iteration 25: Global Deviance = 2001.533
#> GAMLSS-RS iteration 26: Global Deviance = 2001.292
#> GAMLSS-RS iteration 27: Global Deviance = 2001.041
#> GAMLSS-RS iteration 28: Global Deviance = 2000.78
#> GAMLSS-RS iteration 29: Global Deviance = 2000.507
#> GAMLSS-RS iteration 30: Global Deviance = 2000.224
#> GAMLSS-RS iteration 31: Global Deviance = 1999.929
#> GAMLSS-RS iteration 32: Global Deviance = 1999.621
#> GAMLSS-RS iteration 33: Global Deviance = 1999.299
#> GAMLSS-RS iteration 34: Global Deviance = 1998.964
#> GAMLSS-RS iteration 35: Global Deviance = 1998.613
#> GAMLSS-RS iteration 36: Global Deviance = 1998.247
#> GAMLSS-RS iteration 37: Global Deviance = 1997.864
#> GAMLSS-RS iteration 38: Global Deviance = 1997.463
#> GAMLSS-RS iteration 39: Global Deviance = 1997.037
#> GAMLSS-RS iteration 40: Global Deviance = 1996.591
#> GAMLSS-RS iteration 41: Global Deviance = 1996.124
#> GAMLSS-RS iteration 42: Global Deviance = 1995.634
#> GAMLSS-RS iteration 43: Global Deviance = 1995.121
#> GAMLSS-RS iteration 44: Global Deviance = 1994.584
#> GAMLSS-RS iteration 45: Global Deviance = 1994.022
#> GAMLSS-RS iteration 46: Global Deviance = 1993.434
#> GAMLSS-RS iteration 47: Global Deviance = 1992.822
#> GAMLSS-RS iteration 48: Global Deviance = 1992.183
#> GAMLSS-RS iteration 49: Global Deviance = 1991.518
#> GAMLSS-RS iteration 50: Global Deviance = 1990.826
#> GAMLSS-RS iteration 51: Global Deviance = 1990.111
#> GAMLSS-RS iteration 52: Global Deviance = 1989.373
#> GAMLSS-RS iteration 53: Global Deviance = 1988.617
#> GAMLSS-RS iteration 54: Global Deviance = 1987.847
#> GAMLSS-RS iteration 55: Global Deviance = 1987.068
#> GAMLSS-RS iteration 56: Global Deviance = 1986.285
#> GAMLSS-RS iteration 57: Global Deviance = 1985.506
#> GAMLSS-RS iteration 58: Global Deviance = 1984.74
#> GAMLSS-RS iteration 59: Global Deviance = 1983.994
#> GAMLSS-RS iteration 60: Global Deviance = 1983.277
#> GAMLSS-RS iteration 61: Global Deviance = 1982.597
#> GAMLSS-RS iteration 62: Global Deviance = 1981.963
#> GAMLSS-RS iteration 63: Global Deviance = 1981.381
#> GAMLSS-RS iteration 64: Global Deviance = 1980.854
#> GAMLSS-RS iteration 65: Global Deviance = 1980.384
#> GAMLSS-RS iteration 66: Global Deviance = 1979.972
#> GAMLSS-RS iteration 67: Global Deviance = 1979.619
#> GAMLSS-RS iteration 68: Global Deviance = 1979.317
#> GAMLSS-RS iteration 69: Global Deviance = 1979.062
#> GAMLSS-RS iteration 70: Global Deviance = 1978.85
#> GAMLSS-RS iteration 71: Global Deviance = 1978.676
#> GAMLSS-RS iteration 72: Global Deviance = 1978.533
#> GAMLSS-RS iteration 73: Global Deviance = 1978.417
#> GAMLSS-RS iteration 74: Global Deviance = 1978.324
#> GAMLSS-RS iteration 75: Global Deviance = 1978.25
#> GAMLSS-RS iteration 76: Global Deviance = 1978.192
#> GAMLSS-RS iteration 77: Global Deviance = 1978.146
#> GAMLSS-RS iteration 78: Global Deviance = 1978.112
#> GAMLSS-RS iteration 79: Global Deviance = 1978.085
#> GAMLSS-RS iteration 80: Global Deviance = 1978.064
#> GAMLSS-RS iteration 81: Global Deviance = 1978.047
#> GAMLSS-RS iteration 82: Global Deviance = 1978.034
#> GAMLSS-RS iteration 83: Global Deviance = 1978.024
#> GAMLSS-RS iteration 84: Global Deviance = 1978.016
#> GAMLSS-RS iteration 85: Global Deviance = 1978.01
#> GAMLSS-RS iteration 86: Global Deviance = 1978.005
#> GAMLSS-RS iteration 87: Global Deviance = 1978.001
#> GAMLSS-RS iteration 88: Global Deviance = 1977.998
#> GAMLSS-RS iteration 89: Global Deviance = 1977.996
#> GAMLSS-RS iteration 90: Global Deviance = 1977.994
#> GAMLSS-RS iteration 91: Global Deviance = 1977.992
#> GAMLSS-RS iteration 92: Global Deviance = 1977.991
#> GAMLSS-RS iteration 93: Global Deviance = 1977.99
summary(param_mm)
#> ******************************************************************
#> Family:  c("OW", "Odd Weibull")
#>
#> Call:  gamlss(formula = mice ~ 1, sigma.formula = ~1, nu.formula = ~1,
#>     family = myOW, sigma.start = 2, nu.start = 6, control = con.out)
#>
#> Fitting method: RS()
#>
#> ------------------------------------------------------------------
#> Mu link function:  log
#> Mu Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -4.87857    0.01489  -327.6   <2e-16 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> Sigma link function:  log
#> Sigma Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept)   1.8206     0.1349    13.5   <2e-16 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> Nu link function:  log
#> Nu Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept)  -0.2794     0.1640  -1.704     0.09 .
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> ------------------------------------------------------------------
#> No. of observations in the fit:  208
#> Degrees of Freedom for the fit:  3
#>       Residual Deg. of Freedom:  205
#>                       at cycle:  93
#>
#> Global Deviance:     1977.99
#>             AIC:     1983.99
#>             SBC:     1994.003
#> ******************************************************************

Then, we report the results and compare them with those in Cooray (2006)

Parameter MLE (Cooray 2006) GAMLSS
$$\mu$$ 7.61e-03 7.61e-03
$$\sigma$$ 6.2278 6.1754
$$\nu$$ 0.7495 0.7563

# References

Cooray, Kahadawala. 2006. Generalization of the Weibull distribution: The odd Weibull family.” Statistical Modelling 6 (3): 265–77. https://doi.org/10.1191/1471082X06st116oa.
———. 2015. A study of moments and likelihood estimators of the odd Weibull distribution.” Statistical Methodology 26 (September): 72–83. https://doi.org/10.1016/j.stamet.2015.03.003.
Kimball, A. W. 1960. Estimation of Mortality Intensities in Animal Experiments.” Biometrics 16 (4): 505. https://doi.org/10.2307/2527758.
Wang, F. K. 2000. A new model with bath tub-shaped failure rate using an additive Burr XII distribution.” Reliability Engineering and System Safety 70 (3): 305–12. https://doi.org/10.1016/S0951-8320(00)00066-1.