Red Door Analytics

The first thing to read here is the agreement, not the seconds.

On the Weibull frailty fit below, merlin and frailtypack return the same log-likelihood, the same coefficient and the same frailty standard deviation to every digit reported. On the 100,000-row Royston–Parmar fit, merlin and flexsurv return an identical log-likelihood and all three packages agree on the coefficient to four decimal places. Independent implementations, written by different people over more than a decade, converging on the same answer is the useful result. The timings are a practical matter downstream of it.

A note on the packages here

flexsurv, rstpm2, frailtypack, cuRe and flexsurvcure are long-standing, widely used packages, several of them implementing methods merlin also implements. They are the work of researchers whose papers we cite and whose software this field — ourselves included — has relied on for years. merlin was written knowing them, and in places because of them.

Where merlin is faster on these configurations, it is generally because it was written later and narrower: one compiled engine, analytic first and second derivatives, and a single likelihood path rather than a general-purpose one. That is an engineering difference on a particular fit, not a statement about the statistics, and it does not generalise to models or options these packages support and merlin does not.

Measured

Timings

Seconds are wall-clock on a single host, so read them as ratios within a row rather than as absolute performance. Log-likelihoods and coefficients are what each package returned on the same simulated data.

Royston–Parmar, 100,000 rows, 20 covariates, no random effects
PackageMethodSeconds Log-likelihoodCoefficient
merlinanalytic, no RE0.64−83895.670.3024
rstpm2stpm25.06−83896.150.3024
flexsurvspline23.78−83895.670.3024
Weibull frailty, 10,000 rows, 500 clusters
PackageMethodSeconds Log-likelihoodCoefficient
merlinAGHQ0.20−8775.480.2781
merlinLaplace0.47−8775.740.2779
frailtypackWeibull, log-normal frailty38.35−8775.480.2781
Royston–Parmar frailty, 50,000 rows, 2,000 clusters
PackageMethodSeconds Log-likelihoodCoefficient
merlinAGHQ1.66−42080.70.2913
rstpm2frailty, log-normal, Q=711.27−42080.960.2913
coxmeCox frailty5.12not comparable0.2911
Against itself

merlin for Stata, version to version

v3.0.0 rebuilt the likelihood engine with analytic gradients and Hessians. This compares it against v2.4.7 on identical scenarios, each version built from its own source. All ten scenarios are shown, not a selection: the gain ranges from almost nothing on simple fixed-effects fits to a different order of magnitude on the Cox path, which was re-implemented rather than tuned.

merlin for Stata — v3.0.0 against v2.4.7, identical scenarios
ScenarioNv2.4.7 (s) v3.0.0 (s)Same log-lik.
Weibull, 5 covariates5000.2010.178yes
Weibull, fixed effects200,0000.6990.588yes
Royston–Parmar, df 4, fixed effects200,0001.0650.932yes
Log-hazard splines (RCS)10,0001.6130.493yes
Weibull, time-dependent effect10,0001.0370.511yes
Cox, fixed effects10,00019.4610.042yes
Weibull, two-level random effects6,0000.4660.192yes
Gaussian, random intercept + slope3,0001.9310.656yes
Survival prediction, point estimate200,0000.1830.159
Survival prediction with CIs200,0006.6211.692

How these were produced, and what they do not control

  • Generated by Rscript bench/compare.R in the merlin for R repository, under R 4.6.0. The curated snapshot lives at bench/PERF-COMPARISON.md with per-section provenance and commit hashes.
  • One machine, one run per configuration. The benchmark's own notes record that the host thermally throttles under sustained load, which is why comparisons are made within a run and not across days.
  • Comparator package versions are not pinned in the published snapshot. Timings will move with future releases of any package here, including merlin.
  • Each configuration is simulated data chosen to exercise one model class. Nothing here measures a package on the models it was designed for and merlin is not.
  • coxme fits a Cox frailty model, so its log-likelihood is not comparable with the parametric fits and is marked as such.

If a number here does not reproduce for you, we would like to know — tell us.