Create conditional probability tables (CPTs) with gRain

Søren Højsgaard

Create the universe:

yn <- c("yes", "no")
uni <- list(Rain = yn, Sprinkler = yn, GrassWet = yn)

Just a utility:

flatten <- function(x){
    if (length(dim(x))==1) x else   ftable(x, row.vars=1)
} 

Using gRain::cptable

Notice that probabilities are normalized

p.R <- cptable("Rain", levels=uni, values=c(2, 8)) 
p.S_R <- cptable(c("Sprinkler","Rain"), levels = uni, values=c(.01, .99, .4, .6))
p.G_SR <- cptable(~GrassWet:Sprinkler:Rain, levels = uni, values=c(.99, .01, .8, .2, .9, .1, 0, 1)) 

cptlst <- list(p.R, p.S_R, p.G_SR)
lapply(cptlst, flatten)
## [[1]]
## Rain
## yes  no 
## 0.2 0.8 
## 
## [[2]]
##           Rain  yes   no
## Sprinkler               
## yes            0.01 0.40
## no             0.99 0.60
## 
## [[3]]
##          Sprinkler  yes        no     
##          Rain       yes   no  yes   no
## GrassWet                              
## yes                0.99 0.90 0.80 0.00
## no                 0.01 0.10 0.20 1.00

Simplification using Curry

A simpler specification can be obtained by using the Curry function from the functional package (imported in gRain): Some of the arguments are fixed and need not be specified.

cptab <- Curry(cptable, levels=uni)

p.R <- cptab("Rain", values=c(.2, .8)) 
p.S_R <- cptab(c("Sprinkler","Rain"), values=c(.01, .99, .4, .6))
p.G_SR <- cptab(~GrassWet:Sprinkler:Rain, values=c(.99, .01, .8, .2, .9, .1, 0, 1)) 

cptlst <- list(p.R, p.S_R, p.G_SR)
## lapply(cptlst, flatten)

Yet a simpler way using mapply

vpa <- list(~Rain,
            ~Sprinkler|Rain,
            ~GrassWet:Sprinkler:Rain)
cpt <- list(c(.2, .8),
            c(.01, .99, .4, .6),
            c(.99, .01, .8, .2, .9, .1, 0, 1))

cptlst <- mapply(cptab, vpa, cpt)
## lapply(cptlst, flatten)

Using gRbase::arMk - with normalization

Above, cptable creates arrays, and the same can be achieved e.g. with gRbase::arMk.

p.R <- arMk("Rain", levels=uni, values=c(2, 8),
            normalize="first") 
p.S_R <- arMk(c("Sprinkler","Rain"), levels = uni, values=c(.01, .99, .4, .6),
              normalize="first")
p.G_SR <- arMk(~GrassWet:Sprinkler:Rain, levels = uni, values=c(.99, .01, .8, .2, .9, .1, 0, 1),
               normalize="first") 
cptlst <- list(p.R, p.S_R, p.G_SR)

Again, the specification can be simplified using Curry:

arMk2 <- Curry( arMk, levels=uni, normalize="first")
cptlst <- mapply(arMk2, vpa, cpt)

Checking the list of CPTs

So far the list of CPTs is nothing but a list. To check that the CPT's specify a valid network we can check that they all together specify a DAG.

cptl <- compileCPT( cptlst )
cptl
## CPTspec with probabilities:
##  P( Rain )
##  P( Sprinkler | Rain )
##  P( GrassWet | Sprinkler Rain )

On the contrary, below we specify CPTs that do not define a DAG.

uni2 <- list(a=c("yes","no"), b=c("yes", "no"))
p1 <- cptable(~a|b, levels=uni2, values=c(7, 3, 9, 1))
p2 <- cptable(~b|a, levels=uni2, values=c(7, 3, 9, 1))