Error: bad parameters for self coded distribution

59 views
Skip to first unread message

Annie Shapiro

unread,
Oct 24, 2023, 11:40:14 AM10/24/23
to nimble-users
Hi, 

I was wondering if you'd be able to explain what is causing the following error I get when making my model. Below is the error and my code. I've also attached data. 

Defining model Error: bad parameters for distribution dsmooth(k = k, L = L). (No available re-parameterization found.)

f_vec <- nimbleFunction(
   run = function(x = double(0)){

     k <- 6.06 # shape
     r <- 1/.34 # scale = 1/rate
   
     tmp_vec <- rep(NA, x)
     
     for (i in 1:x) {
       tmp_vec[i] <- pgamma(x-(i-1), k, rate = r) - pgamma(x-i, k, rate = r)
     }
     
     return(tmp_vec)
     returnType(double(1))
})

# code my own choose function
my_choose <- nimbleFunction(
  run = function(n = double(0), k = double(0)){
   
    out <- factorial(n)/(factorial(k)*factorial(n-k))
    return(out)
    returnType(double(0))
  })

# code smoother as prior
dsmooth <- nimbleFunction(
  run = function(x = double(1),
                 k = double(0),
                 L = double(0),
                 log = integer(0, default = 0)) {
   
    returnType(double(1))
    tmp_mat <- matrix(NA, nrow = k+1, ncol = L)
   
    for (j in 1:L) {
    for(m in 0:k){
      if ((j+m-k/2) > L){
        tmp_mat[(m+1), j] <- my_choose(k, m)/(2^k)*x[L] }
      else if ((j+m-k/2) < 1){
        tmp_mat[(m+1), j] <- 0 }
      else {
        tmp_mat[(m+1), j] <- my_choose(k, m)/(2^k)*x[j+m-k/2]}
    }}
   
    tmp_vec <- rep(NA, L)
   
    for (c in 1:L) {
      tmp_vec[c] <- sum(tmp_mat[, c])
    }
   
    logProb <- log(tmp_vec)
    if(log) return(logProb)
    else return(exp(logProb))
  })


#====================================================================
# set up MCMC
#====================================================================

modelCode <- nimbleCode({
 
  for (i in 1:L){
    #get mu
    mu[i] <- sum(f_vec(i)*lambda[1:i])
   
    #likelihood
    x[i] ~ dpois(mu[i])
   
    }
 
  #prior
  lambda[1:L] ~ dsmooth(k = k, L = L)
 
   })



modelConsts <- list(L = nrow(dat), k = 4)

modelData <- list(x = dat$Notifications)

initsList <- list(lambda = abs(rnorm(nrow(dat), mean(dat$Notifications), sd(dat$Notifications))))


model <- nimbleModel(modelCode, data = modelData, inits = initsList, constants = modelConsts)
question_data.xlsx

Benjamin R Goldstein

unread,
Oct 24, 2023, 11:51:31 AM10/24/23
to Annie Shapiro, nimble-users
Hi—

I think you want returnType(double(0)) in your function rather than double(1). The former is saying you’ll return a scalar (0=0 dimensions) where the latter is saying you’ll return a vector. 

Haven’t tested this myself though. 

Hope this helps—
Ben

--
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/06da08a5-86cb-4180-9b40-6e18d61195e6n%40googlegroups.com.

Shapiro, Annie

unread,
Oct 24, 2023, 12:12:49 PM10/24/23
to Benjamin R Goldstein, nimble-users
Hi Ben, 

Thanks for your quick response. I think I want my function to return a vector -- I'm thinking of it as a joint distribution, and I believe th. at's how I've coded it (could be wrong, though). Also, I get the same error when I change to double(0). 

Best,
Annie

Wei Zhang

unread,
Oct 24, 2023, 2:46:17 PM10/24/23
to Shapiro, Annie, Benjamin R Goldstein, nimble-users
Hi Annie,

A quick thought: for probability distributions no matter what the dimension of your variable is, the prob density or mass function should return a scalar, i.e. double(0). In your code, tmp_vec is a vector so logProb is a vector, it sounds like you need something like sum(logProb) and then return it. I may be misunderstanding your distribution. Hope it helps. 

Best wishes,
Wei

Chris Paciorek

unread,
Oct 24, 2023, 8:36:27 PM10/24/23
to Wei Zhang, Shapiro, Annie, Benjamin R Goldstein, nimble-users
To follow up on Wei's point, if you wrote the `rsmooth` function to generate from your "smooth" distribution, it would return a vector, but the density function returns the joint density (which is a scalar) of the vector-valued random variable.

Reply all
Reply to author
Forward
0 new messages