gnlrem: Generalized Nonlinear Mixed Effects Regression with TWO...

In swihart/gnlrim: Generalized Non-Linear Random Intercept Models

Generalized Nonlinear Mixed Effects Regression with TWO Random Parameters

Description

gnlrem fits user-specified nonlinear regression equations to one or both parameters of the common one and two parameter distributions. One parameter of the location regression is random with some specified mixing distribution.

Usage

gnlrem(
  y = NULL,
  distribution = "normal",
  mixture = "normal-var",
  random = NULL,
  nest = NULL,
  random2 = NULL,
  nest2 = NULL,
  mu = NULL,
  shape = NULL,
  linear = NULL,
  pmu = NULL,
  pshape = NULL,
  pmix = NULL,
  delta = 1,
  common = FALSE,
  envir = parent.frame(),
  print.level = 0,
  typsize = abs(p),
  ndigit = 10,
  gradtol = 1e-05,
  stepmax = 10 * sqrt(p %*% p),
  steptol = 1e-05,
  iterlim = 100,
  compute_hessian = TRUE,
  compute_kkt = TRUE,
  fscale = 1,
  eps = 1e-04,
  trace = 0,
  maxfun.bobyqa = 10000,
  npt.bobyqa = min(p * 2, p + 2),
  abs.tol.nlminb = 1e-20,
  xf.tol.nlminb = 2.2e-14,
  x.tol.nlminb = 1.5e-08,
  rel.tol.nlminb = 1e-10,
  method = "nlminb",
  ooo = FALSE,
  p_lowb = -Inf,
  p_uppb = Inf,
  points = 5,
  steps = 10,
  int2dmethod = c("romberg_int2d", "cuba")[1],
  tol.pcubature = 1e-05,
  tol.hcubature = 0.001
)

Arguments

y	A response vector of uncensored data, a two column matrix for binomial data, or an object of class, response (created by restovec) or repeated (created by rmna or lvna). If the repeated data object contains more than one response variable, give that object in envir and give the name of the response variable to be used here.
distribution	The distribution for the response: binomial, beta binomial, (removed: double binomial, use repeated::gnlmix()), (removed: mult(iplicative) binomial, use repeated::gnlmix()), Poisson, negative binomial, (removed: double Poisson, use repeated::gnlmix()), (removed: mult(iplicative) Poisson, use repeated::gnlmix()), gamma count, Consul generalized Poisson, logarithmic series, geometric, normal, inverse Gauss, logistic, exponential, gamma, Weibull, extreme value, Cauchy, Pareto, Laplace, Levy, beta, simplex, or two-sided power. (For definitions of distributions, see the corresponding [dpqr]distribution help.) "beta-binomial-TBA" is experimental (see blogpost).
mixture	The mixing distribution for the random parameter: logit-bridge-var, logit-bridge-phi, normal-var (default), normal-phi, Cauchy-scl, Cauchy-phi, stabledist-subgauss-scl, stabledist-subgauss-phi, libstable4u-subgauss-scl, libstable4u-subgauss-phi, libstable4u-subgauss-scl-over-sqrt, libstable4u-subgauss-scl-over-sqrt-phi, logistic, Laplace, inverse Gauss, gamma, inverse gamma, Weibull, beta, simplex, or two-sided power. The first twelve have zero location parameter, the next three have unit location parameter, and the last two have location parameter set to 0.5. cloglog-bridge-delta-free estimates delta (experimental). cloglog-bridge-0-mean makes delta whatever it needs to be to maintain a mean-0 random intercept distribution (experimental). cloglog-bridge-delta-eq-alpha sets delta = alpha for the classic distribution (experimental). betaprime has mean (shp/(shp-1) for shp>1. "beta-HGLM" is experimental and is the same as "beta" but instead it is integrated with intb instead of int1. bivariate-normal-indep is for the 2 random parameters and has correlation fixed at 0. bivariate-normal-corr is for the 2 random parameters and allows correlation between the two random parameters. univariate-normal-corr is for the 2 random parameters and allows correlation between the two random parameters but uses a univariate dnorm for fitting.
random	The name of the random parameter in the mu formula.
nest	The variable classifying observations by the unit upon which they were observed. Ignored if y or envir has class, response or repeated.
random2	The name of the 2nd random parameter in the mu formula.
nest2	The variable classifying observations by the unit corresponding to 'random2'.
mu	A user-specified formula containing named unknown parameters, giving the regression equation for the location parameter. This may contain the keyword, linear referring to a linear part.
shape	A user-specified formula containing named unknown parameters, giving the regression equation for the shape parameter. This may contain the keyword, linear referring to a linear part. If nothing is supplied, this parameter is taken to be constant. This parameter is the logarithm of the usual one.
linear	A formula beginning with ~ in W&R notation, specifying the linear part of the regression function for the location parameter or list of two such expressions for the location and/or shape parameters.
pmu	Vector of initial estimates for the location parameters. These must be supplied in their order of appearance in the formula as a named list.
pshape	Vector of initial estimates for the shape parameters. These must be supplied either in their order of appearance in the expression or in a named list.
pmix	Initial estimate for the untransformed (not logarithm) of the dispersion parameter of the mixing distribution. For mixture="normal-var" and mixture="logit-bridge-var" this is the variance. For mixture="normal-phi", mixture="Cauchy-phi", mixture="stabledist-subgauss-phi", mixture="libstable4u-subgauss-phi", and mixture="logit-bridge-phi" this is the attenuation factor phi. Otherwise it is the scale. The last element must be the scale of the random intercept ('mixture') distribution, if a parameter from mu is needed for the mixture distribution then it must be listed in 'pmix' as well. cloglog-bridge-delta-free pmix will be delta parameter in the LogPS formulation. cloglog-bridge-0-mean pmix will be the variance. cloglog-bridge-delta-eq-alpha pmix will be alpha between 0 and 1. bivariate-normal-indep pmix will have the variance for random and random2. bivariate-normal-corr pmix will have the variance for the random paraters and corr12 as the 3rd element
delta	Scalar or vector giving the unit of measurement (always one for discrete data) for each response value, set to unity by default. For example, if a response is measured to two decimals, delta=0.01. If the response is transformed, this must be multiplied by the Jacobian. The transformation cannot contain unknown parameters. For example, with a log transformation, delta=1/y. (The delta values for the censored response are ignored.)
common	If TRUE, the formulae with unknowns for the location and shape have names in common. All parameter estimates must be supplied in pmu.
envir	Environment in which model formulae are to be interpreted or a data object of class, repeated, tccov, or tvcov; the name of the response variable should be given in y. If y has class repeated, it is used as the environment.
print.level	Arguments for nlm.
typsize	Arguments for nlm.
ndigit	Arguments for nlm.
gradtol	Arguments for nlm.
stepmax	Arguments for nlm.
steptol	Arguments for nlm.
iterlim	Arguments for optimx (itnmax).
compute_hessian	Argument (logical) hessian for optimx. Defaults TRUE. FALSE should only be when ooo=TRUE. To optimize with method=="nlminb" without computing the hessian, do so with the following settings: ooo=TRUE, compute_hessian=FALSE, compute_kkt=FALSE.
compute_kkt	Argument (logical) kkt for optimx (control()). Defaults TRUE. FALSE should only be when ooo=TRUE. To optimize with method=="nlminb" without computing the hessian, do so with the following settings: ooo=TRUE, compute_hessian=FALSE, compute_kkt=FALSE.
fscale	Arguments for nlm.
eps	Precision of the Romberg integration.
trace	Arguments for nlminb.
maxfun.bobyqa	Argument for bobyqa
npt.bobyqa	Argument for bobyqa
abs.tol.nlminb	Argument for nlminb. Default 1e-20 assumes known nonnegative function. Set to 0 if you don't know if function is nonnegative.
xf.tol.nlminb	Argument for nlminb. Default is 2.2e-14
x.tol.nlminb	Argument for nlminb. Default is 1.5e-8
rel.tol.nlminb	Argument for nlminb. Default is 1e-10
method	Arguments for optimx – a string denoting which optimization to use. Accommodate vector of strings, however traditional repeated::gnlmix() output will be computed only for first element if argument ooo=FALSE. Defaults to "nlminb", whereas original gnlmix used "nlm".
ooo	*o*ptimx *o*utput *o*nly. Default is FALSE, which means that the traditional repeated::gnlmix() output for the first element of method will be returned. If TRUE, then only optmix output for the elements of method are returned.
p_lowb	Argument to nlminb / optimx for the method's that are constrained optimization. Lower bounds on parameters. Defaults to -Inf. Can be a vector.
p_uppb	Argument to nlminb / optimx for the method's that are constrained optimization. Upper bounds on parameters. Defaults to Inf. Can be a vector.
points	For the Romberg integration, the number of extrapolation points so that 2*points is the order of integration, by default set to 5; points=2 is Simpson's rule.
steps	For the Romberg integration, the maximum number of steps, by default set to 10.
int2dmethod	EXPERIMENTAL. default "romberg_int2d". Can also be "cuba".
tol.pcubature	default value 1e-5. For the case of mixture "bivariate-subgauss-corr", "bivariate-t-corr", "bivariate-t-varcorr") combined with int2dmethod=="cuba" this is the tolerance that is passed to cuba::pcubature to control the precision of the likelihood. For the case of mixture = "bivariate-subgauss-corr" or "bivariate-t-corr" or "bivariate-t-varcorr" combined with int2dmethod=="cuba", this value may need to be as big as 0.1.
tol.hcubature	default value 0.001. This is the maximum tolerance passed to hcubature which does the integral for the bivariate subgaussian for mixture = "bivariate-subgauss-corr" and the bivariate-t mixture = "bivariate-t-corr" and "bivariate-t-varcorr".

Details

It is recommended that initial estimates for pmu and pshape be obtained from gnlr.

These nonlinear regression models must be supplied as formulae where parameters are unknowns. (See finterp.)

Value

If ooo=TRUE, A list of class gnlm is returned that contains all of the relevant information calculated, including error codes. If ooo=FALSE, then just the optimx output.

Author(s)

Bruce Swihart and J.K. Lindsey

Examples

dose <- c(9,12,4,9,11,10,2,11,12,9,9,9,4,9,11,9,14,7,9,8)
y <- c(8.674419, 11.506066, 11.386742, 27.414532, 12.135699,  4.359469,
      1.900681, 17.425948,  4.503345,  2.691792,  5.731100, 10.534971,
      11.220260,  6.968932,  4.094357, 16.393806, 14.656584,  8.786133,
      20.972267, 17.178012)
id <- rep(1:4, each=5)

gnlrem(y, mu=~a+b*dose+rand, random="rand", nest=id, pmu=c(a=8.7,b=0.25),
       pshape=3.44, pmix=2.3)

## Not run: 
## from repeated::gnlmix
dose <- c(9,12,4,9,11,10,2,11,12,9,9,9,4,9,11,9,14,7,9,8)
#y <- rgamma(20,shape=2+0.3*dose,scale=2)+rep(rnorm(4,0,4),rep(5,4))
y <- c(8.674419, 11.506066, 11.386742, 27.414532, 12.135699,  4.359469,
       1.900681, 17.425948,  4.503345,  2.691792,  5.731100, 10.534971,
      11.220260,  6.968932,  4.094357, 16.393806, 14.656584,  8.786133,
      20.972267, 17.178012)
resp <- restovec(matrix(y, nrow=4, byrow=TRUE), name="y")
reps <- rmna(resp, tvcov=tvctomat(matrix(dose, nrow=4, byrow=TRUE), name="dose"))

# same linear normal model with random normal intercept fitted four ways
# compare with growth::elliptic(reps, model=~dose, preg=c(0,0.6), pre=4)
glmm(y~dose, nest=individuals, data=reps)
gnlmm(reps, mu=~dose, pmu=c(8.7,0.25), psh=3.5, psd=3)
gnlmix(reps, mu=~a+b*dose+rand, random="rand", pmu=c(8.7,0.25),
	pshape=3.44, pmix=2.3)

# gamma model with log link and random normal intercept fitted three ways
glmm(y~dose, family=Gamma(link=log), nest=individuals, data=reps, points=8)
gnlmm(reps, distribution="gamma", mu=~exp(a+b*dose), pmu=c(2,0.03),
	psh=1, psd=0.3)
gnlmix(reps, distribution="gamma", mu=~exp(a+b*dose+rand), random="rand",
	pmu=c(2,0.04), pshape=1, pmix=-2)

# gamma model with log link and random gamma mixtures
gnlmix(reps, distribution="gamma", mixture="gamma",
	mu=~exp(a*rand+b*dose), random="rand", pmu=c(2,0.04),
	pshape=1.24, pmix=3.5)
gnlmix(reps, distribution="gamma", mixture="gamma",
	mu=~exp(a+b*dose)*rand, random="rand", pmu=c(2,0.04),
	pshape=1.24, pmix=2.5)

## End(Not run)

swihart/gnlrim documentation built on Oct. 18, 2023, 8:29 p.m.