exe_path_dev <- "../inst/bin/hisea.exe"
exe_path_inst <- system.file("bin/hisea.exe", package = "RHISEA")
exe_source <- if (file.exists(exe_path_dev)) exe_path_dev else exe_path_inst
if (file.exists(exe_source)) {
file.copy(exe_source, "hisea.exe", overwrite = TRUE)
} else {
stop("Not available")
}
#> [1] TRUE
on.exit({
unlink(list.files(work_dir, pattern = "\\.(std|mix|out|xlsx|rda|exe)$|^fort\\."), force = TRUE)
setwd(old_dir)
}, add = TRUE)
for (i in 1:n_sim) {
cat("\n--- Simulation ", i, " ---\n")
# Activate independent RNG stream for reproducibility
.Random.seed <- streams[[i]]
# Simulate baseline stock data with slightly different lambda values per stock
base_lambda_stock1 <- runif(1, 15, 25)
lambda_shift_stock2 <- runif(1, 5, 15)
stock1_data <- generate_stock_data(100, nv, base_lambda_stock1)
stock2_data <- generate_stock_data(100, nv, base_lambda_stock1 + lambda_shift_stock2)
# Write baseline to .std file
write_std(baseline_std_file, list(stock1_data, stock2_data))
# Generate mixture based on random proportions
current_sim_comp <- as.numeric(gtools::rdirichlet(1, alpha = rep(1, np)))
n_from_stocks <- round(Nmix * current_sim_comp)
if (sum(n_from_stocks) != Nmix) n_from_stocks[1] <- n_from_stocks[1] + (Nmix - sum(n_from_stocks))
mixture_data <- rbind(
stock1_data[sample(nrow(stock1_data), n_from_stocks[1], replace = TRUE), ],
stock2_data[sample(nrow(stock2_data), n_from_stocks[2], replace = TRUE), ]
)
mixture_data <- mixture_data[sample(nrow(mixture_data)), ]
# Write mixture to .mix file
write_mix(mixture_mix_file, mixture_data)
# Run Fortran HISEA executable
system2("./hisea.exe", wait = TRUE, stdout = FALSE)
out_file <- sprintf("hisea%02d.out", i)
if (file.exists("hisea.out")) file.rename("hisea.out", out_file)
# Run R-based estimations with LDA and LDA_MASS
r_result <- run_hisea_all(
type = "ANALYSIS",
np = np,
phi_method = "standard",
nv = nv,
resample_baseline = resample_baseline,
resampled_baseline_sizes = resampled_baseline_sizes,
seed_val = 123456,
nsamps = Nsamps,
Nmix = Nmix,
baseline_input = baseline_std_file,
mix_input = mixture_mix_file,
method_class = "LDA",
stocks_names = stock_labels
)
r_df <- as.data.frame(r_result$mean)
colnames(r_df) <- paste0("Estimator_", 1:5)
r_df$Stock <- stock_labels
r_df$Simulation <- paste0("Sim_", i)
results[[i]] <- r_df
r_result2 <- run_hisea_all(
type = "ANALYSIS",
np = np,
phi_method = "standard",
nv = nv,
resample_baseline = resample_baseline,
resampled_baseline_sizes = resampled_baseline_sizes,
seed_val = 123456,
nsamps = Nsamps,
Nmix = Nmix,
baseline_input = baseline_std_file,
mix_input = mixture_mix_file,
method_class = "LDA_MASS",
stocks_names = stock_labels
)
r_df2 <- as.data.frame(r_result2$mean)
colnames(r_df2) <- paste0("Estimator_", 1:5)
r_df2$Stock <- stock_labels
r_df2$Simulation <- paste0("Sim_", i)
results2[[i]] <- r_df2
# Read Fortran estimators from fort.xx files (10 + estimator index)
fort_estimates <- lapply(1:5, function(j) read_fort_estimator(sprintf("fort.%02d", 10 + j)))
# Store wide and long format comparisons for each estimator
for (j in 1:5) {
wide_df <- data.frame(
Simulation = paste0("Sim_", i),
Stock = stock_labels,
R = r_df[[paste0("Estimator_", j)]],
Mass = r_df2[[paste0("Estimator_", j)]],
Fortran = fort_estimates[[j]]
)
est_wide_comparisons[[j]] <- rbind(est_wide_comparisons[[j]], wide_df)
long_df <- data.frame(
Simulation = rep(paste0("Sim_", i), times = 6),
Estimate = paste0("Estimator_", j),
Stock = rep(stock_labels, times = 3),
Method = rep(c("Fortran", "LDA_R", "LDA_MASS"), each = 2),
Value = c(
fort_estimates[[j]],
r_df[[paste0("Estimator_", j)]],
r_df2[[paste0("Estimator_", j)]]
)
)
est_long_comparisons[[j]] <- rbind(est_long_comparisons[[j]], long_df)
}
}
#>
#> --- Simulation 1 ---
#>
#> --- Simulation 2 ---
#>
#> --- Simulation 3 ---
#>
#> --- Simulation 4 ---
#>
#> --- Simulation 5 ---
#>
#> --- Simulation 6 ---
#>
#> --- Simulation 7 ---
#>
#> --- Simulation 8 ---
#>
#> --- Simulation 9 ---
#>
#> --- Simulation 10 ---
#>
#> --- Simulation 11 ---
#>
#> --- Simulation 12 ---
#>
#> --- Simulation 13 ---
#>
#> --- Simulation 14 ---
#>
#> --- Simulation 15 ---
#>
#> --- Simulation 16 ---
#>
#> --- Simulation 17 ---
#>
#> --- Simulation 18 ---
#>
#> --- Simulation 19 ---
#>
#> --- Simulation 20 ---
#>
#> --- Simulation 21 ---
#>
#> --- Simulation 22 ---
#>
#> --- Simulation 23 ---
#>
#> --- Simulation 24 ---
#>
#> --- Simulation 25 ---
#>
#> --- Simulation 26 ---
#>
#> --- Simulation 27 ---
#>
#> --- Simulation 28 ---
#>
#> --- Simulation 29 ---
#>
#> --- Simulation 30 ---
#>
#> --- Simulation 31 ---
#>
#> --- Simulation 32 ---
#>
#> --- Simulation 33 ---
#>
#> --- Simulation 34 ---
#>
#> --- Simulation 35 ---
#>
#> --- Simulation 36 ---
#>
#> --- Simulation 37 ---
#>
#> --- Simulation 38 ---
#>
#> --- Simulation 39 ---
#>
#> --- Simulation 40 ---
#>
#> --- Simulation 41 ---
#>
#> --- Simulation 42 ---
#>
#> --- Simulation 43 ---
#>
#> --- Simulation 44 ---
#>
#> --- Simulation 45 ---
#>
#> --- Simulation 46 ---
#>
#> --- Simulation 47 ---
#>
#> --- Simulation 48 ---
#>
#> --- Simulation 49 ---
#>
#> --- Simulation 50 ---