Learn how easy it is to sync an existing GitHub or Google Code repo to a SourceForge project! See Demo

Close

Debugging a faulty attempt of ones-tricking

Help
2013-01-29
2013-03-17
  • Hi all,

    I’m currently doing some statistical football modelling in Jags, and I have come across a problem with my implementation of an arbitrary distribution using the ones-trick (I have also tried the zeros-trick, but with the same result). I keep getting the error “invalid parent values” and I have come out short trying to debug it. Could anyone point me in the right direction by either helping me debug my code or pointing out the errors of my doing?

    My implementation in short:

    # Ones trick for implementing the [Dixon&Coles, 1997] mixture model
    
    C <- 1000000
    for(i in 1:noGames) {
            ones[i] <- 1
    }
    rho <- 0.1
    
    delta[i]   <- (
             attack[team[i, 1], timeslice[i, 1]] + defense[team[i, 1], timeslice[i, 1]] -
             attack[team[i, 2], timeslice[i, 2]] - defense[team[i, 2], timeslice[i, 2]]
           ) / 2
    
      log(homeLambda[i])  <-  (
             homeGoalAvg + 
             (
              attack[team[i, 1], timeslice[i, 1]] - 
              defense[team[i, 2], timeslice[i, 2]] -
              gamma * delta[i]
             )
            )
    
      log(awayLambda[i])  <-  (
             awayGoalAvg + 
             (
              attack[team[i, 2], timeslice[i, 2]] - 
              defense[team[i, 1], timeslice[i, 1]] +
              gamma * delta[i]
             )
            )
    
      is0X[i] <- goalsScored[i, 1] == 0
         isX0[i] <- goalsScored[i, 2] == 0
         is1X[i] <- goalsScored[i, 1] == 1
         isX1[i] <- goalsScored[i, 2] == 1
         is00[i] <- is0X[i] * isX0[i]
         is01[i] <- is0X[i] * isX1[i]
         is10[i] <- is1X[i] * isX0[i]
         is11[i] <- is1X[i] * isX1[i]
    
      kappa[i] <- (
          is00[i] * ( 1 + (homeLambda[i] * awayLambda[i] * rho) ) + 
          is01[i] * ( 1 - (homeLambda[i] * rho    ) ) + 
          is10[i] * ( 1 - (awayLambda[i] * rho               ) ) + 
          is11[i] * ( 1 + rho                           ) + 
          1 -    ( is00[i] + is01[i] + is10[i] + is11[i]   )
         )
    
      L[i] <- kappa[i]/C
      ones[i] ~ dbern(L[i])
    
      goalsScored[i, 1] ~ dpois( homeLambda[i] )
      goalsScored[i, 2] ~ dpois( awayLambda[i] )
    

    Any help is greatly appreciated! Thanks.

     
    Last edit: Martyn Plummer 2013-03-11
  • No-one with any experience in getting JAGS to output debug values? JAGS responds with error in parent node ones[x], where x is an arbitrary array index. I need to get a hold of the temporary values of the other vars to be able to debug this problem.

     
  • Luc Coffeng
    Luc Coffeng
    2013-03-07

    It would be much easier to comment on your question if you provide more information on what you want to do. Myself, I'm not familiar with the 'one-trick' or 'zero-trick' you refer to, and I don't know what you are trying to model. More information is welcome.

     
  • Martyn Plummer
    Martyn Plummer
    2013-03-11

    Are you sure that kappa[i] is always positive?

     
  • Thank you both for replying!

    In my test runs kappa[i] is in the interval of <0.4, 1.8> and the code runs perfectly whenever I remove the Bernoulli distribution part. This is really a mind bender for me, and it is a really important part of my masters thesis which is due in just a few months time. I need to come up with a solution fast, so any help is gratefully appreciated.

    Luc: ones- and zero-trick are hacks/workarounds for defining custom distributions in JAGS. What I am trying to do with this is to increase the probability of 0-0 and 1-1 over 0-1 and 1-0 as outcomes of soccer matches. The really short story: Goals scored and conceded by each team in a game is Poisson distributed with mean homeLambda and awayLambda respectively. These lambdas are calculated based on each team's attack and defense strengths. In addition these strengths are modelled by a Brownian Motion (Wiener Process) to evolve randomly over time. (I have also tried to remove this part, but made no real progress towards solving my issue).

     
    Last edit: Thomas R. Andresen 2013-03-12
  • Martyn Plummer
    Martyn Plummer
    2013-03-12

    Try running the model with rho set to zero. After a suitable burn-in, save the parameter values (see below) and use these as initial values for a second run with rho set to a small positive value: you might want to try smaller values of rho before working up to 0.1.

    To save parameters, if you are using the command line interface, use "parameters to myfile.R" (see the user manual). If you are using the rjags package for R, then model$state() returns the current parameter values.

     
  • Thanks Martyn! After setting rho to zero and working my way up to a small positive value (0.001) i managed to get the model running. This leads me to think that I have a error elsewhere in my model (it should work for rho=0.1). I'm grateful for your effort in helping me out.