Documentation 
Fit nonlinear mixedeffects model with stochastic EM algorithm
[BETA,PSI,STATS,B]
= nlmefitsa(X,Y,GROUP,V,MODELFUN,BETA0)
[BETA,PSI,STATS,B]
= nlmefitsa(X,Y,GROUP,V,MODELFUN,BETA0,'Name',Value)
[BETA,PSI,STATS,B] = nlmefitsa(X,Y,GROUP,V,MODELFUN,BETA0) fits a nonlinear mixedeffects regression model and returns estimates of the fixed effects in BETA. By default, nlmefitsa fits a model where each model parameter is the sum of a corresponding fixed and random effect, and the covariance matrix of the random effects is diagonal, i.e., uncorrelated random effects.
The BETA, PSI, and other values this function returns are the result of a random (Monte Carlo) simulation designed to converge to the maximum likelihood estimates of the parameters. Because the results are random, it is advisable to examine the plot of simulation to results to be sure that the simulation has converged. It may also be helpful to run the function multiple times, using multiple starting values, or use the 'Replicates' parameter to perform multiple simulations.
[BETA,PSI,STATS,B] = nlmefitsa(X,Y,GROUP,V,MODELFUN,BETA0,'Name',Value) accepts one or more commaseparated parameter name/value pairs. Specify Name inside single quotes.
Definitions:
In the following list of arguments, the following variable definitions apply:
n — number of observations
h — number of predictor variables
m — number of groups
g — number of groupspecific predictor variables
p — number of parameters
f — number of fixed effects
By default, nlmefitsa fits a model where each model parameter is the sum of a corresponding fixed and random effect. Use the following parameter name/value pairs to fit a model with a different number of or dependence on fixed or random effects. Use at most one parameter name with an 'FE' prefix and one parameter name with an 'RE' prefix. Note that some choices change the way nlmefitsa calls MODELFUN, as described further below.
The default model is equivalent to setting both FEConstDesign and REConstDesign to eye(p), or to setting both FEParamsSelect and REParamsSelect to 1:p.
Additional optional parameter name/value pairs control the iterative algorithm used to maximize the likelihood:
'CovPattern' 
Specifies an rbyr logical or numeric matrix PAT that defines the pattern of the random effects covariance matrix PSI. nlmefitsa computes estimates for the variances along the diagonal of PSI as well as covariances that correspond to nonzeroes in the offdiagonal of PAT. nlmefitsa constrains the remaining covariances, i.e., those corresponding to offdiagonal zeroes in PAT, to be zero. PAT must be a rowcolumn permutation of a block diagonal matrix, and nlmefitsa adds nonzero elements to PAT as needed to produce such a pattern. The default value of PAT is eye(r), corresponding to uncorrelated random effects. Alternatively, specify PAT as a 1byr vector containing values in 1:r. In this case, elements of PAT with equal values define groups of random effects, nlmefitsa estimates covariances only within groups, and constrains covariances across groups to be zero. 
'Cov0' 
Initial value for the covariance matrix PSI. Must be an rbyr positive definite matrix. If empty, the default value depends on the values of BETA0. 
'ComputeStdErrors' 
true to compute standard errors for the coefficient estimates and store them in the output STATS structure, or false (default) to omit this computation. 
'ErrorModel' 
A string specifying the form of the error term. Default is 'constant'. Each model defines the error using a standard normal (Gaussian) variable e, the function value f, and one or two parameters a and b. Choices are
If this parameter is given, the output STATS.errorparam field has the value

'ErrorParameters' 
A scalar or twoelement vector specifying starting values for parameters of the error model. This specifies the a, b, or [a b] values depending on the ErrorModel parameter. 
'LogLikMethod' 
Specifies the method for approximating the loglikelihood. Choices are:

'NBurnIn' 
Number of initial burnin iterations during which the parameter estimates are not recomputed. Default is 5. 
'NChains' 
Number c of "chains" simulated. Default is 1. Setting c>1 causes c simulated coefficient vectors to be computed for each group during each iteration. Default depends on the data, and is chosen to provide about 100 groups across all chains. 
'NIterations' 
Number of iterations. This can be a scalar or a threeelement vector. Controls how many iterations are performed for each of three phases of the algorithm:
Default is [150 150 100]. A scalar is distributed across the three phases in the same proportions as the default. 
'NMCMCIterations' 
Number of Markov Chain Monte Carlo (MCMC) iterations. This can be a scalar or a threeelement vector. Controls how many of three different types of MCMC updates are performed during each phase of the main iteration:
Default is [2 2 2]. A scalar value is treated as a threeelement vector with all elements equal to the scalar. 
'OptimFun' 
Either 'fminsearch' or 'fminunc', specifying the optimization function to be used during the estimation process. Default is 'fminsearch'. Use of 'fminunc' requires Optimization Toolbox™. 
'Options' 
A structure created by a call to statset. nlmefitsa uses the following statset parameters:

'ParamTransform' 
A vector of pvalues specifying a transformation function f() for each of the p parameters: XB = ADESIGN*BETA + BDESIGN*B PHI = f(XB) Each element of the vector must be one of the following integer codes specifying the transformation for the corresponding value of PHI:

'Replicates' 
Number REPS of estimations to perform starting from the starting values in the vector BETA0. If BETA0 is a matrix, REPS must match the number of columns in BETA0. Default is the number of columns in BETA0. 
'Vectorization' 
Determines the possible sizes of the PHI, XFUN, and VFUN input arguments to MODELFUN. Possible values are:
The default for 'Vectorization' is 'SinglePhi'. In all cases, if V is empty, nlmefitsa calls MODELFUN with only two inputs. 
[1] Delyon, B., M. Lavielle, and E. Moulines, Convergence of a stochastic approximation version of the EM algorithm, Annals of Statistics, 27, 94128, 1999.
[2] Mentré, France, and Marc Lavielle, Stochastic EM algorithms in population PKPD analyses, 2008.