tmp/t.R

In gmonette/causalsim: Covariance Matrix of a Causal Graph

#' causalsim: A package to explore linear causal DAGs
#' 
#' The causalsim package uses a matrix containing the 
#' coefficients and unique standard deviations of a causal
#' DAG to generate the marginal covariance matrix and to 
#' calculate the expected value of cofficients of
#' linear models applied to a population generated by
#' the causal DAG.
#' 
#' @docType package
#' @name causalsim
NULL


# putting matrix in strict lower diagonal form if it isn't already
# we'll leave  this to later
zerotolast <- function(mat) {
	abssum <- function(x) sum(abs(x))
	colabssum <- apply(mat, 2, abssum)
	# find the last 0
	lastzero <- max(which(colabssum==0))
	if(lastzero < 1) stop("not acyclic?")
	perm <- seq_len(ncol(mat))
	perm <- c(perm[-lastzero], lastzero)
	mat[perm,perm]
}
z <- matrix(0,3,3, dimnames = list(letters[1:3],letters[1:3]))
z[1,3] <- 1
z
zerotolast(z)


mperm <- function(mat, debug = FALSE) {
 # will switch one direct parent-child pair
 # How do we know to give up??
	ind <- cbind(c(row(mat)), c(col(mat)))
	ind <- cbind(ind, c(mat[ind]))
	ind <- ind[(ind[,2]>ind[,1]) & ind[,3] != 0,,drop=FALSE]
	if(debug) disp(ind)
	if(nrow(ind) > 0) {
		perm <- seq_len(nrow(mat))
	    perm[ind[1,1]] <- ind[1,2]
	    perm[ind[1,2]] <- ind[1,1]
	    mat[perm,perm]
	} else mat
}
msize <- 5
zm <- matrix(0,msize,msize, dimnames = list(letters[1:msize], letters[1:msize]))

zm[col(zm) > row(zm)] <- 1
zm
(zm <- mperm(zm, debug = T))

find_leaves(mat)