#' test 2 graph in mm-way
#'
#' @export
test2graph <- function(graph1, graph2, anchor=c("mean","median"), method=c("ad","ks","kw","ns","steel")){
##############################################################
# preliminary : data names
DNAME = paste(deparse(substitute(graph1)),"and",deparse(substitute(graph2)))
##############################################################
# Check Input & Transform into 'igraph' object if matrix
input1 = graph1
if (!check_network(input1)){
stop("* test2graph : an input graph should be either an 'igraph' object or adjacency matrix.")
}
if (is.matrix(input1)){
if (isSymmetric(input1)){
input1 = igraph::graph_from_adjacency_matrix(input1, mode="undirected")
} else {
input1 = igraph::graph_from_adjacency_matrix(input1, mode="directed")
}
}
input2 = graph2
if (!check_network(input2)){
stop("* test2graph : an input graph should be either an 'igraph' object or adjacency matrix.")
}
if (is.matrix(input2)){
if (isSymmetric(input2)){
input2 = igraph::graph_from_adjacency_matrix(input2, mode="undirected")
} else {
input2 = igraph::graph_from_adjacency_matrix(input2, mode="directed")
}
}
##############################################################
# check condition - (1) both should be connected, and (2) set the anchor selection
if ((!igraph::is_connected(input1))||(!igraph::is_connected(input2))){
stop("* test2graph : we need both input graphs to be connected.")
}
anchor = match.arg(anchor)
##############################################################
# branching to use 'testKgraph'
method = match.arg(method)
if (all(method=="ks")){
# computation : effective resistance whose sqrt is metric
metric1 = sqrt(PNAS::effective(input1))
metric2 = sqrt(PNAS::effective(input2))
# computation : two-sample kolmogorov-smirnov test statistic (DAS)
ids = return_id(metric1, metric2, anchor)
id1 = ids[1]; tgt1 = as.vector(metric1[id1,-id1])
id2 = ids[2]; tgt2 = as.vector(metric2[id2,-id2])
# compute p-value using different methods from twosamples
ksout = suppressWarnings(stats::ks.test(tgt1, tgt2, alternative = "two.sided"))
thestat = as.double(ksout$statistic)
pvalue = as.double(ksout$p.value)
names(thestat) = "Tks"
} else {
listgraph = list()
listgraph[[1]] = input1
listgraph[[2]] = input2
result = testKgraph(listgraph, method=method)
thestat = result$statistic
pvalue = result$p.value
}
##############################################################
# REPORT
hname = "Test for Equality of Two Networks"
Ha = "two networks are not equal"
res = list(statistic=thestat, p.value=pvalue, alternative = Ha, method=hname, data.name = DNAME)
class(res) = "htest"
return(res)
}
# -----------------------------------------------------------------------
#' @keywords internal
#' @noRd
return_id <- function(metric1, metric2, method){
m = nrow(metric1)
n = nrow(metric2)
if (all(method=="random")){
id1 = sample(1:m, 1)
id2 = sample(1:n, 1)
} else if (all(method=="mean")){
id1 = which.min(rowSums(metric1^2))[1]
id2 = which.min(rowSums(metric2^2))[1]
} else if (all(method=="median")){
id1 = which.min(rowSums(metric1))[1]
id2 = which.min(rowSums(metric2))[1]
}
return(c(id1,id2))
}
# # ## test with some named graphs
# # ## scenario 1. star graph of sizes
# library(igraph)
# library(DAS)
# star10c1 = make_star(10)
# star20c1 = make_star(20)
# star20c2 = make_star(20, center=2)
#
#
# # test2graph(star20c1, star20c1)
# metric1 = sqrt(effective(star20c1)); embed1 = DAS::cmds(metric1)$embed
# metric2 = sqrt(effective(star20c2)); embed2 = DAS::cmds(metric2)$embed
#
# id1 = which.min(rowSums(metric1^2))[1]; tgt1 = metric1[id1,-id1]
# id2 = which.min(rowSums(metric2^2))[1]; tgt2 = metric2[id2,-id2]
#
# par(mfrow=c(3,2))
# plot(star20c1, main="star20c1")
# plot(star20c2, main="star20c2")
# plot(embed1[,1], embed1[,2], main="mds of star20c1")
# plot(embed2[,1], embed2[,2], main="mds of star20c2")
# plot(ecdf(tgt1), main="ecdf::star20c1")
# plot(ecdf(tgt2), main="ecdf::star20c2")
# test2graph(star20c1, star20c2, anchor = "mean", method="ad")
# test2graph(star20c1, star20c2, anchor = "mean", method="cvm")
# test2graph(star20c1, star20c2, anchor = "mean", method="kuiper")
# test2graph(star20c1, star20c2, anchor = "mean", method="wass")
# test2graph(star20c1, star20c2, anchor = "mean", method="ks")
# ks.test(tgt1, tgt2)
#
#
# niter = 12345
# pvals = array(0,c(2,niter))
# for (i in 1:niter){
# star20c1 = make_star(20)
# star20c2 = make_star(20, center=2)
# pvals[1,i] = test2graph(star20c1, star20c2, method="ad")$p.value
# pvals[2,i] = test2graph(star20c1, star20c2, method="ks")$p.value
# print(paste("iteration ",i,"/",niter," complete...", sep=""))
# }
#
# allnames = c("ad","cvm","kuiper","wass","ks")
# x11()
# par(mfrow=c(1,2), pty="s")
# for (i in 1:2){
# tgtp = as.vector(pvals[i,])
# alphas = seq(from=0.001, to=0.999, length.out=100)
# errors = rep(0,100)
# for (j in 1:100){
# tgtalpha = alphas[j]
# errors[j] = sum(tgtp <= tgtalpha)/length(tgtp)
# }
# pm = paste("error by ",allnames[i],sep="")
# plot(alphas, errors, main=pm)
# }
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.