*! version 0.1.0 ?????2017

Title

merlin -- user-defined Mata functions

Syntax

Here I describe some of the advanced capabilities of merlin, and hence this is where it is easiest to go wrong.
If this is your first time writing a user-defined function with merlin, then please read the
Description first.

The easiest way to learn the full capabilities is to work through some of the more complex examples, described at
reddooranalytics.se/products/merlin.

Utility functions

real matrix merlin_util_depvar(gml)
returns the dependent variable(s) from the model. If the model is a survival model, with failure(varname) specified, then
the returned matrix has a second column which contains the event indicator variable.

real matrix merlin_util_xzb(gml [, t])
returns the complex linear predictor from the model,
i.e. the XB[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details.

real matrix merlin_util_xzb_deriv(gml [, t])
returns the first derivative with respect to time, d/dt, of the complex linear predictor from the model,
i.e. the dXB[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical differentiation.

real matrix merlin_util_xzb_deriv2(gml [, t])
returns the second derivative with respect to time, d2/dt2, of the complex linear predictor from the model,
i.e. the d2XB[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical differentiation.

real matrix merlin_util_xzb_integ(gml [, t])
returns the integral with respect to time of the complex linear predictor from the model,
i.e. the iXB[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical integration.

real matrix merlin_util_expval(gml [, t])
returns the expected value of the current model, i.e. the EV[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details.

real matrix merlin_util_expval_deriv(gml [, t])
returns the first derivative with respect to time, d/dt, of the expected value of the current model,
i.e. the dEV[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical differentiation.

real matrix merlin_util_xzb_deriv2(gml [, t])
returns the second derivative with respect to time, d2/dt2, of the expected value of the current model,
i.e. the d2EV[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical differentiation.

real matrix merlin_util_xzb_integ(gml [, t])
returns the integral with respect to time of the expected value of the current model,
i.e. the iEV[] element type. If t is specified, it is evaluated at time points contained in t. See
Passing time for details. It is calculated using numerical integration.

real scalar merlin_util_ap(gml, i) returns the ith ancillary
parameter of the current model. Ancillary parameters are specified using the nap(#) option.

real scalar merlin_util_timevar(gml) returns the time variable
of the current model. This corresponds to what was specified in the associated timevar() option.

Indexing other models

All utility functions described above have a corresponding *_mod() function to do the same as *, but
extract the information for another model.
They can be used to build models which depend on other model's linear predictors, for example.

Models are indexed by the order they are specified in the merlin command.

Definition of gml

gml, if it is declared, should be declared transmorphic. gml holds the information about the merlin
optimization problem. Think of it as an M object when writing an evaluator using moptimize.
You should not attempt to alter the contents of gml.

Passing time

If you are passing a user-written function to the hfunction() argument of family(user,...), without
specifying the chfunction(), then you must make sure that your function takes a second argument of time t.
This must then be passed to any other utility functions which may depend on time, and also used as your time variable in the
definition of your hazard function. This is because time will be numerically integrated over in the calculation of the
cumulative hazard function, required in the log likelihood, so merlin must know everything that is dependent on time.

Description

The ability to have user-defined log-likelihood, hazard or cumulative hazard functions, provides substantial flexibility in the models
you can fit with merlin. This functionality means you get the extended, complex linear predictor of an inbuilt merlin
model, for free. You also have the entire library of Mata functions at your disposal, to write as simple or as complex a function
as you so wish. You can then combine your new model with any of the inbuilt models within a multivariate framework, or indeed combine it
with another user-defined model...the possibilities are endless.

For full worked examples, please visit mjcrowther.co.uk/software/merlin.

Writing a merlin log-likelihood function

The structure of a user-defined function, to pass to family(user, llfunction(functionname)} is as follows

real matrix functionname(gml)

{

...

some code

...

logl = some more code

return(logl)

}

Writing a merlin hazard function

The first available structure of a user-defined function, to pass to family(user, hfunction(functionname)} is as follows

real matrix functionname(gml)

{

...

some code

...

haz = some more code

return(haz)

}

which assumes you have also specified the chfunction() option. If you haven't, then you need to use the second structure, as
follows:

real matrix functionname(gml , real colvector t)

{

...

some code

...

haz = some more code

return(haz)

}

Writing a merlin cumulative hazard function

The structure of a user-defined function, to pass to family(user, chfunction(functionname)} is as follows:

real matrix functionname(gml)

{

...

some code

...

chaz = some more code

return(chaz)

}

Examples

For detailed examples, see reddooranalytics.se/products/merlin.

Author

Michael J. Crowther

Red Door Analytics

Stockholm, Sweden

michael@reddooranalytics.se

References

Crowther MJ. Extended multivariate generalised linear and non-linear mixed effects models.
https://arxiv.org/abs/1710.02223

Crowther MJ. merlin - a unified framework for data analysis and methods development in Stata. Stata Journal 2020;20(4):763-784.

Crowther MJ. Multilevel mixed effects parametric survival analysis: Estimation, simulation and application. Stata Journal 2019;19(4):931-949.

Crowther MJ, Lambert PC. Parametric multi-state survival models: flexible modelling allowing transition-specific distributions with
application to estimating clinically useful measures of effect differences. Statistics in Medicine 2017;36(29):4719-4742.

Weibull CE, Lambert PC, Eloranta S, Andersson TM-L, Dickman PW, Crowther MJ. A multi-state model incorporating
estimation of excess hazards and multiple time scales. Statistics in Medicine 2021;40(9):2139-2154.