The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.
Most users arrive with a spreadsheet rather than an R data frame. This vignette covers the file-driven path: from a CSV on disk to fitted models, tables and figures, without writing any modelling code.
Two columns are enough:
| column | meaning |
|---|---|
time |
survival or failure time, strictly positive |
status |
1 if the event was observed, 0 if right-censored |
Covariates go in additional columns. For competing risks, replace
status with cause, coded 0 for censored and 1,
2, … for the causes.
If you would rather start from a working example than a description, write a template and fill it in:
read_survival_data() will guess the time and status
columns when they carry recognisable names, and reports what it
chose:
f <- system.file("extdata", "censored_sample.csv", package = "BetaDanish")
dat <- read_survival_data(f)
head(dat)
#> time status
#> 1 11 1
#> 2 14 1
#> 3 23 1
#> 4 31 1
#> 5 32 1
#> 6 35 1Note that covariates are not retained by default. Ask for them by
name, or use covar_cols = "all":
dat_all <- read_survival_data(f, covar_cols = "all", quiet = TRUE)
names(dat_all)
#> [1] "time" "status" "group"The attached report records what happened, which is worth checking before you model anything:
str(attr(dat_all, "bd_data_report"))
#> List of 13
#> $ file : chr "C:\\Users\\hp\\AppData\\Local\\Temp\\RtmpSsnbRT\\Rinst50843b8c2feb\\BetaDanish\\extdata\\censored_sample.csv"
#> $ time_col : chr "time"
#> $ status_col : chr "status"
#> $ cause_col : NULL
#> $ available_columns: chr [1:3] "time" "status" "group"
#> $ covariates : chr "group"
#> $ dropped_columns : chr(0)
#> $ rows_read : int 91
#> $ rows_kept : int 91
#> $ rows_dropped : int 0
#> $ n_events : int 69
#> $ censoring_prop : num 0.242
#> $ grid_step : num 1One field deserves attention. grid_step is
non-NA when the times look recorded on a coarse grid –
whole days or whole months, say. Rounded times break the assumptions
behind a point-density likelihood, so treat standard errors from such
data as optimistic.
censorColumn-name guessing deliberately ignores names like
censor. In some conventions censor = 1 means
censored; in others it means observed. Guessing wrong
would silently invert every event in your dataset and still produce a
plausible-looking fit, so you must name that column yourself.
bd_analyze_csv() does the whole thing: read, fit,
tabulate, and optionally write results to a directory.
res <- bd_analyze_csv(f, analysis = "univariate", model = "ED",
compare = FALSE, n_starts = 5, seed = 1, quiet = TRUE)
res
#>
#> BetaDanish CSV analysis
#> -----------------------
#> Analysis: univariate
#> Observations: 91 (events: 69, censored: 24.2%)
#> Models fitted: ED
#>
#> model n npar logLik AIC BIC
#> ED 91 3 -476.8013 959.6026 967.1352
#>
#> Warnings: 1 (see $warnings)The tidy tables are in $tables:
res$tables$estimates
#> model parameter estimate std_error lower_95 upper_95
#> ED.1 ED b 0.5297814 0.08547109 0.3622611 0.6973016
#> ED.2 ED c 17.4808483 22.38401058 -26.3910063 61.3527029
#> ED.3 ED k 0.2597284 0.35870994 -0.4433301 0.9627870
res$tables$goodness_of_fit
#> model AIC BIC HQIC AICC KS
#> ED ED 959.6026 967.1352 962.6415 959.8784 0.0814148With model = "both", both the four-parameter Beta-Danish
and the three-parameter Exponentiated Danish submodel are fitted and
compared:
both <- bd_analyze_csv(f, model = "both", compare = FALSE,
n_starts = 5, seed = 2, quiet = TRUE)
#> Likelihood Ratio Test (a = 1 vs a != 1)
#>
#> Model LogLik Chisq Df Pr(>Chisq)
#> 1 Submodel (3-param) -476.8013 NA NA NA
#> 2 Full Model (4-param) -476.8225 0 1 1
both$tables$information_criteria
#> model n npar logLik AIC BIC
#> 1 BD 91 4 -476.8225 961.6450 971.6884
#> 2 ED 91 3 -476.8013 959.6026 967.1352
both$tables$likelihood_ratio_test
#> Model LogLik Chisq Df Pr(>Chisq)
#> 1 Submodel (3-param) -476.8013 NA NA NA
#> 2 Full Model (4-param) -476.8225 0 1 1Read that test alongside the identifiability diagnostics. On many
datasets the four-parameter model fits well but is only weakly
identified, and the likelihood ratio test will not reject the submodel.
See the Identifiability section of ?fit_betadanish, and the
guinea_pig dataset for a contrasting case where the parent
model is well identified.
Nothing is written to disk unless you name a directory. When you do, you get one CSV per table and a PNG per diagnostic figure:
out <- file.path(tempdir(), "bd_results")
saved <- bd_analyze_csv(f, model = "ED", compare = FALSE,
output_dir = out, n_starts = 5, seed = 3, quiet = TRUE)
basename(saved$files)
#> [1] "estimates.csv" "information_criteria.csv"
#> [3] "goodness_of_fit.csv" "data_report.csv"
#> [5] "survival.png" "hazard.png"
#> [7] "density.png" "pp.png"
#> [9] "qq.png"The same entry point handles the regression models. Covariates are taken from the file, either all of them or the ones you name:
g <- system.file("extdata", "covariate_sample.csv", package = "BetaDanish")
# Accelerated failure time
aft <- bd_analyze_csv(g, analysis = "aft",
covariates = c("age", "thickness"))
# Mixture cure model, cure fraction depending on ulceration
cure <- bd_analyze_csv(g, analysis = "cure", cure_formula = ~ ulcer)
# Competing risks
h <- system.file("extdata", "competing_sample.csv", package = "BetaDanish")
cr <- bd_analyze_csv(h, analysis = "competing", cause_col = "cause")Each fit is attempted independently. If one fails, the message is recorded and the rest of the analysis still runs, so a difficult four-parameter fit does not cost you the submodel results beside it:
An empty result means everything succeeded.
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.