library(bninfo)
library(dplyr)
library(ROCR)
data(isachs)
data(tcell_examples)
net <- tcell_examples$net
graphviz.plot(net)

A causal network structure is fit in a unsupervised manner. No functions exist for fitting causal networks that optimize predictions for a single node or a given node.

.data <- select(filter(isachs, INT == 0), -INT)
training <- sample_n(.data, 1000)
test <- sample_n(.data, 1000)
test_labels <- ifelse(test$Mek == "3", 1, 0)
net_unsupervised <- bn.fit(net, training, method = "bayes")

The Tree Augmented Naive Bayes algorithm fits a node like a naive Bayes classifier, but adds some dependence structure. This is a step closer in the right direction.

net_supervised <- tree.bayes(training, "Mek", 
           explanatory = setdiff(names(training), "Mek")) %>%
  bn.fit(training, method = "bayes")
graphviz.plot(net_supervised, shape = "rectangle")

Comparing prediction of Mek == 3 state with an ROC curve:

u_predictions <- rep(NA, nrow(test))
for(r in 1:nrow(test)){
  u_predictions[r] <- cpquery(net_unsupervised, (Mek == "3"), 
           (PKA == test[r, "PKA"]) & (PKC == test[r, "PKC"]) & (Raf == test[r, "Raf"]))
}
s_predictions <- predict(net_supervised, test, 
                         prob = TRUE) %>%
  attr("prob") %>%
  apply(2, function(r) r[3])
prediction(u_predictions, test_labels) %>%
  performance("tpr", "fpr") %>% 
  plot(col = "blue")
prediction(s_predictions, test_labels) %>%
  performance("tpr", "fpr") %>% 
  plot(col = "red", add = TRUE)

Now, showing that interventions also have a role in prediction.

training <- sample_n(.data, 300, replace = TRUE)
test <- rbind(sample_n(filter(.data, Erk == "1"), 100),
              sample_n(filter(.data, Erk != "1"), 100))
test_labels <- ifelse(test$Erk == "3", 1, 0)
naive_net <- naive.bayes(training, "Erk", 
           explanatory = setdiff(names(training), "Erk")) %>%
  bn.fit(training, method = "bayes")
sim_data <- rbn(naive_net, n = 1000)
naive_net_int <- naive.bayes(sim_data[, c("Erk", "PKC", "PKA", "INT")], "Erk",
                             explanatory = c("PKC", "PKA", "INT")) %>%
  bn.fit(sim_data[, c("Erk", "PKC", "PKA", "INT")], method = "bayes")
naive_net_no_int <- naive.bayes(sim_data[, c("Erk", "PKC", "PKA")], "Erk",
                             explanatory = c("PKC", "PKA")) %>%
  bn.fit(sim_data[, c("Erk", "PKC", "PKA")], method = "bayes")
int_predictions <- predict(naive_net_int, test[, c("Erk", "PKC", "PKA", "INT")], 
                           prob = TRUE) %>%
  attr("prob") %>%
  apply(2, function(r) r[3])
no_int_predictions <- predict(naive_net_no_int, test[, c("Erk", "PKC", "PKA")], 
                              prob = TRUE) %>%
  attr("prob") %>%
  apply(2, function(r) r[3])
prediction(int_predictions, test_labels) %>%
  performance("tpr", "fpr") %>% 
  plot(col = "blue", main = "Prediction of Erk activity given PKC and PKA")
prediction(no_int_predictions, test_labels) %>%
  performance("tpr", "fpr") %>% 
  plot(col = "red", add = TRUE)


robertness/bninfo documentation built on May 27, 2019, 10:32 a.m.