Fitting, prediction & simulation
Fitting and inference
A model is fitted from data with fit(Model, @formula(Surv(time, status) ~ ...), df) (or the two-formula variant for GeneralHazard); see Illustrative example for complete calls. The optimizer is seeded automatically from the data.
Fit statistics
GeneralHazardModel <: StatsAPI.StatisticalModel, so a fitted model supports the standard statistical-model accessors. aic, aicc, and bic follow from loglikelihood, dof, and nobs; stderror follows from vcov. coef and vcov are reported on the inference scale [log.(baseline parameters); active regression coefficients], so MvNormal(coef(m), vcov(m)) is a coherent parameter-uncertainty distribution.
For a fitted-model report, coeftable(m) gives the Wald table of covariate effects (with exp(coef) and its confidence interval) and confint(m) the coefficient intervals as a DataFrame; both exclude the baseline distribution parameters, which show reports as the fitted baseline instead.
StatsAPI.loglikelihood — Method
loglikelihood(m::GeneralHazardModel)Maximized log-likelihood of the fitted model, captured from the optimizer at fit time (NaN for models built by the direct, non-optimizing constructor).
StatsAPI.nobs — Method
nobs(m::GeneralHazardModel)Number of observations (event/censoring times) the model was fit to.
StatsAPI.dof — Method
dof(m::GeneralHazardModel)Number of free parameters consumed by the fit: the baseline distribution's parameters plus the regression coefficients that actually enter the hazard for this model's structure (β for PH/AFT, α for AH, both for GH).
StatsAPI.coef — Method
coef(m::GeneralHazardModel)Identified parameters on the scale used for inference: [log.(baseline parameters); active regression coefficients], where the active coefficients are β (PH/AFT), α (AH), or [α; β] (GH). This ordering matches vcov, so MvNormal(coef(m), vcov(m)) is a coherent sampling distribution for the parameter uncertainty (exponentiate the baseline entries to recover the natural-scale baseline parameters).
StatsAPI.vcov — Method
vcov(m::GeneralHazardModel)Covariance of coef: the inverse observed information, i.e. the inverse of the Hessian of the fitted negative log-likelihood at the optimum, restricted to the identified parameters (the inactive coefficient block makes the full Hessian singular). Computed on demand — LBFGS never forms the Hessian during the fit, so there is nothing cached to reuse. stderror follows from this via the generic StatisticalModel method.
StatsAPI.coeftable — Method
coeftable(m::GeneralHazardModel; level::Real=0.95)Wald coefficient table for the covariate effects: estimate, standard error, z, p-value, exp(coef), and the confidence interval for exp(coef) at level level. exp(coef) is a hazard ratio for hazard-level effects and a time/acceleration ratio for time-scale effects. For a General Hazard model the rows are tagged with their role (hazard-level vs time-scale), since a covariate may enter both. Baseline distribution parameters are reported by show, not here: a null of zero is not the relevant hypothesis for a baseline shape/scale.
StatsAPI.confint — Method
confint(m::GeneralHazardModel; level::Real=0.95)Wald confidence intervals for the covariate coefficients at confidence level level. Returns a DataFrame with columns component (the coefficient's role, "hazard-level" or "time-scale"), term, lower, and upper. Baseline distribution parameters are not included (see coeftable).
Brier score
Inverse-probability-of-censoring-weighted Brier score (Graf et al. 1999) and its integrated form work for GeneralHazardModel through the same brier_score(model, ...) / integrated_brier_score(model, ...) API used for Cox. See the Model Evaluation: Brier Score section of the Cox documentation for the mathematical definition and signature list.
Prediction
Once a GeneralHazardModel is fitted (or directly constructed), you can evaluate per-subject cumulative hazards and survival probabilities at user-supplied times. The four hazard structures share the same closed-form expression via the unified representation
\[H(t \,|\, x) = H_0\!\left(t\, c_1(x)\right) c_2(x)\]
where $H_0$ is the cumulative hazard of the baseline distribution and $(c_1, c_2)$ are the method-specific time- and hazard-scale multipliers (c1/c2 in the code). The survival is $S(t \,|\, x) = \exp(-H(t \,|\, x))$.
predict(model, :survival) # length-n vector, each subject at own Tᵢ
predict(model, :expected) # length-n vector of Λᵢ(Tᵢ)
predict(model, :survival, t) # length-n vector at scalar t
predict(model, :expected, t)
predict(model, :survival, ts) # n × length(ts) matrix
predict(model, :expected, ts)The default no-arg form (predict(model) or predict(model, :survival)) evaluates each subject at their own observed time $T_i$, matching the convention used by the Cox interface.
Predict on new data
Each prediction also accepts a newdata::DataFrame argument. The fit's stored formula(s) are re-applied to newdata to rebuild the design matrices $X_1$, $X_2$, so newdata must contain every predictor column referenced in the original @formula(...). For models fit with fit(GHM, formula, df) (one formula), the same formula is stored twice and used for both $X_1$ and $X_2$; for fit(GeneralHazard, formula1, formula2, df) the two are stored separately.
predict(model, :survival, newdata, t) # length-n_new at scalar t
predict(model, :expected, newdata, t)
predict(model, :survival, newdata, ts) # n_new × length(ts) matrix
predict(model, :expected, newdata, ts)Newdata predict requires an explicit time argument — there is no "own time" default for arbitrary new subjects. Models built directly via the positional constructor (without the formula1 / formula2 keyword arguments) do not have stored formulas and will error on newdata predict.
SurvivalModels.predict_expected — Method
predict_expected(m::GeneralHazardModel)
predict_expected(m::GeneralHazardModel, t::Real)
predict_expected(m::GeneralHazardModel, ts::AbstractVector)
predict_expected(m::GeneralHazardModel, newdata::DataFrame, t::Real)
predict_expected(m::GeneralHazardModel, newdata::DataFrame, ts::AbstractVector)Per-subject cumulative hazard $\Lambda_i(t) = H_0(t\, c_{1i})\, c_{2i}$, where $H_0$ is the cumulative hazard of the baseline distribution and $(c_{1i}, c_{2i})$ are the method-specific time- and hazard-scale multipliers (PH, AFT, AH, GH share the same closed form via the unified $H(t \mid x) = H_0(t\, c_1)\, c_2$ representation).
Output shape:
- no time argument → length-
nvector with each subject evaluated at their own observed time $T_i$; t::Real→ length-nvector at the scalar time;ts::AbstractVector→n × length(ts)matrix.
With newdata::DataFrame the design matrices are rebuilt by applying the fit's stored formula(s) — newdata must contain every predictor column referenced in the original @formula(...). Newdata predict requires an explicit time argument (no "own time" default).
SurvivalModels.predict_survival — Method
predict_survival(m::GeneralHazardModel)
predict_survival(m::GeneralHazardModel, t::Real)
predict_survival(m::GeneralHazardModel, ts::AbstractVector)
predict_survival(m::GeneralHazardModel, newdata::DataFrame, t::Real)
predict_survival(m::GeneralHazardModel, newdata::DataFrame, ts::AbstractVector)Per-subject survival probability $S_i(t) = \exp(-\Lambda_i(t))$ derived from predict_expected. Shapes match predict_expected; newdata variants require an explicit time argument.
Simulation
The simulate command (ported from HazReg.jl; formerly simGH, now a deprecated alias) allows one to simulate times to event from the following models:
- General Hazard (GH) model [5] [6].
- Accelerated Failure Time (AFT) model [7].
- Proportional Hazards (PH) model [8].
- Accelerated Hazards (AH) model [9].
SurvivalModels.simulate — Function
simulate(n, method, baseline, X1, X2, α, β; rng=Random.default_rng())
simulate(n, method, baseline, formula[, formula2], df, α, β; rng=Random.default_rng())
simulate(n, model::GeneralHazardModel; rng=Random.default_rng())Simulate n times to event from a general-hazard model with hazard structure method (PHMethod()/AFTMethod()/AHMethod()/GHMethod()) and baseline distribution. Each of the n rows of the design (X1/X2, or built from formula/df via modelcols exactly as in fit) yields one event time by inverting that subject's survival Sᵢ(t)=U, i.e. Tᵢ = quantile(baseline, 1 − (1−Uᵢ)^(1/c2ᵢ)) / c1ᵢ. α are the X2 coefficients and β the X1 coefficients (only the ones the structure uses).
The fitted-model form delegates here using the model's own design/coefficients.
References: