rework getDeps/Parents to work at the node level

nimbleModel/R/MCMC_configuration.R

@@ -28,7 +28,7 @@ samplerConfClass <- R6Class(
       if (name == "crossLevel") {
         control <<- c(
           control,
-          list(dependent_nodes = getNodes(model, getDependencies(model$modelDef, target, self = FALSE, nodesAsChars = FALSE), stochOnly = TRUE, nodesAsChars = FALSE))
+          list(dependent_nodes = getNodes(model, getDependencies(model, target, self = FALSE, nodesAsChars = FALSE), stochOnly = TRUE, nodesAsChars = FALSE))
         )
       } # special case for printing dependents of crossLevel sampler (only)
     },

nimbleModel/R/MCMC_conjugacy.R

@@ -239,7 +239,7 @@ conjugacyRelationshipsClass <- R6Class(
 
       # CHECK: is this sufficiently efficient? Is getting a char representation the best strategy?
       targetNode <- getNodes(model, nodes = nodeRange$toNodeChars(1), nodesAsChars = FALSE)[[1]] # First node as representative.
-      deps <- getNodes(model, getDependencies(model$modelDef, targetNode, self = FALSE, nodesAsChars = FALSE), stochOnly = TRUE, nodesAsChars = FALSE)
+      deps <- getNodes(model, getDependencies(model, targetNode, self = FALSE, nodesAsChars = FALSE), stochOnly = TRUE, nodesAsChars = FALSE)
 
       depTypes <- sapply(deps, function(x) conjugacyObj$checkConjugacyOneDep(model, targetNode, x, restrictLink))
 
@@ -396,8 +396,8 @@ conjugacyClass <- R6Class(
     genSetupFunction = function(dependentCounts, doDependentScreen = FALSE) {
       functionBody <- codeBlockClass()
       functionBody$addCode({
-        calcNodes <- getNodes(model, getDependencies(model$modelDef, target, nodesAsChars = FALSE), nodesAsChars = FALSE)
-        calcNodesDeterm <- getNodes(model, getDependencies(model$modelDef, target, nodesAsChars = FALSE), determOnly = TRUE, nodesAsChars = FALSE)
+        calcNodes <- getNodes(model, getDependencies(model, target, nodesAsChars = FALSE), nodesAsChars = FALSE)
+        calcNodesDeterm <- getNodes(model, getDependencies(model, target, nodesAsChars = FALSE), determOnly = TRUE, nodesAsChars = FALSE)
       })
 
       # if this conjugate sampler is for a multivariate node (i.e., nDim > 0), then we need to determine the size (d)

nimbleModel/R/dataRules.R

@@ -151,7 +151,7 @@ makeRulePieces <- function(elements, varName, all, sequenceThreshold = 0.1) {
 }
 
 
-excludeFromPredictiveRules <- function(modelDef, currentRanges, candidateRules) {
+excludeFromPredictiveRules <- function(model, currentRanges, candidateRules) {
   if (!length(candidateRules)) {
     return(NULL)
   }
@@ -164,8 +164,8 @@ excludeFromPredictiveRules <- function(modelDef, currentRanges, candidateRules)
     } else {
       candidateRules[[varName]] <- NULL
     }
-    parents <- getParents(modelDef, range, nodesAsChars = FALSE)
-    candidateRules <- excludeFromPredictiveRules(modelDef, parents, candidateRules)
+    parents <- getParents(model, range, nodesAsChars = FALSE)
+    candidateRules <- excludeFromPredictiveRules(model, parents, candidateRules)
   }
   return(candidateRules)
 }

nimbleModel/R/graphRules.R

@@ -98,6 +98,7 @@ graphRuleClass <- R6Class(
       applyGraphRule(fromVarRange, self, removeDuplicates = removeDuplicates)
     },
     getFromRange = function() {
+      # Returns maximal extent; so for y[2:5], it returns 1:5; same for y[c(2,5)].
       if (!length(indexSets$fromIndexSlotToSet)) { # no indexing
         varRange <- varRangeClass$new(fromVarName)
       } else {

nimbleModel/R/modelBaseClass.R

@@ -2,6 +2,7 @@
 modelBase_nClass <- nClass(
   classname = "modelBase_nClass",
   Rpublic = list(
+    thisModel = NULL,
     modelDef = NULL,
     dataRules = NULL,
     nondataRules = NULL,
@@ -13,6 +14,8 @@ modelBase_nClass <- nClass(
       if (isTRUE(.GlobalEnv$.debugModelInit)) browser()
       super$initialize(...)
 
+      self$thisModel <- self   # A hack needed because we use `self` as method argument name below.
+
       # TODO: is there a better way to populate declFunNameToIndex in Cpublic?
       declFunNameToIndex <- self$declFunNameToIndex_
 
@@ -111,7 +114,7 @@ modelBase_nClass <- nClass(
           dataRule$rule$apply(dataRule$varName)
         })
       }))
-      self$predictiveRules <- excludeFromPredictiveRules(modelDef, dataRanges, candidateRules)
+      self$predictiveRules <- excludeFromPredictiveRules(self, dataRanges, candidateRules)
 
       # nonpredictive rules
       candidateRules <- unlist(lapply(modelDef$calcRules, function(oneVarRules) {
@@ -365,24 +368,36 @@ modelBase_nClass <- nClass(
         return(expr)
       }
     },
-    getDependencies = function(nodes, self = TRUE, downstream = FALSE, immediateOnly = FALSE,
+    # `self` arg masks the reference to the object.
+    # TODO: perhaps we should rename the arg `includeSelf`, but that is not back compatible.
+    getDependencies = function(nodes, self = TRUE, determOnly = FALSE, stochOnly = FALSE,
+                               includeData = TRUE, dataOnly = FALSE,
+                               includePredictive = nimble::getNimbleOption('getDependenciesIncludesPredictiveNodes'),
+                               predictiveOnly = FALSE, includeRHSonly = FALSE,
+                               downstream = FALSE, immediateOnly = FALSE,
                                nodesAsChars = getNimbleModelOption("nodesAsChars"),
                                returnScalarComponents = FALSE, .sort = FALSE) {
       nimbleModel::getDependencies(
-        modelDef, nodes, self, downstream, immediateOnly,
+        thisModel, nodes, self,
+        determOnly, stochOnly, includeData, dataOnly,
+        includePredictive, predictiveOnly, includeRHSonly,
+        downstream, immediateOnly,
         nodesAsChars, returnScalarComponents, .sort
       )
     },
-    getParents = function(nodes, self = FALSE, upstream = FALSE, immediateOnly = FALSE,
+    getParents = function(nodes, self = FALSE,
+                          determOnly = FALSE, stochOnly = FALSE,
+                          includeData = TRUE, dataOnly = FALSE, includeRHSonly = FALSE,
+                          upstream = FALSE, immediateOnly = FALSE,
                           nodesAsChars = getNimbleModelOption("nodesAsChars"),
                           returnScalarComponents = FALSE, .sort = FALSE) {
       nimbleModel::getParents(
-        modelDef, nodes, self, upstream, immediateOnly,
+        thisModel, nodes, self,
+        determOnly, stochOnly, includeData, dataOnly, includeRHSonly,        
+        upstream, immediateOnly,
         nodesAsChars, returnScalarComponents, .sort
       )
     },
-    # TODO: not working because `nimbleModel::getNodes` needs the model not just modelDef.
-    # Once we integrate modelClass with modelBase_nClass, we should be able to pass `self`.
     getNodes = function(nodes, determOnly = FALSE, stochOnly = FALSE,
                         includeData = TRUE, dataOnly = FALSE,
                         includeRHSonly = FALSE,

nimbleModel/R/modelDef.R

@@ -1024,52 +1024,6 @@ modelDefClass <- R6Class(
 )
 
 
-# Core graph and node querying functions in the model API.
-# These are standalone functions for now, but may become
-# part of model class. That said, more naturally part of modelDef class.
-
-# TODO: move these functions into a new stand-alone code file for user-facing functions?
-
-# Note: `getDependencies` and `getParents` cannot handle `stochOnly` or `determOnly`
-# because a given varRange result for getParents could be partially stochastic and
-# partially deterministic. Instead a user would pass the result through `getNodes()`.
-# Similarly, filtering by RHSonly will be done in `getNodes()`.
-
-# Note: data-related flags not handled as that relates to flags on a model
-# and not part of modelDef.
-
-# TODO: these should presumably take the model not modelDef as the first arg.
-# Once we integrate modelClass with modelBase_nClass, we should be able to
-# pass `self` from the getDeps and getParents methods to these functions.
-
-getDependencies <- function(modelDef, nodes,
-                            self = TRUE,
-                            downstream = FALSE, immediateOnly = FALSE,
-                            nodesAsChars = getNimbleModelOption("nodesAsChars"),
-                            returnScalarComponents = FALSE, .sort = FALSE) {
-  traverseGraph(modelDef$downstreamRules, modelDef$declRules,
-    nodes = nodes,
-    down = TRUE, self = self,
-    follow = downstream, immediateOnly = immediateOnly,
-    nodesAsChars = nodesAsChars, returnScalarComponents = returnScalarComponents,
-    .sort = .sort, modelDef = modelDef
-  )
-}
-
-getParents <- function(modelDef, nodes,
-                       self = FALSE,
-                       upstream = FALSE, immediateOnly = FALSE,
-                       nodesAsChars = getNimbleModelOption("nodesAsChars"),
-                       returnScalarComponents = FALSE, .sort = FALSE) {
-  traverseGraph(modelDef$upstreamRules, modelDef$declRules,
-    nodes = nodes,
-    down = FALSE, self = self,
-    follow = upstream, immediateOnly = immediateOnly,
-    nodesAsChars = nodesAsChars, returnScalarComponents = returnScalarComponents,
-    .sort = .sort, modelDef = modelDef
-  )
-}
-
 
 # Evaluates `if` statements in model code to generate actual model code
 # without any `if` statements. Condition of if statement can use variables

nimbleModel/R/modelFunctions.R

@@ -216,6 +216,46 @@ taggedClass <- R6Class(
   )
 )
 
+getDependencies <- function(model, nodes,
+                            self = TRUE, determOnly = FALSE, stochOnly = FALSE,
+                            includeData = TRUE, dataOnly = FALSE,
+                            includePredictive = nimble::getNimbleOption('getDependenciesIncludesPredictiveNodes'),
+                            predictiveOnly = FALSE, includeRHSonly = FALSE,
+                            downstream = FALSE, immediateOnly = FALSE,
+                            nodesAsChars = getNimbleModelOption("nodesAsChars"),
+                            returnScalarComponents = FALSE, .sort = FALSE) {
+  traverseGraph(model$modelDef$downstreamRules, model$modelDef$declRules,
+    nodes = nodes,
+    down = TRUE, self = self,
+    determOnly = determOnly, stochOnly = stochOnly,
+    includeData = includeData, dataOnly = dataOnly, includePredictive = includePredictive,
+    predictiveOnly = predictiveOnly, includeRHSonly = includeRHSonly,
+    follow = downstream, immediateOnly = immediateOnly,
+    nodesAsChars = nodesAsChars, returnScalarComponents = returnScalarComponents,
+    .sort = .sort, model = model
+  )
+}
+
+getParents <- function(model, nodes,
+                       self = FALSE, determOnly = FALSE, stochOnly = FALSE,
+                       includeData = TRUE, dataOnly = FALSE, includeRHSonly = FALSE,
+                       upstream = FALSE, immediateOnly = FALSE,
+                       nodesAsChars = getNimbleModelOption("nodesAsChars"),
+                       returnScalarComponents = FALSE, .sort = FALSE) {
+  traverseGraph(model$modelDef$upstreamRules, model$modelDef$declRules,
+    nodes = nodes,
+    down = FALSE, self = self,
+    determOnly = determOnly, stochOnly = stochOnly,
+    includeData = includeData, dataOnly = dataOnly, includePredictive = TRUE,
+    predictiveOnly = FALSE, includeRHSonly = includeRHSonly,
+    follow = upstream, immediateOnly = immediateOnly,
+    nodesAsChars = nodesAsChars, returnScalarComponents = returnScalarComponents,
+    .sort = .sort, model = model
+  )
+}
+
+
+
 # This may not optimally aggregate in cases without contiguity - e.g., 2:4 + 6:8 will become a matrix, even though
 # it may be more efficient to leave it as two nodeRanges.
 #' @export
@@ -236,6 +276,7 @@ aggregate_nodes <- function(nodeSet) {
   names(nodeSet) <- declIDs
   nodeIDs <- lapply(nodeSet, \(x) x$getIDs())
   IDsByDecl <- lapply(split(nodeIDs, declIDs), \(x) unique(nimbleModel:::flatten(x)))
+  IDsByDecl <- IDsByDecl[unique(declIDs)]  # Try to keep in order provided.
   nms <- names(IDsByDecl)
   newNodeSet <- lapply(seq_along(IDsByDecl), \(i) {
     if(sum(nms[i] == declIDs) > 1) {
@@ -288,6 +329,16 @@ setdiff_nodes <- function(nodeSet1, nodeSet2) {
   }
   declIDs1 <- sapply(nodeSet1, \(x) x$decl$declRule$ID)
   declIDs2 <- sapply(nodeSet2, \(x) x$decl$declRule$ID)
+
+  nullCases <- sapply(declIDs1, is.null)
+  RHSonly <- nodeSet1[nullCases]
+  nodeSet1 <- nodeSet1[!nullCases]
+  declIDs1 <- unlist(declIDs1[!nullCases])
+
+  nullCases <- sapply(declIDs2, is.null)
+  nodeSet2 <- nodeSet2[!nullCases]
+  declIDs2 <- unlist(declIDs2[!nullCases])
+
   nodeIDs1 <- lapply(nodeSet1, \(x) x$getIDs())
   nodeIDs2 <- lapply(nodeSet2, \(x) x$getIDs())
   excludeNodeIDs <- lapply(split(nodeIDs2, declIDs2), \(x) unique(nimbleModel:::flatten(x)))
@@ -301,7 +352,7 @@ setdiff_nodes <- function(nodeSet1, nodeSet2) {
       )
     }
   }
-  return(newNodeSet1[!sapply(newNodeSet1, is.null)])
+  return(c(newNodeSet1[!sapply(newNodeSet1, is.null)], RHSonly))
 }
 
 #' @export
@@ -922,8 +973,7 @@ splitLatents <- function(model, paramNodes, latentNodes, calcNodes, calcNodesOth
     latentNodes <- margNodes$randomEffectsNodes
     deps <- model$getNodes(model$getDependencies(latentNodes, self = FALSE, nodesAsChars = FALSE), includeData = FALSE, nodesAsChars = FALSE)
     ## By default, we treat "siblings" of latent nodes as latents.
-    ## This attempts to have fixed effects in latents,
-    ## along with random effects.
+    ## This attempts to have fixed effects in latents, along with random effects.
     newLatents <- model$getNodes(model$getParents(deps, nodesAsChars = FALSE), stochOnly = TRUE, includeData = FALSE, nodesAsChars = FALSE)
     paramNodes <- setdiff_nodes(paramNodes, newLatents)
     latentNodes <- aggregate_nodes(c(latentNodes, newLatents))