diff --git a/index.qmd b/index.qmd
index cfb5729..46e81f4 100644
--- a/index.qmd
+++ b/index.qmd
@@ -2919,17 +2919,6 @@ plot_forest_2 <- function(data, intercept = TRUE, MA_mean = TRUE, y_zoom = numer
   print(p)
 }
 
-
-fit_MA_mv <- function(effects_analysis, Z_colname, VZ_colname, estimate_type){
-  Zr <- effects_analysis %>%  pull({{Z_colname}})
-  VZr <- effects_analysis %>%  pull({{VZ_colname}})
-  mod <- ManyEcoEvo::fit_metafor_mv(estimate = Zr, 
-                        variance = VZr, 
-                        estimate_type = estimate_type, 
-                        data = effects_analysis)
-  return(mod)
-}
-
 back_transformed_predictions <- 
   ManyEcoEvo_yi %>% 
   dplyr::mutate(data = 
@@ -2967,7 +2956,10 @@ raw_mod_data_logged <-
 mod_data_logged <- raw_mod_data_logged %>% 
   mutate(MA_mod = 
            map(data, 
-               ~fit_MA_mv(.x, mean_log, std.error_log, "yi")))
+               ~ ManyEcoEvo::fit_MA_mv(effects_analysis = .x, 
+                                      Z_colname = mean_log, 
+                                      VZ_colname = std.error_log, 
+                                      estimate_type = "yi")))
 
 
 plot_data_logged <- mod_data_logged %>% 
diff --git a/supp_mat/SM2_EffectSizeAnalysis.qmd b/supp_mat/SM2_EffectSizeAnalysis.qmd
index d090ccb..8bf5461 100644
--- a/supp_mat/SM2_EffectSizeAnalysis.qmd
+++ b/supp_mat/SM2_EffectSizeAnalysis.qmd
@@ -717,17 +717,6 @@ plot_forest_2 <- function(data, intercept = TRUE, MA_mean = TRUE, y_zoom = numer
   print(p)
 }
 
-# TODO put into R/ and build into package to call!
-fit_MA_mv <- function(effects_analysis, Z_colname, VZ_colname, estimate_type){
-  Zr <- effects_analysis %>%  pull({{Z_colname}})
-  VZr <- effects_analysis %>%  pull({{VZ_colname}})
-  mod <- ManyEcoEvo::fit_metafor_mv(estimate = Zr, 
-                        variance = VZr, 
-                        estimate_type = estimate_type, 
-                        data = effects_analysis)
-  return(mod)
-}
-
 back_transformed_predictions <- 
   ManyEcoEvo_yi %>% 
   dplyr::mutate(data = 
@@ -765,7 +754,11 @@ raw_mod_data_logged <-
 mod_data_logged <- raw_mod_data_logged %>% 
   mutate(MA_mod = 
            map(data, 
-               ~fit_MA_mv(.x, mean_log, std.error_log, "yi")))
+               ~ ManyEcoEvo::fit_MA_mv(effects_analysis = .x,  
+                                      Z_colname = mean_log, 
+                                      VZ_colname = std.error_log,
+                                      estimate_type =  "yi"))
+         )
 
 
 plot_data_logged <- mod_data_logged %>% 
diff --git a/supp_mat/SM3_ExplainingDeviation.qmd b/supp_mat/SM3_ExplainingDeviation.qmd
index c256b65..a86b364 100644
--- a/supp_mat/SM3_ExplainingDeviation.qmd
+++ b/supp_mat/SM3_ExplainingDeviation.qmd
@@ -44,17 +44,6 @@ gt_fmt_yi <- function(gt_tbl, columns, ...) {
 
 
 ```{r define-fns, eval = TRUE}
-#TODO turn into own function and pull out of nested targets function and rm here
-fit_MA_mv <- function(effects_analysis, Z_colname, VZ_colname, estimate_type){
-  Zr <- effects_analysis %>%  pull({{Z_colname}})
-  VZr <- effects_analysis %>%  pull({{VZ_colname}})
-  mod <- ManyEcoEvo::fit_metafor_mv(estimate = Zr, 
-                                    variance = VZr, 
-                                    estimate_type = estimate_type, 
-                                    data = effects_analysis)
-  return(mod)
-}
-
 plot_forest <- function(data, intercept = TRUE, MA_mean = TRUE){
   if(MA_mean == FALSE){
     data <- filter(data, study_id != "overall")
@@ -398,7 +387,7 @@ raw_mod_data_logged <-
 mod_data_logged <- raw_mod_data_logged %>% 
   mutate(MA_mod = 
            map(data, 
-               ~fit_MA_mv(.x, mean_log, std.error_log, "yi")))
+               ~ ManyEcoEvo::fit_MA_mv(.x, mean_log, std.error_log, "yi")))
 
 deviation_models_yi_euc <- 
   raw_mod_data_logged %>%