Problem with indexing inside nimble model

[Athul Sudheesh]

I am trying to implement a bayesian network model and I am sharing my code below. I have highlighted the lines where I am having issues. 

logic for the highlighted section of the code:

Q is a binary matrix (nrow = J, ncol = K) and you can think of Q-Matrix as something similar to an adjacency matrix. Now I would like to use rows from Q matrix  (let's say 1st row) to select elements from another matrix, A (nrow=N, ncol=K) i.e A[1,Q[1,]

RDiBello <- nimbleCode({

# This below code makes sure that this code works for

# both context-BN and non-context BN.

if(is.null(C)){

C = matrix(1,N,K) 

}

# Prior Specification

# Prior for student specific guessing parameter

for (i in 1:N){

for(j in 1:J){

g[i,j] ~ dbeta(mean = 0.5, sd = 0.25)

}

}

# Prior for item specific difficulty parameter

for (j in 1:J){

d[j] ~ dbeta(mean = 0.5, sd = 0.25)

}

# Prior for student's abilities

levels[1:2] <- rep(1,2)

lambda[1:2] ~ ddirch(levels[1:2])

for (i in 1:N){

for (k in 1:K){

A[i,k] ~ dcat(lambda[1:2])

}

}

for (i in 1:N){

for (j in 1:J){

# Q-Matrix specifies which skills are required by

# each item to answer them correctly

q_j[1:K] <- nimbleLogical(Q[j,1:K])

eff_ability[i,j] <- min(A[i,q_j[]])

# d is the item specific difficulty parameter

# g is the student specific guess parameter

# effective ability is the combined effect of skills

theta[i,j] <- eff_ability[i,j] + g[i,j] - d[j]

# z[i,j] is the combined effect of the contextual factors,

# effective ability and the difficulty of the item

z[i,j] <- min(c(theta[i,j], C[i,j]))

# taking the logistic inverse

pi[i,j] <- 1/(1+exp(-(1.7* (z[i,j]))))

# Likelihood

X[i,j] ~ dbern(pi[i,j])

}

}

})

# Nimble flavor of R base's as.logical function

nimbleLogical <- nimbleRcall(function(x = double(1)){}, Rfun = 'as.logical',

returnType = logical(1))

modelConstants <- list(

N = dim(X)[1],

J = dim(X)[2],

K = dim(Q)[2]

)

data <- list(

X = X,

Q = Q,

C = bored

)

initial_values <- NULL

contextModel <- nimbleModel(

code = RDiBello,

data = data,

inits = initial_values,

dimensions = list(q_j = dim(Q)[2]),

constants = modelConstants,

name = "Context"

)

The error I am getting is: 

Error in getSymbolicParentNodesRecurse(x, constNames, indexNames, nimbleFunctionNames, : getSymbolicParentNodesRecurse: only scalar indices are allowed; vector indexing found in A[i, q_j[1:8]] Traceback:

Thanks in advance. 

[Perry de Valpine]

Hi Athul,

I just replied to your previous question about writing covertToLogical as a nimbleFunction. This question looks related, maybe showing how you are trying to use a similar idea in a model.

A key difference between model code and nimbleFunction code is that the model language is declarative (order of statements doesn't matter) while nimbleFunctions are imperative (i.e. a sequence of steps to be done in the order written, what most people are used to in most programming). This means that model code is less flexible than nimbleFunction (or R) code because from model code nimble must construct a fixed set of graphical relationships. All variables in models are double-precision numeric.

One question is whether you are trying to re-use q_j[1:K] within the for-loops. If so, that won't work because of the declarative nature of model code. You can't re-declare q_j[1:K] when iterating over values of j as you would in imperative code.

 Here is code that I think comes closest to what I see written (but I am not sure if it is what you intended):

q_j <- matrix(c(TRUE, TRUE, FALSE, TRUE, FALSE, TRUE), byrow = TRUE, nrow = 2)
initA <- matrix(rnorm(6), nrow=2)

mymin <- nimbleFunction(
  run = function(x = double(1), ind = double(1)) {
    n <- length(x)
    current_min <- Inf
    for(i in 1:n) {
      if(ind[i]) {
        if(x[i] < current_min) current_min <- x[i]
      }
    }
    return(current_min)
    returnType(double())
  }
)

m <- nimbleModel(
  nimbleCode({

    for (i in 1:N){
      for (j in 1:J){

        eff_ability[i,j] <- mymin(A[i, ], q_j[i, ])
      }
    }
  }),
  inits = list(A = initA),
  constants = list(q_j = q_j,
                   N = 2, J = 3)
  )

cm <- compileNimble(m)
cm$calculate()
cm$eff_ability

Note that eff_ability[i, j] will be the same for all i, so I am not sure if this is what you intended, but hopefully it will give you something to work from. I'm not sure if you need only eff_ability[i] or if you need the indices (over which to minimize elements of A) to depend on j.

You could also do it with a ragged array of numeric indices, where NAs indicate the unneeded elements in each row, like this:

q_j_indices <- matrix(c(1, 3, NA, 2, 4, 5), byrow = TRUE, nrow = 2)
initA <- matrix(rnorm(10), nrow=2)

mymin <- nimbleFunction(
  run = function(x = double(1), ind = double(1)) {
    n <- length(ind)
    current_min <- Inf
    for(i in 1:n) {
      if(!is.na(ind[i])) {
        if(x[ind[i]] < current_min) current_min <- x[ind[i]]
      }
    }
    return(current_min)
    returnType(double())
  }
)

m <- nimbleModel(
  nimbleCode({

    for (i in 1:N){
      for (j in 1:J){

        eff_ability[i,j] <- mymin(A[i, ], q_j_indices[i, ])
      }
    }
  }),
  inits = list(A = initA),
  constants = list(q_j_indices = q_j_indices,
                   N = 2, J = 3)
  )

cm <- compileNimble(m)
cm$calculate()
cm$eff_ability

Again eff_ability is redundant across j and you will have to adjust based on what you want.

Good luck with this.

Perry

--
You received this message because you are subscribed to the Google Groups "nimble-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nimble-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/64dc2630-7b3a-4579-8625-556fe808cbb5n%40googlegroups.com.