Skip to content

Deal with non-constant block indices #38

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
paciorek merged 4 commits into from
Jul 13, 2026
Merged

Deal with non-constant block indices #38

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
febf3fc
Clean up printing for nodeRanges.
paciorek Jul 9, 2026
510942d
Handle varying constant block indices,
paciorek Jul 10, 2026
98674e0
Fix some trivial test failures.
paciorek Jul 10, 2026
d4b6f6e
Tweak misbehaving test.
paciorek Jul 13, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
9 changes: 8 additions & 1 deletion nimbleModel/R/indexRuleAll.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -97,7 +97,14 @@ indexRuleAll_setup <- function(toIndexExprList,
unrolledSize <- unrolledIndicesEnv$outputSize

if (length(unrolledResults) > 1) {
allIndices <- do.call("cbind", unrolledResults)
lists <- sapply(unrolledResults, is.list)
if(any(lists)) { # Cases like `y[i,n1[i]:n2[i]]`. We unroll these.
allIndices <- do.call("rbind", lapply(seq_along(unrolledResults[[1]]), \(i)
as.matrix(do.call("expand.grid",
lapply(unrolledResults, \(x) x[[i]])
))
))
} else allIndices <- do.call("cbind", unrolledResults)
} else {
allIndices <- matrix(unrolledResults[[1]], nrow = unrolledSize)
}
Expand Down
7 changes: 4 additions & 3 deletions nimbleModel/R/modelDecl.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,8 @@ modelDeclClass <- R6Class(
context,
sourceLineNumber,
truncated = FALSE,
boundExprs = NULL) {
boundExprs = NULL,
warnNonConstantBlocks = FALSE) {
context <<- context
sourceLineNumber <<- sourceLineNumber
code <<- code
Expand Down Expand Up @@ -69,8 +70,8 @@ modelDeclClass <- R6Class(
targetVarExpr <- targetExpr[[2]]
targetNodeExpr <<- targetExpr
indexedBlocks <- sapply(indexExpr, checkForIndexedIntervals, context)
if (any(indexedBlocks)) {
stop("Non-constant indexing found in `", safeDeparse(indexExpr[[which(indexedBlocks)[1]]]), "` in `", safeDeparse(targetNodeExpr), "`. Block indices must be constant.")
if (any(indexedBlocks) && warnNonConstantBlocks) {
messageIfVerbose("Non-constant indexing found in `", safeDeparse(indexExpr[[which(indexedBlocks)[1]]]), "` in `", safeDeparse(targetNodeExpr), "`. Graph queries will incorrectly show elements of multivariate nodes as individual nodes but other functionality should operate correctly.")
}
} else {
# There is a transformation, possibly with a subscript.
Expand Down
3 changes: 2 additions & 1 deletion nimbleModel/R/modelDef.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -192,7 +192,8 @@ modelDefClass <- R6Class(
declInfo[[iAns]] <<- modelDeclClass$new(
code[[i]],
contexts[[contextID]],
lineNumber
lineNumber,
warnNonConstantBlocks = TRUE
)
}
if (code[[i]][[1]] == "for") {
Expand Down
4 changes: 2 additions & 2 deletions nimbleModel/R/nodeRules.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -497,7 +497,7 @@ nodeRangeClass <- R6Class(

indexVars <- paste0("idx", seq_len(sum(boolExternalIndexRanges & !multiSlot)))
exprs[boolExternalIndexRanges & !multiSlot] <- lapply(indexVars, as.name)
forText <- paste(indexVars, "in", sv)
forText <- paste0("`", indexVars, "`", " in ", sv)
result <- safeDeparse(do.call("call",
c(
list("[", nm),
Expand All @@ -512,7 +512,7 @@ nodeRangeClass <- R6Class(
forText <- "`"
}
result <- paste0("`", result, forText)
}
} else result <- paste0("`", result, "`")
return(result)
},

Expand Down
10 changes: 9 additions & 1 deletion nimbleModel/R/originalIndexingRules.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,11 +38,19 @@ originalIndexingRuleClass <- R6Class(

# Produces a varRange, though it's not really a range for a variable
# but rather a range for the indices.
# (2023-06-10, commit 30ede6)
# Do not remove duplicates because in generation of `calcRange`s there
# can be cases where we need duplicated values in order to have correct
# number of logProbs.
# (2026-07-09) actually we need to remove duplicates for y[i,n1[i]:n2[i]] case
# as that has one 'node' per scalar element.
# It appears that the 2023-06-10 thinking was incorrect as it seems to have
# been based on wanting `calculate()` to return as many elements as were
# computed in any given deterministic calculation, e.g., y[1:2] <- foo()
# returning 2 elements, but we actually don't return deterministic results anyway.
# See test in line 1302 of test-nodeRules.R in commit 30ede6.
apply = function(fromVarRange) {
graphRule$apply(fromVarRange, removeDuplicates = FALSE)
graphRule$apply(fromVarRange, removeDuplicates = TRUE)
}
)
)
16 changes: 0 additions & 16 deletions nimbleModel/tests/testthat/test-modelDefClass.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -282,19 +282,3 @@ test_that("detection of duplicated declarations", {

})

test_that("detection of non-constant block indices", {
code <- quote({
for(i in 1:2) {
x[i:5] ~ dmnorm(z[i:5], Q[i:5,i:5])
}
})
expect_error(m <- modelClass$new(code), "Non-constant indexing")

code <- quote({
for(i in 1:2) {
x[i,(i+1):5] ~ dmnorm(z[(i+1):5], Q[(i+1):5,(i+1):5])
}
})
expect_error(m <- modelClass$new(code), "Non-constant indexing")

})
8 changes: 3 additions & 5 deletions nimbleModel/tests/testthat/test-modelGraph.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1548,16 +1548,14 @@ test_that("basic hierarchical models", {
})

test_that("mixed-length block dependences", {
# See issue 31 regarding getting scalar elements of nodes in such cases.
code <- quote({
for(i in 1:3)
y[i, n1[i]:n2[i]] ~ dmulti(p[n1[i]:n2[i]], 10)
})
# Issue #31
expect_error({
modelDef <- modelDefClass$new(code, constants = list(n1 = c(3,1,2), n2 = c(6,3,2)));
result <- getDependencies(modelDef, 'p[2]');
modelDef <- modelDefClass$new(code, constants = list(n1 = c(3,1,2), n2 = c(6,3,2)))
result <- getDependencies(modelDef, 'p[2]')
expect_equal(result[[1]]$indexRanges, list(newIndexRange(matrix(c(2,2,2,3,1,2,3,2), ncol = 2))))
})
})


Expand Down
177 changes: 162 additions & 15 deletions nimbleModel/tests/testthat/test-nimbleModel.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -659,7 +659,7 @@ test_that("non-sequential indexing cases", {
nr <- m$getNodes()[[1]]
expect_true(inherits(nr$indexRanges[[1]], "indexRangeMatrixClass"))
expect_identical(nr$numExternalIndexRanges, 1L)
expect_identical(nr$toChar(), "`y[idx1]`, for idx1 in c(2, 3, 5)")
expect_identical(nr$toChar(), "`y[idx1]`, for `idx1` in c(2, 3, 5)")

code <- nimbleCode({
y[c(2,3,5)] ~ dmnorm(mu[1:3], pr[1:3,1:3])
Expand All @@ -686,22 +686,23 @@ test_that("non-sequential indexing cases", {
nr <- m$getNodes()[[2]]
expect_true(inherits(nr$indexRanges[[1]], "indexRangeMatrixClass"))
expect_identical(nr$numExternalIndexRanges, 0L)
expect_identical(nr$toChar(), "y[c(2, 3, 5)]")
expect_identical(nr$toChar(), "`y[c(2, 3, 5)]`")

if(FALSE) { # No operator def for nimC: This was part of the call: y[i = nimC(2, 3, 5)]
code <- nimbleCode({
# Cannot compile. No operator def for nimC: This was part of the call: y[i = nimC(2, 3, 5)]
code <- nimbleCode({
y[c(2,3,5)] <- x[1:3] + 1
})
mclass <- nimbleModel(code, inits = list(x = 1:3),
returnClass = TRUE)
m <- mclass$new()
cmclass <- nCompile(mclass)
cm <- cmclass$new()
m$calculate()
cm$calculate()
expect_identical(m$y, c(NA,2,3,NA,4))
expect_identical(cm$y, c(NA,2,3,NA,4))
}
})
mclass <- nimbleModel(code, inits = list(x = 1:3),
returnClass = TRUE)
m <- mclass$new()
m$calculate()
expect_identical(m$y, c(NA,2,3,NA,4))
expect_error({
cmclass <- nCompile(mclass)
cm <- cmclass$new()
cm$calculate()
expect_identical(cm$y, c(NA,2,3,NA,4))
})

})

Expand Down Expand Up @@ -1406,6 +1407,152 @@ test_that("isData", {

})

test_that("non-constant block indexing", {
# Changing dimensions.
code <- quote({
for(i in 1:3)
y[i, n1[i]:n2[i]] ~ dmnorm(mu[n1[i]:n2[i]], pr[n1[i]:n2[i],n1[i]:n2[i]])
})
model <- nimbleModel(code, data=list(y=matrix(rnorm(7*6),7)), constants = list(mu=rep(0,6), pr=diag(6),n1 = c(3,1,2), n2 = c(6,3,3)))
expect_identical(model$getNodes('y',nodesAsChars=TRUE),
c(paste0("y[1, ", 3:6, "]"),
paste0("y[2, ", 1:3, "]"),
paste0("y[3, ", 2:3, "]")))

truth <- dmnorm_chol(model$y[1,3:6],rep(0,4),diag(4), log=TRUE)+
dmnorm_chol(model$y[2,1:3],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[3,2:3],rep(0,2),diag(2), log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth) # Formerly incorrect as discussed in issue 31.

# All same dimensions, constant precision.
code <- quote({
for(i in 1:3)
y[i, n1[i]:n2[i]] ~ dmnorm(mu[n1[i]:n2[i]], pr[1:3,1:3])
})
model <- nimbleModel(code, data=list(y=matrix(rnorm(7*6),7)), constants = list(mu=rep(0,6), pr=diag(6),n1 = c(3,1,2), n2 = c(5,3,4)))
expect_identical(model$getNodes('y',nodesAsChars=TRUE),
c(paste0("y[1, ", 3:5, "]"),
paste0("y[2, ", 1:3, "]"),
paste0("y[3, ", 2:4, "]")))
truth <- dmnorm_chol(model$y[1,3:5],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[2,1:3],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[3,2:4],rep(0,3),diag(3), log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)

# Various crazy cases.
code <- quote({
for(i in 1:3)
y[n1[i]:n2[i], k[i], i] ~ dmnorm(mu[n1[i]:n2[i]], pr[1:3,1:3])
})
model <- nimbleModel(code, data=list(y=array(rnorm(5*5*3),c(5,5,3))), constants = list(k = c(2,4,5), mu=rep(0,6), pr=diag(6), n1 = c(3,1,2), n2 = c(5,3,4)))
expect_identical(sort(model$getNodes('y',nodesAsChars=TRUE)),
sort(c(paste0("y[", 3:5, ", 2, 1]"),
paste0("y[", 1:3, ", 4, 2]"),
paste0("y[", 2:4, ", 5, 3]"))))
truth <- dmnorm_chol(model$y[3:5,2,1],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[1:3,4,2],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[2:4,5,3],rep(0,3),diag(3), log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)

code <- quote({
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i]] ~ dmnorm(mu[n1[i]:n2[i]], pr[n1[i]:n2[i],n1[i]:n2[i]])
})
model <- nimbleModel(code, data=list(y=array(rnorm(7*3*2),c(3,2,7))), constants = list(mu=rep(0,6), pr=diag(6),n1 = c(3,1,2), n2 = c(6,3,3)))
expect_identical(sort(model$getNodes('y',nodesAsChars=TRUE)),
sort(c(paste0("y[1, 1, ", 3:6, "]"),
paste0("y[1, 2, ", 3:6, "]"),
paste0("y[2, 1, ", 1:3, "]"),
paste0("y[2, 2, ", 1:3, "]"),
paste0("y[3, 1, ", 2:3, "]"),
paste0("y[3, 2, ", 2:3, "]"))))
truth <- dmnorm_chol(model$y[1,1,3:6],rep(0,3),diag(4), log=TRUE)+
dmnorm_chol(model$y[1,2,3:6],rep(0,3),diag(4), log=TRUE)+
dmnorm_chol(model$y[2,1,1:3],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[2,2,1:3],rep(0,3),diag(3), log=TRUE)+
dmnorm_chol(model$y[3,1,2:3],rep(0,3),diag(2), log=TRUE)+
dmnorm_chol(model$y[3,2,2:3],rep(0,3),diag(2), log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)

code <- quote({
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n1[i]:n2[i]] ~ dwish(pr[1:4,1:4], df = 10)
})
y <- array(0, c(3,2,6,6))
n1 <- c(3,1,2); n2 <- c(6,4,5)
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n1[i]:n2[i]] <- rwish_chol(1, diag(4), 10)
model <- nimbleModel(code, data=list(y=y), constants =list( pr=diag(6),n1 = n1,n2=n2))
expect_length(model$getNodes('y',nodesAsChars=TRUE), 6*4*4)
truth <- dwish_chol(model$y[1,1,3:6,3:6],diag(4), 10, log=TRUE)+
dwish_chol(model$y[1,2,3:6,3:6],diag(4), 10, log=TRUE)+
dwish_chol(model$y[2,1,1:4,1:4],diag(4), 10, log=TRUE)+
dwish_chol(model$y[2,2,1:4,1:4],diag(4), 10, log=TRUE)+
dwish_chol(model$y[3,1,2:5,2:5],diag(4), 10, log=TRUE)+
dwish_chol(model$y[3,2,2:5,2:5],diag(4),10, log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)

code <- quote({
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n1[i]:n2[i]] ~ dwish(pr[n1[i]:n2[i],n1[i]:n2[i]], df = 10)
})
y <- array(0, c(3,2,6,6))
n1 <- c(3,1,2); n2 <- c(6,4,5)
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n1[i]:n2[i]] <- rwish_chol(1, diag(4), 10, 10)
model <- nimbleModel(code, data=list(y=y), constants =list( pr=diag(6),n1 = n1,n2=n2))
expect_length(model$getNodes('y',nodesAsChars=TRUE), 6*4*4)
truth <- dwish_chol(model$y[1,1,3:6,3:6],diag(4), 10, log=TRUE)+
dwish_chol(model$y[1,2,3:6,3:6],diag(4), 10, log=TRUE)+
dwish_chol(model$y[2,1,1:4,1:4],diag(4), 10, log=TRUE)+
dwish_chol(model$y[2,2,1:4,1:4],diag(4), 10, log=TRUE)+
dwish_chol(model$y[3,1,2:5,2:5],diag(4), 10, log=TRUE)+
dwish_chol(model$y[3,2,2:5,2:5],diag(4),10, log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)

code <- quote({
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n3[j]:n4[j]] ~ dwish(pr[n1[i]:n2[i],n3[j]:n4[j]], df = 10)
})
y <- array(0, c(3,2,6,7))
n1 <- c(3,1,2); n2 <- c(6,4,5)
n3 <- c(2,4); n4 <- c(5,7)
for(i in 1:3)
for(j in 1:2)
y[i,j,n1[i]:n2[i],n3[j]:n4[j]] <- rwish_chol(1, diag(4), 10, 10)
model <- nimbleModel(code, data=list(y=y), constants =list( pr=diag(7),n1 = n1,n2=n2,n3=n3,n4=n4))
expect_length(model$getNodes('y',nodesAsChars=TRUE), 6*4*4)

## SS case
code <- quote({
for(i in 2:4)
y[i, n1[i]:n2[i]] ~ dmnorm(y[i-1,n1[i]:n2[i]], pr[1:3,1:3])
})
model <- nimbleModel(code, data=list(y=matrix(rnorm(4*5),4)), constants = list(pr=diag(3), n1 = c(1,3,1,2), n2 = c(1,5,3,4)))
expect_identical(sort(model$getNodes('y',nodesAsChars=TRUE)),
sort(c(paste0("y[2, ", 3:5, "]"),
paste0("y[3, ", 1:3, "]"),
paste0("y[4, ", 2:4, "]"))))

truth <- dmnorm_chol(model$y[2,3:5],model$y[1,3:5],diag(3), log=TRUE)+
dmnorm_chol(model$y[3,1:3],model$y[2,1:3],diag(3), log=TRUE)+
dmnorm_chol(model$y[4,2:4],model$y[3,2:4],diag(3), log=TRUE)
model$calculate()
expect_identical(model$getLogProb('y'), truth)
})




Expand Down
6 changes: 3 additions & 3 deletions nimbleModel/tests/testthat/test-nodeRules.R
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1358,17 +1358,17 @@ test_that("nodeRange::print works correctly", {
LHSrule <- nodeRuleClass$new(LHS, 1, context_i)

expect_identical(LHSrule$apply()$toChar(),
"`mu[1:3, idx1]`, for idx1 in 2:8")
"`mu[1:3, idx1]`, for `idx1` in 2:8")

LHS <- quote(mu[1:3,j,i])
LHSrule <- nodeRuleClass$new(LHS, 1, context_ij)
expect_identical(LHSrule$apply()$toChar(),
"`mu[1:3, idx1, idx2]`, for idx1 in 3:5, idx2 in 2:8")
"`mu[1:3, idx1, idx2]`, for `idx1` in 3:5, `idx2` in 2:8")

LHS <- quote(mu[1:3,j,i-1])
LHSrule <- nodeRuleClass$new(LHS, 1, context_ij)
expect_identical(LHSrule$apply()$toChar(),
"`mu[1:3, idx1, idx2]`, for idx1 in 3:5, idx2 in 1:7")
"`mu[1:3, idx1, idx2]`, for `idx1` in 3:5, `idx2` in 1:7")

})

Expand Down
Loading