Do-it-yourself toolkit for Lambert W \times F distribution

Description

IMPORTANT: This toolkit functionality is still under active development; function names, arguments, return values, etc. may change.

This do-it-yourself Lambert W \times F toolkit implements the flexible input/output framework of Lambert W \times F random variables (see References). Using a modular approach, it allows users to create their own Lambert W \times 'MyFavoriteDistribution' RVs. See Details below.

If the distribution you inted to use is not already implemented (get_distnames), then you can create it:

create input:

use create_LambertW_input with your favorite distribution,

create output:

pass it as an input argument to create_LambertW_output,

use output:

use Rs standard functionality for distributions such as random number generation (rY), pdf (dY) and cdf (pY), quantile function (qY), etc. for this newly generated Lambert W \times 'MyFavoriteDistribution'.

create_LambertW_output converts the input LambertW_input representing random variable X sim F_X to the Lambert W \times F_X output.

Usage

create_LambertW_input(distname = NULL, beta, input.u = list(beta2tau = NULL,
  d = NULL, p = NULL, r = NULL, q = NULL, distname = "MyFavoriteDistribution",
  is.non.negative = FALSE))

create_LambertW_output(LambertW.input = NULL, theta = NULL,
  distname = LambertW.input$distname)

Arguments

Details

create_LambertW_output takes an object of class LambertW_input and creates a class LambertW_output for standard distributions as well as the user-defined distribution. This LambertW_output represents the RV Y sim Lambert W \times 'MyFavoriteDistribution' with all its properties and R functionality, such as random number generation (rY), pdf (dY) and cdf (pY), etc.

create_LambertW_input allows users to define their own Lambert W\times F distribution by supplying the necessary functions about the input random variable U and oldsymbol β. Here U is the zero mean and/or unit variance version of X sim F_X(x mid oldsymbol β) (see References).

The argument input.u must be a list containing all of the following:

beta2tau

R function of (beta): converts oldsymbol β to τ for the user defined distribution

distname

optional; users can specify the name of their input distribution. By default it's called "MyFavoriteDistribution". The distribution name will be used in plots and summaries of the Lambert W\times F input (and output) object.

is.non.negative

logical; users should specify whether the distribution is for non-negative random variables or not. This will help for plotting and theoretical quantile computation.

d

R function of (u, beta): probability density function (pdf) of U,

p

R function of (u, beta): cumulative distribution function (cdf) of U,

q

R function of (p, beta): quantile function of U,

r

R function (n, beta): random number generator for U,

Value

create_LambertW_output returns a list of class LambertW_output with values that are (for the most part) functions themselves (see Examples):

Author(s)

Georg M. Goerg

Examples

# create a Gaussian N(1, 2) input
Gauss.input <- create_LambertW_input("normal", beta = c(1, 2))

# create a heavy-tailed version of a normal
# gamma = 0, alpha = 1 are set by default; beta comes from input
params <- list(delta = c(0.3)) 
LW.Gauss <- create_LambertW_output(LambertW.input = Gauss.input, 
                                   theta = params)
LW.Gauss

op <- par(no.readonly = TRUE)
par(mfrow = c(2, 1), mar = c(3, 3, 2, 1))
curve(LW.Gauss$d(x, params), -7, 10, col = "red")
# parameter will get detected automatically from the input
curve(LW.Gauss$d(x), -7, 10, col = "blue") # same in blue; 

# compare to the input case (i.e. set delta = 0)
params.0 <- params 
params.0$delta <- 0

# to evaluate the RV at a different parameter value, 
# it is necessary to pass the new parameter
curve(LW.Gauss$d(x, params.0), -7, 10, add = TRUE, col = 1) #' par(op)

curve(LW.Gauss$p(x, params), -7, 10, col = "red")
curve(LW.Gauss$p(x, params.0), -7, 10, add = TRUE, col = 1)

test_normality(LW.Gauss$r(n = 100), add.legend = FALSE)

## generate a positively skewed version of a shifted, scaled t_3
t.input <- create_LambertW_input("t", beta = c(2, 1, 3))
t.input
params <- list(gamma = 0.05) # skew it
LW.t <- create_LambertW_output(LambertW.input = t.input, theta = params)
LW.t

plot(t.input$d, -7, 11, col = 1)
plot(LW.t$d, -7, 11, col = 2, add = TRUE)
abline(v = t.input$beta["location"], lty = 2)

# draw samples from the skewed t_3
yy <- LW.t$r(n = 100)
test_normality(yy)

### create a skewed exponential distribution
exp.input <- create_LambertW_input("exp", beta = 1)
plot(exp.input)
params <- list(gamma = 0.2)
LW.exp <- create_LambertW_output(exp.input, theta = params)
plot(LW.exp)

# create a heavy-tail exponential distribution
params <- list(delta = 0.2)
LW.exp <- create_LambertW_output(exp.input, theta = params)
plot(LW.exp)

# create a skewed chi-square distribution with 5 df
chi.input <- create_LambertW_input("chisq", beta = 5)
plot(chi.input)
params <- list(gamma = sqrt(2)*0.2)
LW.chi <- create_LambertW_output(chi.input, theta = params)
plot(LW.chi)


# a demo on how a user-defined U input needs to look like
user.tmp <- list(d = function(u, beta) dnorm(u),
                 r = function(n, beta) rnorm(n),
                 p = function(u, beta) pnorm(u),
                 q = function(p, beta) qnorm(p),
                 beta2tau = function(beta) {
                   c(mu_x = beta[1], sigma_x = beta[2], 
                     gamma = 0, alpha = 1, delta = 0)
                   },
                 distname = "MyNormal",
                 is.non.negative = FALSE)
my.input <- create_LambertW_input(input.u = user.tmp, beta = c(0, 1))
my.input
plot(my.input)

Results

R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

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(LambertW)
Loading required package: MASS
Loading required package: ggplot2
This is 'LambertW' version 0.6.4.  Please see the NEWS file and citation("LambertW").

> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/LambertW/LambertW-toolkit.Rd_%03d_medium.png", width=480, height=480)
> ### Name: LambertW-toolkit
> ### Title: Do-it-yourself toolkit for Lambert W \times F distribution
> ### Aliases: LambertW-toolkit create_LambertW_input create_LambertW_output
> ### Keywords: datagen distribution models univar
> 
> ### ** Examples
> 
> 
> # create a Gaussian N(1, 2) input
> Gauss.input <- create_LambertW_input("normal", beta = c(1, 2))
> 
> # create a heavy-tailed version of a normal
> # gamma = 0, alpha = 1 are set by default; beta comes from input
> params <- list(delta = c(0.3)) 
> LW.Gauss <- create_LambertW_output(LambertW.input = Gauss.input, 
+                                    theta = params)
> LW.Gauss
 Input distribution: normal
 Output distribution: heavy-tail (one parameter) Lambert W x normal(1,2)
 with (input) parameters: mu=1, sigma=2
 and transformation parameters: delta = 0.3
> 
> op <- par(no.readonly = TRUE)
> par(mfrow = c(2, 1), mar = c(3, 3, 2, 1))
> curve(LW.Gauss$d(x, params), -7, 10, col = "red")
> # parameter will get detected automatically from the input
> curve(LW.Gauss$d(x), -7, 10, col = "blue") # same in blue; 
> 
> # compare to the input case (i.e. set delta = 0)
> params.0 <- params 
> params.0$delta <- 0
> 
> # to evaluate the RV at a different parameter value, 
> # it is necessary to pass the new parameter
> curve(LW.Gauss$d(x, params.0), -7, 10, add = TRUE, col = 1) #' par(op)
> 
> curve(LW.Gauss$p(x, params), -7, 10, col = "red")
> curve(LW.Gauss$p(x, params.0), -7, 10, add = TRUE, col = 1)
> 
> test_normality(LW.Gauss$r(n = 100), add.legend = FALSE)
$seed
[1] 619124

$shapiro.wilk

	Shapiro-Wilk normality test

data:  data.test
W = 0.87957, p-value = 1.74e-07


$shapiro.francia

	Shapiro-Francia normality test

data:  data.test
W = 0.87304, p-value = 6.609e-07


$anderson.darling

	Anderson-Darling normality test

data:  data
A = 2.31, p-value = 6.824e-06


> 
> ## generate a positively skewed version of a shifted, scaled t_3
> t.input <- create_LambertW_input("t", beta = c(2, 1, 3))
> t.input
 Input distribution: t
 with parameters: location=2, scale=1, df=3
> params <- list(gamma = 0.05) # skew it
> LW.t <- create_LambertW_output(LambertW.input = t.input, theta = params)
> LW.t
 Input distribution: t
 Output distribution: skewed Lambert W x t(2,1,3)
 with (input) parameters: location=2, scale=1, df=3
 and transformation parameters: gamma=0.05
> 
> plot(t.input$d, -7, 11, col = 1)
> plot(LW.t$d, -7, 11, col = 2, add = TRUE)
> abline(v = t.input$beta["location"], lty = 2)
> 
> # draw samples from the skewed t_3
> yy <- LW.t$r(n = 100)
> test_normality(yy)
$seed
[1] 671065

$shapiro.wilk

	Shapiro-Wilk normality test

data:  data.test
W = 0.86909, p-value = 6.441e-08


$shapiro.francia

	Shapiro-Francia normality test

data:  data.test
W = 0.85404, p-value = 1.593e-07


$anderson.darling

	Anderson-Darling normality test

data:  data
A = 2.7087, p-value = 7.149e-07


> 
> ### create a skewed exponential distribution
> exp.input <- create_LambertW_input("exp", beta = 1)
> plot(exp.input)
> params <- list(gamma = 0.2)
> LW.exp <- create_LambertW_output(exp.input, theta = params)
> plot(LW.exp)
> 
> # create a heavy-tail exponential distribution
> params <- list(delta = 0.2)
> LW.exp <- create_LambertW_output(exp.input, theta = params)
> plot(LW.exp)
> 
> # create a skewed chi-square distribution with 5 df
> chi.input <- create_LambertW_input("chisq", beta = 5)
> plot(chi.input)
> params <- list(gamma = sqrt(2)*0.2)
> LW.chi <- create_LambertW_output(chi.input, theta = params)
> plot(LW.chi)
> 
> 
> # a demo on how a user-defined U input needs to look like
> user.tmp <- list(d = function(u, beta) dnorm(u),
+                  r = function(n, beta) rnorm(n),
+                  p = function(u, beta) pnorm(u),
+                  q = function(p, beta) qnorm(p),
+                  beta2tau = function(beta) {
+                    c(mu_x = beta[1], sigma_x = beta[2], 
+                      gamma = 0, alpha = 1, delta = 0)
+                    },
+                  distname = "MyNormal",
+                  is.non.negative = FALSE)
> my.input <- create_LambertW_input(input.u = user.tmp, beta = c(0, 1))
> my.input
 Note: This is a user-defined Lambert W x F distribution. 
 * * * * * * * * 
 Input distribution: MyNormal
 with parameters: =0, =1
> plot(my.input)
> 
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>