plot_b02nested.R
In rbiogeme: Interface to 'Biogeme'

#!/usr/bin/env Rscript

# Nested-logit estimation using native sampling of alternatives.
# The entire utility, partition, MEV sample, and nest specification is here.

library(rbiogeme)

# prepare_sampling_example() is defined in example_utils.R. It supplies only
# explicit input data and output paths, not model equations.
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"))

prepared <- prepare_sampling_example(
  commandArgs(trailingOnly = TRUE),
  default_model = "nested_downtown_20"
)
alternatives <- prepared$alternatives
observations <- prepared$observations

# Main sampling uses 20 alternatives from downtown and its complement. The
# MEV terms use all 33 Asian alternatives, exactly as in native b02nested.
all_alternatives <- sort(unique(as.integer(alternatives$ID)))
asian <- all_alternatives[alternatives$Asian[match(all_alternatives, alternatives$ID)] == 1]
downtown <- all_alternatives[alternatives$downtown[match(all_alternatives, alternatives$ID)] == 1]
partition <- biogeme_sampling_partition(
  segments = list(downtown, setdiff(all_alternatives, downtown)),
  sample_sizes = sampling_segment_sizes(20, 2),
  full_set = all_alternatives
)
mev_partition <- biogeme_sampling_partition(
  segments = list(asian),
  sample_sizes = 33,
  full_set = asian
)

log_dist <- cross_variable(
  "log_dist",
  log(((variable("user_lat") - variable("rest_lat"))^2 +
    (variable("user_lon") - variable("rest_lon"))^2)^0.5)
)

# Complete utility specification, with native parameter names.
beta_rating <- biogeme_beta("beta_rating", start = 0)
beta_price <- biogeme_beta("beta_price", start = 0)
beta_chinese <- biogeme_beta("beta_chinese", start = 0)
beta_japanese <- biogeme_beta("beta_japanese", start = 0)
beta_korean <- biogeme_beta("beta_korean", start = 0)
beta_indian <- biogeme_beta("beta_indian", start = 0)
beta_french <- biogeme_beta("beta_french", start = 0)
beta_mexican <- biogeme_beta("beta_mexican", start = 0)
beta_lebanese <- biogeme_beta("beta_lebanese", start = 0)
beta_ethiopian <- biogeme_beta("beta_ethiopian", start = 0)
beta_log_dist <- biogeme_beta("beta_log_dist", start = 0)
utility <- beta_rating * variable("rating") +
  beta_price * variable("price") +
  beta_chinese * variable("category_Chinese") +
  beta_japanese * variable("category_Japanese") +
  beta_korean * variable("category_Korean") +
  beta_indian * variable("category_Indian") +
  beta_french * variable("category_French") +
  beta_mexican * variable("category_Mexican") +
  beta_lebanese * variable("category_Lebanese") +
  beta_ethiopian * variable("category_Ethiopian") +
  beta_log_dist * variable("log_dist")

# The Asian nest contains all Asian alternatives; other alternatives are
# treated as native trivial nests by Biogeme.
mu_asian <- biogeme_beta("mu_asian", start = 1, lower = 1)
nests <- nested_nests(
  choice_set = all_alternatives,
  nests = list(nested_nest(mu_asian, asian, name = "asian"))
)

control <- biogeme_control(
    output_directory = prepared$output,
  generate_html = FALSE,
  generate_yaml = FALSE,
  save_iterations = FALSE
)
model <- sampled_alternatives_model(
  alternatives = alternatives,
  individuals = observations,
  choice_column = "nested_0",
  id_column = "ID",
  utility = utility,
  combined_variables = list(log_dist),
  partition = partition,
  mev_partition = mev_partition,
  mev_sample_sizes = 33,
  nests = nests,
  biogeme_file_name = file.path(prepared$output, "nested_downtown_20.dat"),
  model_type = "nested",
  control = control
)

fit <- estimate_sampled_alternatives(
  model,
  model_name = "nested_downtown_20",
  control = control
)
print(fit)
comparison <- compare_sampling_parameters(fit, sampling_true_parameters())
print(comparison$data)
print(comparison$message)

invisible(fit)