inst/examples/swissmetro/plot_b01c_logit.R

In rbiogeme: Interface to 'Biogeme'

#!/usr/bin/env Rscript

# b01c. Quick estimation of a multinomial logit model
#
# This example uses the same model specification as b01a, but delegates to
# native BIOGEME.quick_estimate(). Quick estimation is useful when estimated
# parameter values and lightweight statistics are needed without the usual
# derivative-based post-estimation calculations.

library(rbiogeme)

# The shared helper contains command-line parsing and data preparation. The
# complete model specification remains in this script.
script_path <- commandArgs(trailingOnly = FALSE)
script_path <- sub("^--file=", "", script_path[startsWith(script_path, "--file=")][[1L]])
source(file.path(dirname(normalizePath(script_path)), "example_utils.R"))

build_b01c_model <- function(database) {
  # Each biogeme_beta() call creates a symbolic native Biogeme parameter.
  # The Swissmetro ASC is fixed at zero for identification.
  asc_car <- biogeme_beta("asc_car", start = 0)
  asc_train <- biogeme_beta("asc_train", start = 0)
  asc_sm <- biogeme_beta("asc_sm", start = 0, fixed = TRUE)
  b_time <- biogeme_beta("b_time", start = 0)
  b_cost <- biogeme_beta("b_cost", start = 0)

  # variable() creates symbolic data nodes. R arithmetic combines these
  # nodes into the utility expressions; it does not evaluate the data here.
  logit_model(
    database = database,
    choice = "CHOICE",
    utilities = list(
      `1` = asc_train + b_time * variable("TRAIN_TT_SCALED") +
        b_cost * variable("TRAIN_COST_SCALED"),
      `2` = asc_sm + b_time * variable("SM_TT_SCALED") +
        b_cost * variable("SM_COST_SCALED"),
      `3` = asc_car + b_time * variable("CAR_TT_SCALED") +
        b_cost * variable("CAR_CO_SCALED")
    ),
    availability = list(
      `1` = variable("TRAIN_AV_SP"),
      `2` = variable("SM_AV"),
      `3` = variable("CAR_AV_SP")
    )
  )
}

# prepare_swissmetro_example() is defined in example_utils.R. It parses the
# command line, validates the data/Python paths, configures the bridge, reads
# the data, and creates a fresh output directory.
prepared <- prepare_swissmetro_example(
  commandArgs(trailingOnly = TRUE),
  default_model = "b01c_logit"
)

# quick_estimate() does not load a YAML result, but explicitly remove the
# exact b01c result/iteration files so a reused output directory cannot make
# this example appear to succeed from stale files.
stale_files <- c("b01c_logit.yaml", "__b01c_logit.iter")
stale_files <- file.path(prepared$output, stale_files)
stale_files <- stale_files[file.exists(stale_files)]
if (length(stale_files) > 0L) unlink(stale_files, force = TRUE)

database <- swissmetro_data(prepared$data)
model <- build_b01c_model(database)
control <- biogeme_control(
    output_directory = prepared$output,
  model_name = "b01c_logit",
  generate_html = FALSE,
  generate_yaml = FALSE,
  save_iterations = FALSE
)

# The complete symbolic model is compiled once, then native Biogeme performs
# quick estimation. No R callback runs inside the likelihood or optimizer.
fit <- quick_estimate(model, model_name = "b01c_logit", control = control)

# quick_estimate() intentionally does not write YAML automatically. Save the
# returned native result explicitly, matching the Python example's
# results.dump_yaml_file(filename=...).
print(summary(fit))
print(coef(fit))
print("General statistics")
print(biogeme_general_statistics(fit))
save_results(fit, filename = file.path(prepared$output, "b01c_logit.yaml"))

invisible(fit)