Competing risks in merlin
Two cause-specific models in one command — the second set of brackets — and the arithmetic that shows they are the same two models you would have fitted separately.
Fitting your first model in merlin fitted one model and checked it against streg, which left the syntax looking like a lot of brackets for no return. This is where the return starts: a second set of brackets is a second model, fitted at the same time as the first.
The data
hypoxia.dta is the competing-risks dataset from Stata's own documentation: 109 women treated for cervical cancer, followed to the first recurrence of disease. A recurrence is either pelvic or distant, never both, and most patients had neither.
webuse hypoxia, clear
gen byte pelvic = failtype==1
gen byte distant = failtype==2. tab failtype
Failure |
type: 1 if |
pelrec, 2 |
if disrec & |
not pelrec, |
0 otherwise | Freq. Percent Cum.
------------+-----------------------------------
0 | 59 54.13 54.13
1 | 33 30.28 84.40
2 | 17 15.60 100.00
------------+-----------------------------------
Total | 109 100.00So 33 pelvic recurrences, 17 distant, and 59 patients censored without either. dftime is the time to whichever came first.
One cause at a time
A cause-specific hazard model asks about one cause and treats the competing cause as censoring: a patient whose disease came back distantly contributes follow-up time to the pelvic model, and then leaves it without a pelvic event. That is what failure(pelvic) does here, because pelvic is 0 for those patients.
. merlin (dftime ifp tumsize, family(weibull, failure(pelvic)))
Fitting full model:
Iteration 0: Log likelihood = -353.129
Iteration 1: Log likelihood = -68.116343
Iteration 2: Log likelihood = -47.721332
Iteration 3: Log likelihood = -41.385181
Iteration 4: Log likelihood = -41.374054
Iteration 5: Log likelihood = -41.37405
Fixed effects regression model Number of obs = 109
Log likelihood = -41.37405
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
dftime: |
ifp | .0392511 .0197033 1.99 0.046 .0006334 .0778689
tumsize | .4020289 .0893376 4.50 0.000 .2269304 .5771274
_cons | -4.569662 .7529348 -6.07 0.000 -6.045387 -3.093937
log(gamma) | -1.228312 .155388 -7.90 0.000 -1.532867 -.9237574
------------------------------------------------------------------------------Both causes in one command
Now the second set of brackets. Each pair carries its own response, its own covariates and its own family, and merlin fits them together.
. merlin (dftime ifp tumsize, family(weibull, failure(pelvic))) ///
> (dftime ifp tumsize, family(weibull, failure(distant)))
Fitting full model:
Iteration 0: Log likelihood = -706.258
Iteration 1: Log likelihood = -136.30085
Iteration 2: Log likelihood = -112.17419
Iteration 3: Log likelihood = -105.4927
Iteration 4: Log likelihood = -105.46881
Iteration 5: Log likelihood = -105.46879
Iteration 6: Log likelihood = -105.46879
Fixed effects regression model Number of obs = 109
Log likelihood = -105.46879
------------------------------------------------------------------------------
| Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
dftime: |
ifp | .0392511 .0197033 1.99 0.046 .0006334 .0778689
tumsize | .4020289 .0893376 4.50 0.000 .2269304 .5771274
_cons | -4.569663 .7529349 -6.07 0.000 -6.045388 -3.093937
log(gamma) | -1.228312 .155388 -7.90 0.000 -1.532867 -.9237574
-------------+----------------------------------------------------------------
dftime: |
ifp | .0766309 .0274645 2.79 0.005 .0228015 .1304602
tumsize | .0910814 .1507554 0.60 0.546 -.2043938 .3865565
_cons | -4.870947 1.152572 -4.23 0.000 -7.129947 -2.611946
log(gamma) | -.1708737 .2106839 -0.81 0.417 -.5838066 .2420592
------------------------------------------------------------------------------The two blocks appear in the order you wrote them. Both are labelled dftime, because merlin labels an equation by its response variable and both causes are measured on the same clock — the failure indicator is what separates them.
Is it really the same model?
It should be. With nothing shared between the two sets of brackets, the joint likelihood is the product of the two separate ones, so the log-likelihoods should add up.
pelvic alone = -41.374050
distant alone = -64.094744
their sum = -105.468794
both together = -105.468794They do. Across all eight parameters the largest disagreement between the joint fit and the two separate fits is 3.5e-07, which is the optimiser's convergence tolerance rather than a difference between the models. Fitting the two causes together has not changed either of them.
Then why write it this way?
Because nothing is shared yet. The moment the two causes share a parameter — the same coefficient constrained across both, or a random effect that carries a patient's frailty into each — the joint likelihood stops factorising and two separate commands can no longer produce it. The second set of brackets is what makes that possible to write down; this tutorial is the check that it costs nothing when you do not need it.
MethodCompeting risks