在 JAGS/BUGS 中建模混合模型

机器算法验证 混合模式 马尔可夫链蒙特卡罗 锯齿
2022-03-25 19:03:18

我目前正在 JAGS 中实现足球结果预测模型。实际上,我已经实现了几个,但我已经遇到了迄今为​​止最困难的挑战:Rue & Salvesen 在他们的论文“Prediction and retrospective analysis of football matches in a League”中描述的一个模型。他们的模型使用混合模型来截断以 5 个进球后攻击/防守强度为条件的泊松分布。他们还采用了 Dixon & Coles (1997) 的一项法则,以增加低分比赛中 0-0 和 1-1 结果的概率。

我的问题如下,我正在尝试实现混合模型:

πg1(xA,B,yA,B|λA,B(x),λA,B(y))=κ(xA,B,yA,B|λA,B(x),λA,B(y))Po(xA,B|λA,B(x))Po(yA,B|λA,B(y))
在哪里xA,B表示主队在 A 队和 B 队之间的比赛中的进球数,并且log(λA,B(x))表示球队的实力。我试图通过使用 0-ones 技巧在 JAGS 中实现这两个定律,但到目前为止没有运气(error: illegal parent values)。到目前为止,我的 JAGS 模型:

data {
    C <- 10000

    for(i in 1:noGames) {
        zeros[i] <- 0
    }

    homeGoalAvg <- 0.395
    awayGoalAvg <- 0.098

    rho <- 0.1
}

model {

    ### Time model - Brownian motion
    tau ~ dgamma(10, 0.1)
    precision ~ dgamma(0.1, 1)

    for(t in 1:noTeams) {
        attack[t, 1] ~ dnorm(0, precision)
        defence[t, 1] ~ dnorm(0, precision)

        for(s in 2:noTimeslices) {
            attack[t, s] ~ dnorm(attack[t, (s-1)], (tau * precision) / 
                                         (abs(days[t,s]-days[t,s-1])))
            defence[t, s] ~ dnorm(defence[t, (s-1)], (tau * precision) / 
                                          (abs(days[t,s]-days[t,s-1])))
        }
    }

    ### Goal model
    gamma ~ dunif(0, 0.1)

    for(i in 1:noGames) {

        delta[i]            <-  (
                                attack[team[i, 1], timeslice[i, 1]] + 
                                defence[team[i, 1], timeslice[i, 1]] -
                                attack[team[i, 2], timeslice[i, 2]] - 
                                defence[team[i, 2], timeslice[i, 2]]
                            ) / 2

        log(homeLambda[i])  <-  (
                                    homeGoalAvg + 
                                    (
                                        attack[team[i, 1], timeslice[i, 1]] - 
                                        defence[team[i, 2], timeslice[i, 2]] -
                                        gamma * delta[i]
                                    )
                                )

        log(awayLambda[i])  <-  (
                                    awayGoalAvg + 
                                    (
                                        attack[team[i, 2], timeslice[i, 2]] - 
                                        defence[team[i, 1], timeslice[i, 1]] +
                                        gamma * delta[i]
                                    )
                                )

        goalsScored[i, 1] ~ dpois( homeLambda[i] )
        goalsScored[i, 2] ~ dpois( awayLambda[i] )

        is0X[i] <- ifelse(goalsScored[i, 1]==0, 1, 0)
        isX0[i] <- ifelse(goalsScored[i, 2]==0, 1, 0)
        is1X[i] <- ifelse(goalsScored[i, 1]==1, 1, 0)
        isX1[i] <- ifelse(goalsScored[i, 2]==1, 1, 0)
        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]     )
                    )

        # This does not work!
        zeros[i] ~ dpois(-log(kappa[i]) + C)
    }

}
2个回答

有时 ifelse 不起作用。反而

is0X[i] <- ifelse(goalsScored[i, 1]==0, 1, 0)

你应该试试

is0X[i] <- goalsScored[i, 1]==0

targetsScored[i, 1]==0 如果为真则返回 1,如果为假则返回 0

我不认为你可以zeros[i] ~ dpois(-log(kappa[i]) + C)在模型构造内部定义。

尝试将代码修改为zeros ~ dpois(-log(kappa[i]) + C)(取出'[i]')。

定义模型后,重新定义 zeors 处的数据:

data$zero=0

试试这是否有效。

有关更多信息,请参阅JAGS 的过零技巧:随机查找根