updates to 1.1.0 ...

... to accompany accepted paper at Fish and Fisheries

.Rbuildignore

@@ -3,6 +3,6 @@
 ^\.github$
 ^cran-comments\.md$
 ^CRAN-SUBMISSION$
-
+^data-raw$
 
 

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: phylosem
 Type: Package
 Title: Phylogenetic Structural Equation Model
-Version: 1.0.2
-Date: 2023-09-22
+Version: 1.1.0
+Date: 2023-10-03
 Authors@R: 
     c(person(given = "James",
              family = "Thorson",
@@ -28,6 +28,7 @@ Suggests:
     TreeTools,
     Rphylopars,
     phylolm,
+    fishtree,
     phyr,
     knitr,
     rmarkdown,
@@ -45,7 +46,7 @@ Description: Applies phylogenetic comparative methods (PCM) and phylogenetic tra
     allow users to easily read all of the documentation and source code.  PCM using SEM 
     includes phylogenetic linear models and structural equation models as nested submodels, 
     but also allows imputation of missing values.  Features and comparison with other packages
-    are described in Thorson and van der Bijl (In revisions). 
+    are described in Thorson and van der Bijl (In press). 
 License: GPL-3
 Encoding: UTF-8
 RoxygenNote: 7.2.3

NAMESPACE

@@ -1,12 +1,13 @@
 # Generated by roxygen2: do not edit by hand
 
-S3method(AIC,phylosem)
 S3method(average,compare_phylosem)
 S3method(best,compare_phylosem)
 S3method(choice,compare_phylosem)
 S3method(coef,phylosem)
+S3method(logLik,phylosem)
 S3method(print,phylosem)
 S3method(summary,phylosem)
+S3method(vcov,phylosem)
 export(TMBAIC)
 export(as_fitted_DAG)
 export(as_phylo4d)
@@ -33,10 +34,12 @@ importFrom(sem,sem)
 importFrom(sem,specifyEquations)
 importFrom(sem,specifyModel)
 importFrom(stats,AIC)
+importFrom(stats,logLik)
 importFrom(stats,na.omit)
 importFrom(stats,nlminb)
 importFrom(stats,optimHess)
 importFrom(stats,plogis)
 importFrom(stats,pnorm)
 importFrom(stats,rnorm)
+importFrom(stats,vcov)
 useDynLib(phylosem, .registration = TRUE)

NEWS.md

@@ -1,3 +1,10 @@
+# phylosem 1.1.0
+
+* Adding new S3 generic functions: logLik and vcov
+* Adding vignette "fisheries" associated with paper
+  "Trees for fishes: The neglected role for phylogenetic comparative
+  methods in fisheries science" accepted at Fish and Fisheries
+
 # phylosem 1.0.2
 
 * Fixed bug arising when data are missing for observations in the earliest tips and for the first inputted variable

R/data.R

@@ -0,0 +1,14 @@
+
+#' Fisheries natural mortality example
+#'
+#' Data used to demonstrate phylogenetic comparative methods for fisheries science. 
+#' Specifically a copy of the Then et al. database \doi{10.1093/icesjms/fsu136} 
+#' using file "Mlifehist_ver1.0.csv" accessed from 
+#' \url{https://www.vims.edu/research/departments/fisheries/programs/mort_db/index.php}
+#'
+#' @name Mlifehist_ver1_0
+#' @docType data
+#' @usage data(Mlifehist_ver1_0)
+#' @keywords data
+NULL
+

R/phylosem.R

@@ -72,9 +72,9 @@
 #' **Introducing the package, its features, and comparison with other software
 #' (to cite when using phylosem):**
 #'
-#' Thorson, J. T., & van der Bijl, W. (In revision). phylosem: A fast and simple
+#' Thorson, J. T., & van der Bijl, W. (In press). phylosem: A fast and simple
 #' R package for phylogenetic inference and trait imputation using phylogenetic
-#' structural equation models.
+#' structural equation models. Journal of Evolutionary Biology.
 #'
 #' *Statistical methods for phylogenetic structural equation models*
 #'
@@ -347,12 +347,12 @@ function( sem,
                         random = tmb_inputs$random,
                         DLL = "phylosem" )
   if(quiet==FALSE) list_parameters(obj)
-  results = list( data=data,
-                  SEM_model=SEM_model,
-                  obj=obj,
-                  call=match.call(),
-                  tree=tree,
-                  tmb_inputs=tmb_inputs )
+  results = list( "data" = data,
+                  "SEM_model" = SEM_model,
+                  "obj" = obj,
+                  "call" = match.call(),
+                  "tree" = tree,
+                  "tmb_inputs" = tmb_inputs )
 
   # Export stuff
   if( run_model==FALSE ){
@@ -421,19 +421,59 @@ coef.phylosem = function( object, standardized=FALSE, ... ){
   return( data.frame(Path=object$SEM_model[,1], Parameter=object$SEM_model[,2], Estimate=SEM_params ) )
 }
 
-#' Calculate AIC
+#' Extract Variance-Covariance Matrix
 #'
-#' @inheritParams TMBAIC
+#' extract the covariance of fixed effects, or both fixed and random effects.
 #'
-#' @title Calculate Akaike Information Criterion from marginal likelihood
-#'
-#' @param object Output from \code{\link{phylosem}}
-#' @param ... Not used
-#' @return Akaike Information Criterion
-#' @method AIC phylosem
+#' @param object output from \code{phylosem}
+#' @param which whether to extract the covariance among fixed effects, random effects, or both
+#' @param ... ignored, for method compatibility
+#' @importFrom stats vcov
+#' @method vcov phylosem
 #' @export
-AIC.phylosem = function( object, ..., k = 2 ){
-  return( TMBAIC(object$opt, ..., k=k) )
+vcov.phylosem <-
+function( object,
+          which = c("fixed", "random", "both"),
+          ...) {
+
+  which = match.arg(which)
+
+  if( which=="fixed" ){
+    V = object$opt$SD$cov.fixed
+    if(is.null(V)){
+      warning("Please re-run `phylosem` with `getsd=TRUE`, or confirm that the model is converged")
+    }
+  }
+  if( which=="random" ){
+    V = solve(object$obj$env$spHess(random=TRUE))
+  }
+  if( which=="both" ){
+    H = object$opt$SD$jointPrecision
+    if(is.null(H)){
+      warning("Please re-run `phylosem` with `getsd=TRUE` and `getJointPrecision=TRUE`, or confirm that the model is converged")
+      V = NULL
+    }else{
+      V = solve(H)
+    }
+  }
+
+  return( V )
+}
+
+# Extract the (marginal) log-likelihood of a phylosem model
+#
+# @return object of class \code{logLik} with attributes
+#   \item{val}{log-likelihood}
+#   \item{df}{number of parameters}
+#' @importFrom stats logLik
+#' @export
+logLik.phylosem <- function(object, ...) {
+  val = -1 * object$opt$objective
+  df = length( object$opt$par )
+  out = structure( val,
+             df = df,
+             class = "logLik")
+  return(out)
 }
 
 #' summarize phylosem

data-raw/Mlifehist_ver1.0.R

@@ -0,0 +1,10 @@
+
+data_dir =  R'(C:\Users\James.Thorson\Desktop\Git\phylosem\data-raw)'
+
+# Load and format
+Mlifehist_ver1_0 = read.csv( file.path(data_dir,"Mlifehist_ver1.0.csv") )
+Mlifehist_ver1_0 = Mlifehist_ver1_0[,-match(c("GrowthRef","Mref","tmaxRef","SpecificLocation","CommonName"),names(Mlifehist_ver1_0))]
+
+# Export
+setwd( R'(C:\Users\James.Thorson\Desktop\Git\phylosem)' )
+usethis::use_data( Mlifehist_ver1_0, overwrite=TRUE )