Update v61_pedigree_model_fitting.Rmd

Author: smasongarrison

vignettes/articles/v61_pedigree_model_fitting.Rmd

@@ -186,6 +186,7 @@ multi_model_ce <- buildPedigreeMx(
 ```{r, include = FALSE}
 fitted_multi_ace <- NULL # ensure variable exists for later use
 fitted_multi_mace <- NULL # ensure variable exists for later use
+fitted_multi_ce <- NULL # ensure variable exists for later use
 ```
 ```{r, eval = run_models}
 fitted_multi_ace <- mxRun(multi_model_ace,
@@ -199,15 +200,33 @@ fitted_multi_ce <- mxRun(multi_model_ce, intervals= TRUE)
 saveRDS(fitted_multi_ce, "inst/extdata/fitted_multi_ce.rds")
 ```
 
+```{r, eval = FALSE, include = FALSE}
+tt_ace <- mxTryHard(fitted_multi_ace, intervals = TRUE)
+saveRDS(tt_ace, "inst/extdata/fitted_multi_ace_tryhard.rds")
+tt_mace <- mxTryHard(fitted_multi_mace, intervals = TRUE)
+saveRDS(tt_mace, "inst/extdata/fitted_multi_mace_tryhard.rds")
+
+tt_ce <- mxTryHard(fitted_multi_ce, intervals = TRUE)
+saveRDS(tt_ce, "inst/extdata/fitted_multi_ce_tryhard.rds")
+```
+
+Note that fitting these models can take some time, especially with many families and large pedigrees. The `mxTryHard()` function can be used to attempt to find better-fitting solutions if the initial optimization does not converge well. In practice, you may want to experiment with different starting values or optimization settings to improve convergence.
+
 ```
 fitted_multi_ace <- mxRun(multi_model_ace)
 fitted_multi_mace <- mxRun(multi_model_mace)
+fitted_multi_ce <- mxRun(multi_model_ce)
 ```
-```{r, include = FALSE, eval = file.exists(rds_file_ace) && !run_models && file.exists(rds_file_mace)}
+
+```{r, include = FALSE, eval = file.exists(rds_file_ace) && !run_models && file.exists(rds_file_mace) && file.exists(file.path(rds_dir, "fitted_multi_ce.rds"))}
 fitted_multi_ace <- readRDS(rds_file_ace)
 fitted_multi_mace <- readRDS(rds_file_mace)
 fitted_multi_ce <- readRDS(file.path(rds_dir, "fitted_multi_ce.rds"))
 ```
+
+
+
+
 ```{r, eval = has_openmx && !is.null(fitted_multi_ace)}
 summary(fitted_multi_ace,verbose=T)
 summary(fitted_multi_ace)$CI
@@ -246,17 +265,20 @@ cat("Common nuclear  (Vcn):", fitted_multi_mace$ModelOne$Vcn$values/total_var_ma
 cat("Mitochondrial (Vmt):", fitted_multi_mace$ModelOne$Vmt$values/total_var_mace, "\n")
 cat("Unique environ. (Ver):", fitted_multi_mace$ModelOne$Ver$values/total_var_mace, "\n")
 ```
-Now we can compare the ACE and MACE models to see if including the mitochondrial component changes our estimates of heritability and environmental contributions. This can provide insights into the role of mitochondrial inheritance in the trait of interest (LRS) in these squirrels.
+Now we can compare the ACE and MACE models to see if including the mitochondrial component changes our estimates of heritability and environmental contributions. This can provide insights into the role of mitochondrial inheritance in the trait of interest (LRS) in these squirrels. In the original paper, the authors found that a maternally inherited component (which could be due to mitochondria or maternal effects) explained a significant portion of the variance in LRS,. 
 
-```{r, eval = has_openmx && !is.null(fitted_multi_ace) && !is.null(fitted_multi_mace)}
+```{r, eval = has_openmx && !is.null(fitted_multi_ace) && !is.null(fitted_multi_mace) && !is.null(fitted_multi_ce)}
 
 mxCompare(fitted_multi_mace, fitted_multi_ace)
+mxCompare(fitted_multi_ace, fitted_multi_ce)
 
-tt_ace <- mxTryHard(fitted_multi_ace, intervals = TRUE)
-saveRDS(tt_ace, "inst/extdata/fitted_multi_ace_tryhard.rds")
-tt_mace <- mxTryHard(fitted_multi_mace, intervals = TRUE)
-saveRDS(tt_mace, "inst/extdata/fitted_multi_mace_tryhard.rds")
 
-tt_ce <- mxTryHard(fitted_multi_ce, intervals = TRUE)
-saveRDS(tt_ce, "inst/extdata/fitted_multi_ce_tryhard.rds")
 ```
+
+However, as you can see when we compare the ACE and MACE models, the inclusion of the mitochondrial component does not substantially change the estimates of additive genetic and common nuclear environmental variance, suggesting that the mitochondrial component may not be a major contributor to LRS in this dataset.
+
+
+
+
+
+