diff --git a/tests/testthat/helper-integration-tests-setup.R b/tests/testthat/helper-integration-tests-setup.R
index cd0f3210..4281d767 100644
--- a/tests/testthat/helper-integration-tests-setup.R
+++ b/tests/testthat/helper-integration-tests-setup.R
@@ -329,12 +329,6 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
   # of how that is done. Other sections of the code below leave defaults in
   # place as appropriate.
 
-  # set up logR_sd
-  # logR_sd is NOT logged. It needs to enter the model logged b/c the exp() is
-  # taken before the likelihood calculation
-  recruitment$log_sigma_recruit$value <- log(om_input$logR_sd)
-  recruitment$log_sigma_recruit$is_random_effect <- FALSE
-  recruitment$log_sigma_recruit$estimated <- FALSE
   # set up log_rzero (equilibrium recruitment)
   recruitment$log_rzero$value <- log(om_input$R0)
   recruitment$log_rzero$is_random_effect <- FALSE
@@ -344,11 +338,26 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
   recruitment$logit_steep$is_random_effect <- FALSE
   recruitment$logit_steep$estimated <- FALSE
   # turn on estimation of deviations
-  recruitment$estimate_log_devs <- TRUE
   # recruit deviations should enter the model in normal space.
   # The log is taken in the likelihood calculations
   # alternative setting: recruitment$log_devs <- rep(0, length(om_input$logR.resid))
-  recruitment$log_devs <- om_input$logR.resid[-1]
+  recruitment$log_devs <- methods::new(ParameterVector, om_input$logR.resid[-1], om_input$nyr-1)
+
+  recruitment_distribution <- new(TMBDnormDistribution)
+  # set up logR_sd using the normal log_sd parameter
+  # logR_sd is NOT logged. It needs to enter the model logged b/c the exp() is
+  # taken before the likelihood calculation
+  recruitment_distribution$log_sd <- new(ParameterVector, 1)
+  recruitment_distribution$log_sd[1]$value <- log(om_input$logR_sd)
+  recruitment_distribution$log_sd[1]$estimated <- FALSE
+  recruitment_distribution$x <- new(ParameterVector, om_input$nyr)
+  recruitment_distribution$expected_values <- new(ParameterVector, om_input$nyr)
+  for (i in 1:om_input$nyr) {
+    recruitment_distribution$x[i]$value <- 0
+    recruitment_distribution$expected_values[i]$value <- 0
+  }
+  recruitment_distribution$set_distribution_links("random_effects", recruitment$log_devs$get_id())
+  recruitment$estimate_log_devs <- TRUE
 
   # Data
   catch <- em_input$L.obs$fleet1
@@ -396,23 +405,29 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
   fishing_fleet <- new(Fleet)
   fishing_fleet$nages <- om_input$nages
   fishing_fleet$nyears <- om_input$nyr
-  fishing_fleet$log_Fmort <- log(om_output$f)
-  fishing_fleet$estimate_F <- TRUE
-  fishing_fleet$random_F <- FALSE
+  fishing_fleet$log_Fmort <- methods::new(ParameterVector, log(om_output$f), om_input$nyr)
+  fishing_fleet$log_Fmort$set_all_estimable(TRUE)
   fishing_fleet$log_q <- log(1.0)
   fishing_fleet$estimate_q <- FALSE
   fishing_fleet$random_q <- FALSE
-  fishing_fleet$log_obs_error <- rep(log(sqrt(log(em_input$cv.L$fleet1^2 + 1))), om_input$nyr)
-  fishing_fleet$estimate_obs_error <- FALSE
-  # Modules are linked together using module IDs
-  # Each module has a get_id() function that returns the unique ID for that module
-  # Each fleet uses the module IDs to link up the correct module to the correct fleet
-  # Note: Likelihoods not yet set up as a stand-alone modules, so no get_id()
-  fishing_fleet$SetAgeCompLikelihood(1)
-  fishing_fleet$SetIndexLikelihood(1)
   fishing_fleet$SetSelectivity(fishing_fleet_selectivity$get_id())
-  fishing_fleet$SetObservedIndexData(fishing_fleet_index$get_id())
-  fishing_fleet$SetObservedAgeCompData(fishing_fleet_age_comp$get_id())
+
+  # Set up fishery index data using the lognormal
+  fishing_fleet_index_distribution <- methods::new(TMBDlnormDistribution)
+  # lognormal observation error transformed on the log scale
+  fishing_fleet_index_distribution$log_logsd <- new(ParameterVector, om_input$nyr)
+  for (y in 1:om_input$nyr) {
+    fishing_fleet_index_distribution$log_logsd[y]$value <- log(sqrt(log(em_input$cv.L$fleet1^2 + 1)))
+  }
+  fishing_fleet_index_distribution$log_logsd$set_all_estimable(FALSE)
+  # Set Data using the IDs from the modules defined above
+  fishing_fleet_index_distribution$set_observed_data(fishing_fleet_index$get_id())
+  fishing_fleet_index_distribution$set_distribution_links("data", fishing_fleet$log_expected_index$get_id())
+
+  # Set up fishery age composition data using the multinomial
+  fishing_fleet_agecomp_distribution <- methods::new(TMBDmultinomDistribution)
+  fishing_fleet_agecomp_distribution$set_observed_data(fishing_fleet_age_comp$get_id())
+  fishing_fleet_agecomp_distribution$set_distribution_links("data", fishing_fleet$proportion_catch_numbers_at_age$get_id())
 
   # Create the survey fleet
   survey_fleet_selectivity <- new(LogisticSelectivity)
@@ -429,25 +444,36 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
   survey_fleet$is_survey <- TRUE
   survey_fleet$nages <- om_input$nages
   survey_fleet$nyears <- om_input$nyr
-  survey_fleet$estimate_F <- FALSE
-  survey_fleet$random_F <- FALSE
   survey_fleet$log_q <- log(om_output$survey_q$survey1)
   survey_fleet$estimate_q <- TRUE
   survey_fleet$random_q <- FALSE
-  survey_fleet$log_obs_error <- rep(log(sqrt(log(em_input$cv.survey$survey1^2 + 1))), om_input$nyr)
-  survey_fleet$estimate_obs_error <- FALSE
-  survey_fleet$SetAgeCompLikelihood(1)
-  survey_fleet$SetIndexLikelihood(1)
   survey_fleet$SetSelectivity(survey_fleet_selectivity$get_id())
-  survey_fleet$SetObservedIndexData(survey_fleet_index$get_id())
-  survey_fleet$SetObservedAgeCompData(survey_fleet_age_comp$get_id())
+
+  # Set up survey index data using the lognormal
+  survey_fleet_index_distribution <- methods::new(TMBDlnormDistribution)
+  # lognormal observation error transformed on the log scale
+  # sd = sqrt(log(cv^2 + 1)), sd is log transformed
+  survey_fleet_index_distribution$log_logsd <- new(ParameterVector, om_input$nyr)
+  for (y in 1:om_input$nyr) {
+    survey_fleet_index_distribution$log_logsd[y]$value <- log(sqrt(log(em_input$cv.survey$survey1^2 + 1)))
+  }
+  survey_fleet_index_distribution$log_logsd$set_all_estimable(FALSE)
+  # Set Data using the IDs from the modules defined above
+  survey_fleet_index_distribution$set_observed_data(survey_fleet_index$get_id())
+  survey_fleet_index_distribution$set_distribution_links("data", survey_fleet$log_expected_index$get_id())
+
+  # Age composition data
+
+  survey_fleet_agecomp_distribution <- methods::new(TMBDmultinomDistribution)
+  survey_fleet_agecomp_distribution$set_observed_data(survey_fleet_age_comp$get_id())
+  survey_fleet_agecomp_distribution$set_distribution_links("data", survey_fleet$proportion_catch_numbers_at_age$get_id())
 
   # Population
   population <- new(Population)
-  population$log_M <- rep(log(om_input$M.age[1]), om_input$nyr * om_input$nages)
-  population$estimate_M <- FALSE
-  population$log_init_naa <- log(om_output$N.age[1, ])
-  population$estimate_init_naa <- TRUE
+  population$log_M <- methods::new(ParameterVector, rep(log(om_input$M.age[1]), om_input$nyr * om_input$nages), om_input$nyr * om_input$nages)
+  population$log_M$set_all_estimable(FALSE)
+  population$log_init_naa <- methods::new(ParameterVector, log(om_output$N.age[1, ]), om_input$nages)
+  population$log_init_naa$set_all_estimable(TRUE)
   population$nages <- om_input$nages
   population$ages <- om_input$ages
   population$nfleets <- sum(om_input$fleet_num, om_input$survey_num)
@@ -461,6 +487,7 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
   CreateTMBModel()
   # Create parameter list from Rcpp modules
   parameters <- list(p = get_fixed())
+
   input <- list()
   input$parameters <- parameters
   input$version <- "Model Comparison Project example"
@@ -468,14 +495,7 @@ setup_and_run_FIMS_with_wrappers <- function(iter_id,
 
   clear()
   # Return the results as a list
-  return(list(
-    parameters = fit$input$parameters,
-    obj = fit$obj,
-    opt = fit$opt,
-    report = fit$rep,
-    sdr_report = fit$sd,
-    fit = fit
-  ))
+  return(fit)
 
 
 }
diff --git a/tests/testthat/test-integration-fims-estimation-with-wrappers.R b/tests/testthat/test-integration-fims-estimation-with-wrappers.R
index da9c6f6c..d44abc18 100644
--- a/tests/testthat/test-integration-fims-estimation-with-wrappers.R
+++ b/tests/testthat/test-integration-fims-estimation-with-wrappers.R
@@ -11,20 +11,12 @@ test_that("deterministic test of fims", {
     estimation_mode = FALSE
   )
 
-  # Set up TMB's computational graph
-  obj <- result$obj
-
-  # Calculate standard errors
-  sdr <- TMB::sdreport(obj)
-  sdr_report <- summary(sdr, "report")
-  sdr_fixed <- summary(sdr, "fixed")
-
   # Call report using deterministic parameter values
   # obj$report() requires parameter list to avoid errors
-  report <- obj$report(obj$par)
+  report <- result$rep
 
   # Compare log(R0) to true value
-  fims_logR0 <- sdr_fixed[1, "Estimate"]
+  fims_logR0 <- as.numeric(result$obj$par[1])
   expect_gt(fims_logR0, 0.0)
   expect_equal(fims_logR0, log(om_input_list[[iter_id]]$R0))
 
@@ -154,11 +146,12 @@ test_that("estimation test of fims using fit_fims()", {
     estimation_mode = TRUE
   )
 
+  result$sd
   # Compare FIMS results with model comparison project OM values
   validate_fims(
-    report = result$report,
-    sdr = TMB::sdreport(result$obj),
-    sdr_report = result$sdr_report,
+    report = result$rep,
+    sdr = result$sd,
+    sdr_report = result$sd,
     om_input = om_input_list[[iter_id]],
     om_output = om_output_list[[iter_id]],
     em_input = em_input_list[[iter_id]],
diff --git a/tests/testthat/test-integration-fims-estimation-without-wrappers.R b/tests/testthat/test-integration-fims-estimation-without-wrappers.R
index 733ac926..569907eb 100644
--- a/tests/testthat/test-integration-fims-estimation-without-wrappers.R
+++ b/tests/testthat/test-integration-fims-estimation-without-wrappers.R
@@ -144,7 +144,7 @@ test_that("deterministic test of fims", {
 test_that("nll test of fims", {
   iter_id <- 1
 
-  result <- setup_and_run_FIMS(
+  result <- setup_and_run_FIMS_without_wrappers(
     iter_id = iter_id,
     om_input_list = om_input_list,
     om_output_list = om_output_list,
diff --git a/tests/testthat/test-parallel-with-snowfall-with-wrappers.R b/tests/testthat/test-parallel-with-snowfall-with-wrappers.R
index 100eae9c..70a5bb5f 100644
--- a/tests/testthat/test-parallel-with-snowfall-with-wrappers.R
+++ b/tests/testthat/test-parallel-with-snowfall-with-wrappers.R
@@ -55,14 +55,14 @@ test_that("Run FIMS in parallel using {snowfall}", {
   # Compare parameters in results:
   # Verify that the results from both runs are equivalent.
   expect_setequal(
-    unname(unlist(lapply(results_parallel, `[[`, "parameters"))),
-    unname(unlist(lapply(estimation_results_serial, `[[`, "parameters")))
+    unname(unlist(lapply(results_parallel, `[[`, "parList"))),
+    unname(unlist(lapply(estimation_results_serial, `[[`, "parList")))
   )
 
   # Compare sdr_report values in results:
   # Verify that the results from both runs are equivalent.
   expect_setequal(
-    unlist(lapply(results_parallel, `[[`, "sdr_report")),
-    unlist(lapply(estimation_results_serial, `[[`, "sdr_report"))
+    unname(unlist(lapply(results_parallel, `[[`, "sd"))),
+    unname(unlist(lapply(estimation_results_serial, `[[`, "sd")))
   )
 })