diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml
index 8b21a5c..2de6cbc 100644
--- a/.github/workflows/R-CMD-check.yaml
+++ b/.github/workflows/R-CMD-check.yaml
@@ -129,6 +129,15 @@ jobs:
             dependencies = FALSE,
             upgrade = "never"
           )
+          # Ensure testthat is exactly 3.1.8
+          required_testthat <- package_version("3.1.8")
+          testthat_ok <- requireNamespace("testthat", quietly = TRUE) &&
+            packageVersion("testthat") == required_testthat
+
+          if (!testthat_ok) {
+            install.packages('ellipsis')
+            install.packages('https://cran.r-project.org/src/contrib/Archive/testthat/testthat_3.1.8.tar.gz', repos = NULL, type = 'source')
+          }
         shell: Rscript {0}
 
       - name: Check
diff --git a/tests/testthat/test-estim_param.R b/tests/testthat/test-estim_param.R
index 0e9f746..42296f8 100644
--- a/tests/testthat/test-estim_param.R
+++ b/tests/testthat/test-estim_param.R
@@ -145,18 +145,10 @@ nlo <- lapply(nlo, function(x) {
 })
 
 test_that("Test Vignette simple_case", {
-  expect_equal(nlo_new[[1]]$x0, c(0.0023862966, 118.3769), tolerance = 1e-7)
-  expect_equal(nlo_new[[2]]$x0, c(0.0008777006, 290.9086), tolerance = 1e-7)
-  expect_equal(
-    sapply(nlo_new, "[[", "solution"),
-    sapply(nlo, "[[", "solution"),
-    tolerance = 1e-4
-  )
-  expect_equal(
-    sapply(nlo_new, "[[", "objective"),
-    sapply(nlo, "[[", "objective"),
-    tolerance = 1e-4
-  )
+  expect_equal(as.numeric(res$init_values[1,]), c(0.0023862966, 118.3769), tolerance = 1e-7)
+  expect_equal(as.numeric(res$init_values[2,]), c(0.0008777006, 290.9086), tolerance = 1e-7)
+  expect_equivalent(t(res$est_values), sapply(nlo, "[[", "solution"), tolerance = 1e-4)
+  expect_equal(res$crit_values, sapply(nlo, "[[", "objective"), tolerance = 1e-4)
   expect_true(file.exists(file.path(data_dir, "EstimatedVSinit.pdf")))
   expect_true(file.exists(file.path(data_dir, "ValuesVSit.pdf")))
   expect_true(file.exists(file.path(data_dir, "ValuesVSit_2D.pdf")))
@@ -274,41 +266,16 @@ source(file.path(tmpdir, "Parameter_estimation_Specific_and_Varietal.R"))
 
 ## load the results
 load(file.path(data_dir, "optim_results.Rdata"))
-nlo_new <- lapply(res$nlo, function(x) {
-  x$call <- NULL
-  x
-}) # remove "call" since it may change between code versions ...
-load(system.file(
-  file.path("extdata", "ResSpecVar_2rep4it", stics_version),
-  "optim_results.Rdata",
-  package = "CroptimizR"
-))
-nlo <- lapply(res$nlo, function(x) {
-  x$call <- NULL
-  x
-}) # remove "call" since it may change between code versions ...
+res_new <- res
+rm(res)
+load(system.file(file.path("extdata", "ResSpecVar_2rep4it", stics_version), "optim_results.Rdata", package = "CroptimizR"))
+res_ref <- res
 
 test_that("Test Vignette specific and varietal", {
-  expect_equal(
-    nlo_new[[1]]$x0,
-    c(0.001386297, 293.3769, 299.3398),
-    tolerance = 1e-7
-  )
-  expect_equal(
-    nlo_new[[2]]$x0,
-    c(0.001877701, 115.9086, 162.9456),
-    tolerance = 1e-7
-  )
-  expect_equal(
-    sapply(nlo_new, "[[", "solution"),
-    sapply(nlo, "[[", "solution"),
-    tolerance = 1e-4
-  )
-  expect_equal(
-    sapply(nlo_new, "[[", "objective"),
-    sapply(nlo, "[[", "objective"),
-    tolerance = 1e-4
-  )
+  expect_equal(as.numeric(res_new$init_values[1,]), c(0.001386297, 293.3769, 299.3398), tolerance = 1e-7)
+  expect_equal(as.numeric(res_new$init_values[2,]), c(0.001877701, 115.9086, 162.9456), tolerance = 1e-7)
+  expect_equal(res_new$est_values, res$est_values, tolerance = 1e-4)
+  expect_equal(res_new$crit_values, res$crit_values, tolerance = 1e-4)
   expect_true(file.exists(file.path(data_dir, "EstimatedVSinit.pdf")))
   expect_true(file.exists(file.path(data_dir, "ValuesVSit.pdf")))
   expect_true(file.exists(file.path(data_dir, "ValuesVSit_2D.pdf")))
@@ -796,13 +763,14 @@ test_that("Test rotation", {
 #   expect_true(file.exists(file.path(out_dir,"marginalPlots.pdf")))
 # })
 
-# Test Vignette AgMIP phase III protocol
-# --------------------------------------
+
+# Test Vignette Parameter selection
+# ---------------------------------
 
 tmpdir <- normalizePath(tempdir(), winslash = "/", mustWork = FALSE)
 vignette_rmd <- file.path(tmpdir, "AgMIP_Calibration_protocol.Rmd")
 download.file(
-  "https://raw.github.com/SticsRPacks/CroptimizR/main/vignettes/AgMIP_Calibration_protocol.Rmd",
+  "https://raw.github.com/SticsRPacks/CroptimizR/main/vignettes/Parameter_selection.Rmd",
   vignette_rmd
 )
 
@@ -874,8 +842,7 @@ if (Sys.getenv("CI") != "") {
 ## generate the R script
 knitr::purl(
   input = vignette_rmd,
-  output = file.path(tmpdir, "AgMIP_Calibration_protocol.R"),
-  documentation = 2
+  output = file.path(tmpdir, "Parameter_selection.R"), documentation = 2
 )
 
 ## Seems that optim_options and optim_results.Rdata are not overwritten => try to remove them before run
@@ -885,23 +852,14 @@ if (file.exists(file.path(data_dir, "optim_results.Rdata"))) {
 rm(optim_options, param_info)
 
 ## run it
-source(file.path(tmpdir, "AgMIP_Calibration_protocol.R"))
+source(file.path(tmpdir, "Parameter_selection.R"))
 
 ## load the results
 load(file.path(data_dir, "optim_results.Rdata"))
-nlo_new <- lapply(res$nlo, function(x) {
-  x$call <- NULL
-  x
-}) # remove "call" since it may change between code versions ...
-load(system.file(
-  file.path("extdata", "ResAgmipPheno_2rep4it", stics_version),
-  "optim_results.Rdata",
-  package = "CroptimizR"
-))
-nlo <- lapply(res$nlo, function(x) {
-  x$call <- NULL
-  x
-}) # remove "call" since it may change between code versions ...
+res_new <- res
+rm(res)
+load(system.file(file.path("extdata", "ResAgmipPheno_2rep4it", stics_version), "optim_results.Rdata", package = "CroptimizR"))
+res_ref <- res
 
 test_that("Test Vignette AgMIP Phase III protocol", {
   expect_equal(optim_options$nb_rep, 2)
@@ -916,18 +874,10 @@ test_that("Test Vignette AgMIP Phase III protocol", {
       tdmax = c(30.67603, 27.26267)
     )
   )
-  expect_equal(nlo_new[[1]]$x0, c(138.2656, 506.3340), tolerance = 1e-4)
-  expect_equal(nlo_new[[2]]$x0, c(446.0060, 639.4217), tolerance = 1e-4)
-  expect_equal(
-    sapply(nlo_new, "[[", "solution"),
-    sapply(nlo, "[[", "solution"),
-    tolerance = 1e-4
-  )
-  expect_equal(
-    sapply(nlo_new, "[[", "objective"),
-    sapply(nlo, "[[", "objective"),
-    tolerance = 1e-4
-  )
+  expect_equal(as.numeric(res_new$init_values[1,]), c(138.2656, 506.3340), tolerance = 1e-4)
+  expect_equal(as.numeric(res_new$init_values[2,]), c(446.0060, 639.4217), tolerance = 1e-4)
+  expect_equal(res_new$est_values, res$est_values, tolerance = 1e-4)
+  expect_equal(res_new$crit_values, res$crit_values, tolerance = 1e-4)
   expect_true(file.exists(file.path(data_dir, "param_selection_steps.csv")))
   for (i in 1:4) {
     expect_true(file.exists(file.path(