tests/regression_tests/hojsgaard_model_tests/triangle_model/triangle-model_mrf.R

In npetraco/CRFutil: Handy Utility Functions for Use with CRF and gRbase packages

library(CRFutil)

# Load samps generated by true model:
fpth <- "/home/npetraco/codes/R/CRFutil/tests/regression_tests/hojsgaard_model_tests/triangle_model/triangle_data/"
setwd(fpth)
load(paste0(fpth,"triangle_samp.RData"))
head(samps)

# Instantiate an empty model for the Bayesian fit:
# Fully connected
grphf <- ~1:2+1:3+2:3
adj <- ug(grphf, result="matrix")
# Check the graph:
gp <- ug(grphf, result = "graph")
dev.off()
iplot(gp)

f0 <- function(y){ as.numeric(c((y==1),(y==2)))}
n.states <- 2


mrf.m <- make.crf(adj, n.states)
mrf.m <- make.features(mrf.m)
mrf.m <- make.par(mrf.m, 6)
mrf.m$node.par[1,1,] <- 1
mrf.m$node.par[2,1,] <- 2
mrf.m$node.par[3,1,] <- 3

mrf.m$edge.par[[1]][1,1,1] <- 4
mrf.m$edge.par[[1]][2,2,1] <- 4
mrf.m$edge.par[[2]][1,1,1] <- 5
mrf.m$edge.par[[2]][2,2,1] <- 5
mrf.m$edge.par[[3]][1,1,1] <- 6
mrf.m$edge.par[[3]][2,2,1] <- 6
#----------

# Fit model to samples from the known model and obtain an estimate for w:
mrf.m$par
mrf.m$node.pot
train.mrf(mrf.m, samps, nll=mrf.exact.nll, infer.method = infer.exact)
mrf.m$par.stat
mrf.m$par
# 0.1211475  0.8901192 -0.7510040 -0.4105264  0.6291559  1.0374155
mrf.m$node.pot

pot.info2 <- make.gRbase.potentials(mrf.m, node.names = c("X.1","X.2","X.3"), state.nmes = c("1","2"))
#pot.info2

gR.dist.info2    <- distribution.from.potentials(pot.info2$node.potentials, pot.info2$edge.potentials)
#logZ2            <- gR.dist.info2$logZ
joint.dist.info2 <- as.data.frame(as.table(gR.dist.info2$state.probs))
joint.dist.info2
mrf.dist.info <- joint.dist.info2[,c(2,3,1,4)]
colnames(mrf.dist.info)[4] <- c("MRF.Freq")
mrf.dist.info
save(mrf.dist.info, file = paste0(fpth,"triangle_mrf_dist.RData"))