Last data update: 2014.03.03

R: Locally D-optimal designs for Weibull model
ldweibullR Documentation

Locally D-optimal designs for Weibull model

Description

Finds Locally D-optimal designs for Weibull regression model which is defined as E(y) = a-bexp(-λ*x^h) with Var(y) = σ^2, where a, b, λ, h and σ are unknown parameters.

Usage


ldweibull(a, b, lambda, h, lb, ub, user.points = NULL, user.weights = NULL, 
...,  n.restarts = 1, n.sim = 1, tol = 1e-8, prec = 53, rseed = NULL) 

Arguments

a

initial value for paremeter a.

b

initial value for paremeter b.

lambda

initial value for paremeter λ.

h

initial value for paremeter h.

lb

lower bound of design interval, must be greater than 0. Value 0 for lower bound is not allowed, instead of 0 a small value such as 10^-10 can be used.

ub

upper bound of design interval.

user.points

(optional) vector of user design points which calculation of its D-efficiency is aimed. Each element of user.points must be within the design interval.

user.weights

(optional) vector of weights which its elements correspond to user.points elements. The sum of weights should be 1; otherwise they will be normalized.

...

(optional) additional parameters will be passed to function curve.

prec

(optional) a number, the maximal precision to be used for D-efficiency calculation, in bite. Must be at least 2 (default 53), see 'Details'.

n.restarts

(optional optimization parameter) number of solver restarts required in optimization process (default 1), see 'Details'.

n.sim

(optional optimization parameter) number of random parameters to generator for every restart of solver in optimization process (default 1), see 'Details'.

tol

(optional optimization parameter) relative tolerance on feasibility and optimality in optimization process (default 1e-8).

rseed

(optional optimization parameter) a seed to initiate the random number generator, else system time will be used.

Details

While D-efficiency is NaN, an increase in prec can be beneficial to achieve a numeric value, however, it can slow down the calculation speed.

Values of n.restarts and n.sim should be chosen according to the length of design interval.

Value

plot of derivative function, see 'Note'.

a list containing the following values:

points

obtained design points

weights

corresponding weights to the obtained design points

det.value

value of Fisher information matrix determinant at the obtained design

user.eff

D-efficeincy of user design, if user.design and user.weights are not NULL.

Note

To verify optimality of obtained design, derivate function (symmetry of Frechet derivative with respect to the x-axis) will be plotted on the design interval. Based on the equivalence theorem (Kiefer, 1974), a design is optimal if and only if its derivative function are equal or less than 0 on the design interval. The equality must be achieved just at the obtained points.

Author(s)

Ehsan Masoudi, Majid Sarmad and Hooshang Talebi

References

Masoudi, E., Sarmad, M. and Talebi, H. 2012, An Almost General Code in R to Find Optimal Design, In Proceedings of the 1st ISM International Statistical Conference 2012, 292-297.

Dette, H., Pepelyshev, A. (2008), Efficient Experimental Designs for Sigmoidal Growth Models, Statistical Planning and Inference, 138, 2-17.

Kiefer, J. C. 1974, General equivalence theory for optimum designs (approximate theory), Ann. Statist., 2, 849-879.

See Also

cfisher, cfderiv and eff.

Examples

ldweibull(a = 1, b = 1, lambda = 2, h = 1, lb = 10^-10, ub =3) 
# $points: 0.0000000001 0.1713914120 0.8002692550 3.0000000000

## usage of n.sim and n.restars:
# Various responses for different rseed

ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, rseed = 1) 
# $points: 0.0010000  0.2991952  5.2428039 19.0000000

ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, rseed = 19) 
# $points: 0.001000  1.217404  3.566328 19.000000

ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, n.sim = 10, n.restarts = 10) 
# (valid respone) $points: 0.0010000, 0.1205858,  0.5544623, 19.0000000

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(LDOD)
Loading required package: Rsolnp
Loading required package: Rmpfr
Loading required package: gmp

Attaching package: 'gmp'

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

    %*%, apply, crossprod, matrix, tcrossprod

C code of R package 'Rmpfr': GMP using 64 bits per limb


Attaching package: 'Rmpfr'

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

    dbinom, dnorm, dpois, pnorm

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

    cbind, pmax, pmin, rbind

> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/LDOD/ldweibull.Rd_%03d_medium.png", width=480, height=480)
> ### Name: ldweibull
> ### Title: Locally D-optimal designs for Weibull model
> ### Aliases: ldweibull
> ### Keywords: optimal design Weibull equivalence theorem
> 
> ### ** Examples
> 
> ldweibull(a = 1, b = 1, lambda = 2, h = 1, lb = 10^-10, ub =3) 

Iter: 1 fn: 11.2682	 Pars:  0.80027 0.17139
Iter: 2 fn: 11.2682	 Pars:  0.80027 0.17139
solnp--> Completed in 2 iterations
$points
[1] 0.0000000001 0.1713913339 0.8002685655 3.0000000000

$weights
[1] 0.25 0.25 0.25 0.25

$det.value
[1] 1.277295e-05

> # $points: 0.0000000001 0.1713914120 0.8002692550 3.0000000000
> 
> ## usage of n.sim and n.restars:
> # Various responses for different rseed
> 
> ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, rseed = 1) 

Iter: 1 fn: 33.6254	 Pars:  0.29403 5.24512
Iter: 2 fn: 33.6254	 Pars:  0.29403 5.24512
solnp--> Completed in 2 iterations
$points
[1]  0.0010000  0.2940269  5.2451226 19.0000000

$weights
[1] 0.25 0.25 0.25 0.25

$det.value
[1] 2.492623e-15

> # $points: 0.0010000  0.2991952  5.2428039 19.0000000
> 
> ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, rseed = 19) 

Iter: 1 fn: 28.9843	 Pars:  1.21740 3.56633
Iter: 2 fn: 28.9843	 Pars:  1.21740 3.56633
solnp--> Completed in 2 iterations
$points
[1]  0.001000  1.217404  3.566328 19.000000

$weights
[1] 0.25 0.25 0.25 0.25

$det.value
[1] 2.583807e-13

> # $points: 0.001000  1.217404  3.566328 19.000000
> 
> ldweibull(a = 1, b = 1, lambda = 3, h = 1, lb = 0.001, ub = 19, n.sim = 10, n.restarts = 10) 

Iter: 1 fn: 12.0279	 Pars:  0.12059 0.55446
Iter: 2 fn: 12.0279	 Pars:  0.12059 0.55446
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.12059 0.55446
Iter: 2 fn: 12.0279	 Pars:  0.12059 0.55446
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.55446 0.12059
Iter: 2 fn: 12.0279	 Pars:  0.55446 0.12059
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.12058 0.55446
Iter: 2 fn: 12.0279	 Pars:  0.12059 0.55446
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.55447 0.12058
Iter: 2 fn: 12.0279	 Pars:  0.55446 0.12059
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.55446 0.12059
Iter: 2 fn: 12.0279	 Pars:  0.55446 0.12059
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.12059 0.55446
Iter: 2 fn: 12.0279	 Pars:  0.12059 0.55446
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.55446 0.12059
Iter: 2 fn: 12.0279	 Pars:  0.55446 0.12059
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.12059 0.55446
Iter: 2 fn: 12.0279	 Pars:  0.12059 0.55446
solnp--> Completed in 2 iterations

Iter: 1 fn: 12.0279	 Pars:  0.55446 0.12059
Iter: 2 fn: 12.0279	 Pars:  0.55446 0.12059
solnp--> Completed in 2 iterations
$points
[1]  0.0010000  0.1205858  0.5544623 19.0000000

$weights
[1] 0.25 0.25 0.25 0.25

$det.value
[1] 5.975097e-06

> # (valid respone) $points: 0.0010000, 0.1205858,  0.5544623, 19.0000000
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>