From 626a6b4790efa9c989be8b10a1ba3250b4346275 Mon Sep 17 00:00:00 2001
From: Jim Thorson <James.T.Thorson@gmail.com>
Date: Wed, 5 Feb 2025 17:21:04 -0800
Subject: [PATCH] add `total_effect(.)` to vignette

---
 R/utility.R            | 22 ++++++++++++----------
 man/total_effect.Rd    |  8 ++++----
 vignettes/vignette.Rmd | 13 +++++++++++++
 3 files changed, 29 insertions(+), 14 deletions(-)

diff --git a/R/utility.R b/R/utility.R
index 8235d2a..3f1390c 100644
--- a/R/utility.R
+++ b/R/utility.R
@@ -134,41 +134,43 @@ function( object,
 #' from 0 (simultaneous effects) to a user-specified maximum lag.
 #'
 #' @param object Output from \code{\link{dsem}}
-#' @param n_times Number of lags over which to calculate total effects
+#' @param n_lags Number of lags over which to calculate total effects
 #'
 #' @details
 #' Total effects are calculated from the Leontief matrix \eqn{\mathbf{(I-P)^{-1}}
 #'
 #' @return
-#' A data frame listing the time-lag (deltaT), variable that is undergoing some 
+#' A data frame listing the time-lag (lag), variable that is undergoing some 
 #' exogenous change (from), and the variable being impacted (to), along with the 
-#' total effect including direct and indirect pathways
+#' total effect (effect) including direct and indirect pathways
 #'
 #' @export
 total_effect <-
 function( object,
-          n_times = 4 ){
+          n_lags = 4 ){
 
-  #library(Matrix)
+  # Extract path matrix
   Z = object$internal$tsdata
   P_kk = make_matrices( 
     beta_p = object$internal$parhat$beta,
     model = object$sem_full,
-    times = seq_len(n_times),
+    times = seq_len(n_lags),
     variables = colnames(Z)
   )$P_kk            
 
-  #
+  # Define innovations
   delta_kj = kronecker( Diagonal(n=ncol(Z)), 
-                        sparseMatrix(i=1, j=1, x=1, dims=c(n_times,1)) )
+                        sparseMatrix(i=1, j=1, x=1, dims=c(n_lags,1)) )
   IminusRho_kk = Diagonal(n=nrow(P_kk)) - P_kk
+  
+  # Calculate total effect using sparse matrices
   Impact_kj = solve( IminusRho_kk, delta_kj )
   
   # Make into data frame
-  out = expand.grid( "deltaT" = seq_len(n_times)-1, 
+  out = expand.grid( "lag" = seq_len(n_lags)-1, 
                      "to" = colnames(Z), 
                      "from" = colnames(Z) )
-  out$total_effect = as.vector(Impact_kj)
+  out$effect = as.vector(Impact_kj)
   return(out)
 }
 
diff --git a/man/total_effect.Rd b/man/total_effect.Rd
index a49b945..47f93aa 100644
--- a/man/total_effect.Rd
+++ b/man/total_effect.Rd
@@ -4,17 +4,17 @@
 \alias{total_effect}
 \title{Calculate total effects}
 \usage{
-total_effect(object, n_times = 4)
+total_effect(object, n_lags = 4)
 }
 \arguments{
 \item{object}{Output from \code{\link{dsem}}}
 
-\item{n_times}{Number of lags over which to calculate total effects}
+\item{n_lags}{Number of lags over which to calculate total effects}
 }
 \value{
-A data frame listing the time-lag (deltaT), variable that is undergoing some 
+A data frame listing the time-lag (lag), variable that is undergoing some 
 exogenous change (from), and the variable being impacted (to), along with the 
-total effect including direct and indirect pathways
+total effect (effect) including direct and indirect pathways
 }
 \description{
 Calculate a data frame of total effects, representing the 
diff --git a/vignettes/vignette.Rmd b/vignettes/vignette.Rmd
index 97443fc..2023318 100644
--- a/vignettes/vignette.Rmd
+++ b/vignettes/vignette.Rmd
@@ -323,6 +323,19 @@ ggarrange( p1, p2, p3,
            ncol = 1, nrow = 3)
 ```
 
+We can then plot the total effects, which shows that effects propagate through time due to both interactions and density dependence:
+
+```{r, echo=TRUE, message=FALSE, fig.width=5, fig.height=7, warning=FALSE}
+# Calculate total effects
+effect = total_effect( fit )
+
+# Plot total effect
+ggplot( effect) + 
+  geom_bar( aes(lag, effect, fill=lag), stat='identity', col='black', position='dodge' ) + 
+  facet_grid( from ~ to  )
+```
+
+
 Results again show that `dsem` can estimate parameters for a vector autoregressive model (VAM), and it exactly matches results from `vars`, using `dynlm`, or using `MARSS`.
 
 ## Multi-causal ecosystem synthesis