# NOT RUN {
## Simulate the true DAG
set.seed(123)
p <- 7
myDAG <- randomDAG(p, prob = 0.2) ## true DAG
myCPDAG <- dag2cpdag(myDAG) ## true CPDAG
myPDAG <- addBgKnowledge(myCPDAG,2,3) ## true PDAG with background knowledge 2 -> 3
covTrue <- trueCov(myDAG) ## true covariance matrix
## simulate Gaussian data from the true DAG
n <- 10000
dat <- rmvDAG(n, myDAG)
## estimate CPDAG and PDAG -- see help(pc)
suffStat <- list(C = cor(dat), n = n)
pc.fit <- pc(suffStat, indepTest = gaussCItest, p=p, alpha = 0.01)
pc.fit.pdag <- addBgKnowledge(pc.fit@graph,2,3)
if (require(Rgraphviz)) {
## plot the true and estimated graphs
par(mfrow = c(1,3))
plot(myDAG, main = "True DAG")
plot(pc.fit, main = "Estimated CPDAG")
plot(pc.fit.pdag, main = "Estimated CPDAG")
}
## Supppose that we know the true CPDAG and covariance matrix
(l.ida.cpdag <- ida(2,5, covTrue, myCPDAG, method = "local", type = "cpdag"))
(g.ida.cpdag <- ida(2,5, covTrue, myCPDAG, method = "global", type = "cpdag"))
## The global and local method produce the same unique values.
stopifnot(all.equal(sort(unique(g.ida.cpdag)),
sort(unique(l.ida.cpdag))))
## Supppose that we know the true PDAG and covariance matrix
(l.ida.pdag <- ida(2,5, covTrue, myPDAG, method = "local", type = "pdag"))
(g.ida.pdag <- ida(2,5, covTrue, myPDAG, method = "global", type = "pdag"))
## The global and local method produce the same unique values.
stopifnot(all.equal(sort(unique(g.ida.pdag)),
sort(unique(l.ida.pdag))))
## From the true DAG, we can compute the true causal effect of 2 on 5
(ce.2.5 <- causalEffect(myDAG, 5, 2))
## Indeed, this value is contained in the values found by both the
## local and global method
## When working with data, we have to use the estimated CPDAG and
## the sample covariance matrix
(l.ida.cpdag <- ida(2,5, cov(dat), pc.fit@graph, method = "local", type = "cpdag"))
(g.ida.cpdag <- ida(2,5, cov(dat), pc.fit@graph, method = "global", type = "cpdag"))
## The unique values of the local and the global method are still identical,
stopifnot(all.equal(sort(unique(g.ida.cpdag)), sort(unique(l.ida.cpdag))))
## and the true causal effect is contained in both sets, up to a small
## estimation error (0.868 vs. 0.857)
stopifnot(all.equal(ce.2.5, max(l.ida.cpdag), tolerance = 0.02))
## When working with data, we have to use the estimated PDAG and
## the sample covariance matrix
(l.ida.pdag <- ida(2,5, cov(dat), pc.fit.pdag, method = "local", type = "pdag"))
(g.ida.pdag <- ida(2,5, cov(dat), pc.fit.pdag, method = "global", type = "pdag"))
## The unique values of the local and the global method are still identical,
stopifnot(all.equal(sort(unique(g.ida.pdag)), sort(unique(l.ida.pdag))))
## and the true causal effect is contained in both sets, up to a small
## estimation error (0.868 vs. 0.857)
stopifnot(all.equal(ce.2.5, max(l.ida.pdag), tolerance = 0.02))
# }
Run the code above in your browser using DataCamp Workspace