## Programs to apply net coin analysis
# Image is a files vector with length and order equal to nrow(nodes). Place as nodes field
# Batch
netCoin <- function(nodes = NULL, links = NULL, tree = NULL,
community = NULL, layout = NULL,
name = NULL, label = NULL, group = NULL, groupText = FALSE,
labelSize = NULL, size = NULL, color = NULL, shape = NULL,
border = NULL, legend = NULL, sort = NULL, decreasing = FALSE,
ntext = NULL, info = NULL, image = NULL, imageNames = NULL,
centrality = NULL,
nodeBipolar = FALSE, nodeFilter = NULL, degreeFilter = NULL,
lwidth = NULL, lweight = NULL, lcolor = NULL, ltext = NULL,
intensity = NULL, linkBipolar = FALSE, linkFilter = NULL,
repulsion = 25, distance = 10, zoom = 1,
fixed = showCoordinates, limits = NULL,
main = NULL, note = NULL, showCoordinates = FALSE, showArrows = FALSE,
showLegend = TRUE, frequencies = FALSE, showAxes = FALSE,
axesLabels = NULL, scenarios = NULL, help = NULL, helpOn = FALSE,
mode = c("network","heatmap"), roundedItems = FALSE, controls = 1:8,
cex = 1, background = NULL, defaultColor = "#1f77b4",
language = c("en","es","ca"), dir = NULL)
{
if(is.null(links) && is.null(nodes)){
stop("You must explicit a nodes or links data frame.")
}
if(inherits(nodes, 'netCoin')){
stop("Using netCoin function to change netCoin object attributes is deprecated. Use addNetCoin instead.")
}
name <- nameByLanguage(name,language,nodes)
if(!is.null(nodes)){
if (all(inherits(nodes,c("tbl_df","tbl","data.frame"),TRUE))) nodes<-as.data.frame(nodes) # convert haven objects
}
if(!is.null(links)){
if (all(inherits(links,c("tbl_df","tbl","data.frame"),TRUE))) links<-as.data.frame(links) # convert haven objects
}
color <- setAttrByValueKey("color",color,nodes)
shape <- setAttrByValueKey("shape",shape,nodes)
lcolor <- setAttrByValueKey("lcolor",lcolor,links)
net <- network_rd3(nodes = nodes, links = links, tree = tree,
community = community, layout = layout,
name = name, label = label, group = group, groupText = groupText,
labelSize = labelSize, size = size, color = color, shape = shape,
border = border, legend = legend,
sort = sort, decreasing = decreasing, ntext = ntext, info = info,
image = image, imageNames = imageNames,
nodeBipolar = nodeBipolar, nodeFilter = nodeFilter, degreeFilter = degreeFilter,
source = "Source", target = "Target",
lwidth = lwidth, lweight = lweight, lcolor = lcolor, ltext = ltext,
intensity = intensity, linkBipolar = linkBipolar, linkFilter = linkFilter,
repulsion = repulsion, distance = distance, zoom = zoom,
fixed = fixed, limits = limits,
main = main, note = note, showCoordinates = showCoordinates, showArrows = showArrows,
showLegend = showLegend, frequencies = frequencies, showAxes = showAxes,
axesLabels = axesLabels, scenarios = scenarios, help = help, helpOn = helpOn,
mode = mode, roundedItems = roundedItems, controls = controls, cex = cex,
background = background, defaultColor = defaultColor,
language = language, dir = dir)
class(net) <- c("netCoin",class(net))
if(!is.null(centrality)){
columns <- calCentr(net, centrality)$nodes
for(col in setdiff(colnames(columns),c("nodes","degree"))){
net$nodes[[col]] <- columns[[col]]
}
}
return(net)
}
setAttrByValueKey <- function(name,item,items){
if(is.list(item) && !is.data.frame(item)){
checkedlist <- list()
for(k in names(item)){
if(!k %in% colnames(items) || !(is.character(items[[k]]) || is.factor(items[[k]]))){
warning(paste0(name,": the names in the list must match character columns of the items, but '",k,"' doesn't"))
}else{
if(!is.character(item[[k]]) || is.null(names(item[[k]]))){
warning(paste0(name,": each item in the list must be a named character vector describing value-",name,", but '",k,"' doesn't"))
}else{
checkedlist[[k]] <- unname(item[[k]][items[[k]]])
}
}
}
if(length(checkedlist)){
item <- as.data.frame(checkedlist)
}else{
item <- NULL
}
}
return(item)
}
addNetCoin <- function(x, ...){
arguments <- list(...)
for(n in c("nodes","links","tree")){
if(!(n %in% names(arguments))){
arguments[[n]] <- x[[n]]
}
}
options <- x$options
getOpt <- function(opt,item=opt){
if(item %in% names(arguments)){
return(arguments[[item]])
}else{
if(!is.null(options[[opt]])){
return(options[[opt]])
}else{
return(NULL)
}
}
}
attributes <- c(
"name" = "nodeName",
"cex" = "cex",
"distance" = "distance",
"repulsion" = "repulsion",
"zoom" = "zoom",
"scenarios" = "scenarios",
"limits" = "limits",
"main" = "main",
"note" = "note",
"help" = "help",
"background" = "background",
"language" = "language",
"nodeBipolar" = "nodeBipolar",
"linkBipolar" = "linkBipolar",
"helpOn" = "helpOn",
"frequencies" = "frequencies",
"defaultColor" = "defaultColor",
"controls" = "controls",
"mode" = "mode",
"axesLabels" = "axesLabels",
"fixed" = "fixed",
"showCoordinates" = "showCoordinates",
"showArrows" = "showArrows",
"showLegend" = "showLegend",
"showAxes" = "showAxes",
"roundedItems" = "roundedItems",
"label" = "nodeLabel",
"labelSize" = "nodeLabelSize",
"group" = "nodeGroup",
"groupText" = "groupText",
"size" = "nodeSize",
"color" = "nodeColor",
"shape" = "nodeShape",
"border" = "nodeBorder",
"legend" = "nodeLegend",
"ntext" = "nodeText",
"info" = "nodeInfo",
"sort" = "nodeOrder",
"decreasing" = "decreasing",
"image" = "imageItems",
"imageNames" = "imageNames",
"lwidth" = "linkWidth",
"lweight" = "linkWeight",
"lcolor" = "linkColor",
"ltext" = "linkText",
"intensity" = "linkIntensity"
)
for(item in names(attributes)){
arguments[[item]] <- getOpt(attributes[[item]],item)
}
return(do.call(netCoin,arguments))
}
# Program to apply nets to correlations
netCorr<-function(variables, weight=NULL, pairwise=FALSE,
minimum=-Inf, maximum=Inf, sort=FALSE, decreasing=TRUE,
frequency=FALSE, means=TRUE,
method=c("pearson", "kendall", "spearman"), criteria="p", Bonferroni=FALSE,
minL=0, maxL=Inf,
sortL=NULL, decreasingL=TRUE,
igraph=FALSE, ...)
{
arguments <- list(...)
arguments$name <- nameByLanguage(arguments$name,arguments$language,arguments$nodes)
method <- r.method(method[1])
if(!("size" %in% names(arguments))) arguments$size <- "mean"
if(!("lwidth" %in% names(arguments))) arguments$lwidth <- "value"
if(!("lweight" %in% names(arguments))) arguments$lweight <- "value"
if(!pairwise) variables<-na.omit(variables)
if (inherits(weight,"character")) variables <- variables[, setdiff(names(variables), weight)]
cases<-nrow(variables)
if(exists("linkBipolar", arguments) && arguments$linkBipolar) xl <- 2
else xl <-1
if (criteria=="p" & maxL==Inf) maxL<-.5*xl
if (criteria=="p" & Bonferroni) maxL<-maxL/choose(cases,2)
if (criteria %in% c("pearson", "kendall", "spearman")) criteria <- "value"
statistics <-data.frame(name=colnames(variables),
mean=round(apply(variables,2,mean, na.rm=TRUE),2),
std=round(sqrt(apply(variables,2,var, na.rm=TRUE)),2),
min=apply(variables,2,min, na.rm=TRUE),
max=apply(variables,2,max, na.rm=TRUE))
colnames(statistics)[1] <- arguments$name
if(!is.null(arguments$nodes)) arguments$nodes <- merge(statistics, arguments$nodes, by=arguments$name, all.x=TRUE, sort=FALSE)
else arguments$nodes <- statistics
if (pairwise) use <- "pairwise.complete.obs"
else use <- "complete.obs"
R<-cor(variables[,arguments$nodes[,2]>=minimum & arguments$nodes[,2]<=maximum],method=method[1], use=use)
E<-edgeList(R, "shape", min=-1, max=1, directed=FALSE, diagonal=FALSE)
E$z<-E$value/sqrt(1-E$value^2)*sqrt(cases-2)
if(exists("linkBipolar", arguments) && arguments$linkBipolar) E$p <- (1-pt(abs(E$z),40000-2))*2
else xl <- E$p<-1-pt(E$z,cases-2)
E<-E[E[[criteria]]>=minL & E[[criteria]]<=maxL,]
if (!is.null(sortL)) E<-E[order((-1*decreasingL+!decreasingL)*E[[sortL]]),]
arguments$links <- E
if(exists("layout", arguments) && is.character(arguments$layout) && tolower(substr(arguments$layout,1,2))=="pc") arguments$layout <- layoutPCA(R)
xNx <- do.call(netCoin,arguments)
if (igraph) return(toIgraph(xNx))
else return(xNx)
}
# Program to apply evolving nets to correlations.
d_netCorr <- function(variables, nodes= NULL, weight=NULL,
pairwise=FALSE, minimum=-Inf, maximum=Inf,
frequency=FALSE, means=TRUE,
method=c("pearson", "kendall", "spearman"), criteria="value", Bonferroni=FALSE,
minL=0, maxL=Inf,
sortL=NULL, decreasingL=TRUE,
factorial=c("null", "pc", "nf", "vf", "of"),
components=TRUE, backcomponents=FALSE, sequence=seq(.20, 1, .01),
textFilter=c(1, .99), speed=50, dir=NULL, ...)
{
arguments <- list(variables= variables, nodes=nodes, weight=weight,
pairwise=pairwise, minimum=minimum, maximum=maximum,
frequency=frequency, means=means,
method=method, criteria=criteria, Bonferroni=Bonferroni,
minL=minL, maxL=maxL, sortL=sortL, decreasingL=decreasingL, ...)
if(exists("layout", arguments) & !identical(factorial, c("null", "pc","nf","rf","of"))) warning("Argument factorial is incompatible with layout")
if(exists("layout", arguments)) C <- arguments$layout else C <- layoutFact(variables, factorial)
G <- do.call (netCorr, arguments)
g <- list()
g$lineplots <- c("component", "Degree", "closeness", "betweenness", "eigen", "ratio" )
g$mode <- "frame"
g$speed <- speed
if(!is.null(dir)){
g$dir <- dir
}
ch=0
if(length(textFilter)<2) textFilter[2]<- .99
for(I in sequence) {
G$links$text <- ifelse(abs(G$links$value)>=min(textFilter[1], I+textFilter[2]), sprintf("%.2f", G$links$value), "")
H <- addNetCoin(G, linkFilter = paste0("value>",I), lwidth="value",
ltext="text", size="mean", layout=C,
main=paste0("Correlation: ", I))
if(exists("hidden", H$links)) H$links <- H$links[!H$links$hidden, ]
else H$links$hidden <- FALSE
enodes <- unique(c(H$links$Source,H$links$Target))
H$nodes <- H$nodes[H$nodes[[H$options$nodeName]] %in% enodes,]
comps <- igraph::components(toIgraph(H))
cc <- sum(comps$csize>1)
if(((cc > ch | !components) | (cc < ch & components & backcomponents)) & nrow(H$links)>0) {
central <- calCentr(H)
H$nodes <- cbind(H$nodes, central$nodes)
H$nodes <- cbind(H$nodes, compon(comps$membership))
names(H$nodes)[ncol(H$nodes)] <- "component"
names(H$nodes)[names(H$nodes)=="degree"] <- "Degree"
H$links$ratio <- H$links$value/I
g[[paste0("C",I)]] <- H
}
ch <- cc
}
multi <- do.call(multigraphCreate, g)
return(multi)
}
# Complete netCoin from an incidences matrix
allNet<-function(incidences, weight = NULL, subsample = FALSE, pairwise = FALSE,
minimum=1, maximum = nrow(incidences), sort = FALSE, decreasing = TRUE,
frequency = FALSE, percentages = TRUE,
procedures = "Haberman", criteria = "Z", Bonferroni = FALSE,
support = -Inf, minL = -Inf, maxL = Inf,
directed = FALSE, diagonal = FALSE, sortL = NULL, decreasingL = TRUE,
igraph = FALSE, dir=NULL, ...)
{
arguments <- list(...)
arguments$dir<-dir
if((criteria=="Z" | criteria=="hyp") & maxL==Inf) maxL=.5
if(!("language" %in% names(arguments))) arguments$language <- "en"
arguments$name <- nameByLanguage(arguments$name,arguments$language,arguments$nodes)
if (!("size" %in% names(arguments)))
if(percentages)
arguments$size <- "%"
if (!("level" %in% names(arguments))) level<-.95 else level <-arguments$level
if (!pairwise) incidences<-na.omit(incidences)
if (inherits(weight, "character")) incidences <- incidences[, setdiff(names(incidences), weight)]
incidences <- incidences[,colSums(incidences)>0]
if (all(is.na(incidences) | incidences==0 | incidences==1)) {
C<-coin(incidences, minimum, maximum, sort, decreasing, weight=weight, subsample=subsample, pairwise = pairwise)
if(exists("size",arguments))if(arguments$size=="frequency")frequency=TRUE
O<-asNodes(C,frequency,percentages,arguments$language)
names(O)[1]<-arguments$name
if (is.null(arguments$nodes)){
if(any(sapply(incidences,function(X) {"label" %in% names(attributes(X))}))) {
label <- "label"
if(arguments$language %in% c("es","ca")){
label <- "etiqueta"
}
O[[label]] <- "NULL"
O[[label]] <- ifelse(sapply(incidences, attr, "label")=="NULL", O[[arguments$name]], sapply(incidences, attr, "label"))
arguments$label <- label
}
arguments$nodes<-O
}else{
nodesOrder<-as.character(arguments$nodes[[arguments$name]])
arguments$nodes<-merge(O,arguments$nodes[,setdiff(names(arguments$nodes),frequencyList),drop=FALSE],by.x=arguments$name,by.y=arguments$name,all.y=TRUE, sort=FALSE)
row.names(arguments$nodes)<-arguments$nodes[[arguments$name]]
arguments$nodes<-arguments$nodes[nodesOrder,]
}
procedures<-union(procedures,unlist(arguments[c("lwidth","lweight","lcolor","ltext")]))
arguments$links<-edgeList(C, procedures, criteria, level, Bonferroni, minL, maxL, support,
directed, diagonal, sortL, decreasingL)
for(lattr in c("lwidth","lweight","lcolor","ltext"))
if(!is.null(arguments[[lattr]])) arguments[[lattr]]<-i.method(c.method(arguments[[lattr]]))
if(is.character(arguments$layout)){
if(tolower(substr(arguments$layout,1,2))=="mc")arguments$layout<-layoutMCA(incidences)
else if(tolower(substr(arguments$layout,1,2))=="pc")arguments$layout<-layoutPCA(C)
}
arguments$scenarios <- attr(C,"n")
xNx <- do.call(netCoin,arguments)
if (igraph) return(toIgraph(xNx))
else return(xNx)
}
else warning("Input is not a dichotomous matrix of incidences")
}
# surCoin is a wrapper to build a netCoin object from an original non-dichotomized data.frame. See below dichotomize()
surCoin<-function(data,variables=names(data), commonlabel=NULL,
dichotomies=NULL, valueDicho=1, metric=NULL, exogenous=NULL,
weight=NULL, subsample=FALSE, pairwise=FALSE,
minimum=1, maximum=nrow(data), sort=FALSE, decreasing=TRUE,
frequency=FALSE, percentages=TRUE,
procedures="Haberman", criteria="Z", Bonferroni=FALSE,
support=-Inf, minL=-Inf, maxL=Inf,
directed=FALSE, diagonal=FALSE, sortL=NULL, decreasingL=TRUE,
igraph=FALSE, coin=FALSE, dir=NULL, ...)
{
arguments <- list(...)
if((criteria=="Z" | criteria=="hyp") & maxL==Inf) maxL=.5
varOrder <- variables # To order variables later before coin
#Check methods. No necessary because edgeList call these routines.
#procedures<-i.method(c.method(procedures))
#criteria<-i.method(c.method(criteria))
procedures<-union(procedures,unlist(arguments[c("lwidth","lweight","lcolor","ltext")]))
#Metrics
if(!is.null(metric)) {
variables<-setdiff(variables,metric)
procedures<-intersect(union(procedures,"Z"),c("Pearson","Haberman","Z"))
criteria<-intersect(union(criteria,"Z"),c("Pearson","Haberman","Z"))[1]
}
#Names
if(!("language" %in% names(arguments))) arguments$language <- "en"
nodes <- arguments$nodes
if (inherits(nodes,"tbl_df")) nodes<-as.data.frame(nodes)
name <- arguments$name <- nameByLanguage(arguments$name,arguments$language,arguments$nodes)
if (!("level" %in% names(arguments))) level<-.95 else level <-arguments$level
#Data.frame
if (all(inherits(data,c("tbl_df","tbl","data.frame"),TRUE))) data<-as.data.frame(data) # convert haven objects
if (inherits(weight,"character")) variables <- setdiff(variables,weight)
#ignore constant dichotomies
if(!is.null(dichotomies)){
cons <- character()
for(i in dichotomies){
if(length(unique(data[,i]))==1){
cons <- c(cons,i)
}
}
dichotomies <- setdiff(dichotomies,cons)
if(!length(dichotomies)){
dichotomies <- NULL
}
}
allvar<-union(union(metric,dichotomies),variables)
if (!pairwise & inherits(weight,"character")) {
if (!is.null(weight)) weight <- data[rowSums(is.na(data[,allvar]))<1,weight]
data <- data[complete.cases(data[,allvar]),]
}
if(!is.null(weight)) {
if(inherits(weight,"character")){
weight<-data[,weight]
data<-data[,allvar]
}
else{
if(length(weight)!=dim(data)[1]) stop("Weights have not the correct dimensions")
if (pairwise) data <- cbind(data[,allvar],weight)[,1:length(data[,allvar])]
else data <- na.omit(cbind(data[,allvar],weight))[,1:length(data[,allvar])]
}
}
else data<-data[,allvar]
data[,variables]<-as_factor(data[,variables])
arguments$scenarios <- sum(rowSums(!is.na(data))>0) # Number of scenarios
#Size
if(!("size" %in% names(arguments)))
if(percentages)
arguments$size <- "%"
#Dichotomies
if(!is.null(dichotomies)){
dichos<-dicho(data, dichotomies, valueDicho, newlabel = FALSE)
variables<-setdiff(variables,dichotomies)
}
#Dichotomize
if (length(variables)>0){
incidences<-dichotomize(data, variables, "", min=minimum, length=0, values=NULL, sparse=FALSE, add=FALSE, sort=sort, nas=NULL)
if(!is.null(dichotomies)) incidences<-cbind(dichos,incidences)
}
else if(exists("dichos")) incidences<-dichos
#Nodes filter
if (!is.null(nodes)) {
nonAmong <-setdiff(as.character(nodes[[name]]),names(incidences))
nodeList <-setdiff(as.character(nodes[[name]]),nonAmong)
incidences<-incidences[,nodeList, drop=FALSE]
nonAmongM <- setdiff(metric, as.character(nodes[[name]]))
metric <- setdiff(metric,nonAmongM)
if (length(metric)==0) metric <-NULL
if (length(nonAmong)>0)
warning(paste0(toString(nonAmong)," is/are not present in the data frame."))
# nodes <- nodes[as.character(nodes[[name]]) %in% union(names(incidences),metric),]
}
#Nodes elaboration
if (!exists("incidences") | ncol(incidences)<2) stop("There are no more than 1 qualitative category. Try netCorr.")
if (all(is.na(incidences) | incidences==0 | incidences==1)) {
incidences <- incidences[,names(incidences)[order(match(sub(":.*","",names(incidences)),varOrder))], drop=FALSE]
C<-coin(incidences, minimum, maximum, sort, decreasing, weight=weight, subsample=subsample, pairwise = pairwise)
if(coin) return(C)
O<-asNodes(C,frequency,percentages,arguments$language)
names(O)[1]<-name
O$variable <- sub(":.*","",O[,name])
if(!is.null(nodes)) {
O<-merge(O,nodes[,setdiff(names(nodes),frequencyList),drop=FALSE],by.x=name,by.y=name,all.x=TRUE)
}else {
if (!is.null(commonlabel)) { # Preserve the prename (variable) of a node if specified in commonlabel
arguments$label<-getByLanguage(labelList,arguments$language)
provlabels<-as.character(O[[name]])
O[[arguments$label]]<-ifelse(substr(O[[name]],1,regexpr('\\:',O[[name]])-1) %in% commonlabel,provlabels,substr(O[[name]],regexpr('\\:',O[[name]])+1,1000000L))
}
}
#Links elaboration
E<-edgeList(C, procedures, criteria, level, Bonferroni, minL, maxL, support,
directed, diagonal, sortL, decreasingL, pairwise)
for(lattr in c("lwidth","lweight","lcolor","ltext"))
if(!is.null(arguments[[lattr]])) arguments[[lattr]]<-i.method(c.method(arguments[[lattr]]))
if(!is.null(arguments$layout)) {
layout2 <- layout <- arguments$layout
if (inherits(layout,"matrix") && is.null(metric)){
if (!is.null(nodes)){
if(nrow(layout)==nrow(nodes)){
Oxy <- matrix(NA,nrow(O),2)
rownames(Oxy) <- as.character(O[,name])
rownames(layout) <- as.character(nodes[,name])
layoutnames <- intersect(rownames(Oxy),rownames(layout))
Oxy[layoutnames,] <- layout[layoutnames,]
arguments$layout <- Oxy
} else warning("layout must have a coordinate per node")
} else warning("layout must be applied to the nodes variable")
} else {
if(is.character(layout)){
if(tolower(substr(layout,1,2))=="mc")arguments$layout<-layoutMCA(incidences)
else if(tolower(substr(layout,1,2))=="pc")arguments$layout<-layoutPCA(C)
}
else if(!is.null(metric)) arguments$layout<-NULL # There is metric information and not MCA or PCA
}
}
if(!is.null(metric)) {
#Metric nodes elaboration
if(percentages) O$mean<-O$`%`/100
O$min<-0
O$max<-1
means<-sapply(na.omit(data[,metric, drop=F]),mean)
mins<-sapply(na.omit(data[,metric, drop=F]),min)
maxs<-sapply(na.omit(data[,metric, drop=F]),max)
P<-(means-mins)/(maxs-mins)*100
O2<-data.frame(name=names(means),mean=means,min=mins,max=maxs,P=P,variable=names(means))
colnames(O2)[1] <- name
colnames(O2)[5] <- "%"
if(!is.null(nodes)){
O2 <- O2[as.character(O2[[name]]) %in% as.character(nodes[[name]]),] #nodes filter
for(col in as.character(O2[[name]]))
O2[as.character(O2[[name]])==col,colnames(nodes)] <- nodes[as.character(nodes[[name]])==col,]
}
O<-rbind.all.columns(O,O2)
#Metric links elaboration
methods<-union(procedures,criteria)
if (pairwise) R <- corrp(data[,metric, drop=F], cbind(incidences, data[,metric, drop=F]), weight=weight)
else R <- corr(data[,metric, drop=F], cbind(incidences,data[,metric, drop=F]), weight=weight)
if (nrow(R)==1) row.names(R)<-metric
allvar<-union(as.character(nodes[[name]]),c(names(incidences),metric))
order1<-intersect(allvar,rownames(R))
order2<-intersect(allvar,colnames(R))
R<-R[order1,order2, drop=F]
Pearson<-mats2edges(R)
colnames(Pearson)[3]<-"Pearson"
D<-as.data.frame(Pearson)
if("Haberman" %in% methods){
H<-R*sqrt(nrow(incidences))
Haberman<-mats2edges(H)
D<-cbind(D,Haberman[,3]); colnames(D)[length(D)]<-"Haberman"
}
if("Z" %in% methods) {
t<-R/sqrt((1-pmin(1L,R))/(nrow(incidences)-2))
Z <- mats2edges(1-pt(t,nrow(incidences)-2))
D<-cbind(D,Z[,3]); colnames(D)[length(D)]<-"Z"
}
D<-D[,c("Source","Target",methods)]
D<-D[D[criteria] > minL & D[criteria] < maxL,]
colnames(D)<-sub("^Z$","p(Z)",colnames(D))
if(is.null(E))E<-D
else E<-rbind.all.columns(E,D)
#Layout
if (inherits(layout,"matrix")){
if (!is.null(nodes)){
if(nrow(layout2)==nrow(nodes)){
Oxy <- matrix(NA,nrow(O),2)
rownames(Oxy) <- as.character(O[,name])
rownames(layout2) <- as.character(nodes[,name])
layoutnames <- intersect(rownames(Oxy),rownames(layout2))
Oxy[layoutnames,] <- layout2[layoutnames,]
layout2 <- Oxy
} else warning("layout must have a coordinate per node")
} else warning("layout must be applied to the nodes variable")
arguments$layout <- layout2
}
}
if (!is.null(exogenous)) {
exogenous2<-intersect(exogenous,c(metric,dichotomies))
E$chaine<-ifelse(((substr(E$Source,1,regexpr("\\:",E$Source)-1) %in% exogenous) |
(E$Source %in% exogenous2)) &
((substr(E$Target,1,regexpr("\\:",E$Target)-1) %in% exogenous) |
(E$Target %in% exogenous2)),"No","Yes")
arguments$linkFilter<-paste(ifelse(is.null(arguments$linkFilter),"",paste(arguments$linkFilter,"&")),"chaine=='Yes'")
}
if("showArrows" %in% names(arguments$options) & exists("nodes")) E<-orderEdges(E,nodes[[name]])
if(exists("ltext",arguments)) {
if(toupper(arguments$ltext) == "Z") arguments$ltext <- "p(Z)"
if(arguments$ltext =="Fisher") arguments$ltext <- "p(Fisher)"
}
if(!is.null(dir)){
arguments$dir <- dir
}
arguments$nodes <- O
arguments$links <- E
xNx <- do.call(netCoin,arguments)
if (igraph) {
return(toIgraph(xNx))
} else {
return(xNx)
}
} else warning("Input is not a dichotomous matrix of incidences")
}
# surScat is a wrapper to build a netCoin object from an original non-dichotomized data.frame and see frequencies.
surScat <- function(data, variables=names(data), active=variables, type= c("mca", "pca"), nclusters=2, maxN=2000, ...) {
if(type[1]=="mca") {
B <- as.data.frame(droplevels(as_factor(na.omit(data[,variables]))))
b <- B[,active]
m <- as.matrix(dichotomize(b,variables=names(b), sort=F, add=F, nas=NULL))
cc <- layoutMca(m, rows=T)
}
else {
B <- na.omit(data[,variables])
b <- as.data.frame(sapply(B[,active], as.numeric))
factors <- setdiff(variables, active)
B[, factors] <- as.data.frame(droplevels(as_factor(B[,factors])))
B[, active] <- b
m <- prcomp(b, center = TRUE, scale. = TRUE)
cc <- m$x[,1:2]
}
for(i in nclusters) {
G <- stats::kmeans(cc, centers=i)
g <- paste0("Grupos(",sprintf(paste0("%0",length(levels(G$cluster)),"d"),i),")")
B[[g]] <- paste0("Grupo: ",sprintf(paste0("%0",length(levels(G$cluster)),"d"),G$cluster))
}
arguments <- list(...)
arguments$name <- nameByLanguage(NULL,arguments$language,NULL)
if(is.character(rownames(B))) B[[arguments$name]] <- rownames(B)
else B[[arguments$name]] <- sprintf(paste0("%0",nchar(nrow(B)),"d"),as.numeric(rownames(B)))
B <- B[, c(active, setdiff(names(B), active))]
if(nrow(B)>maxN) {
set.seed(2020)
rcases <- sample(1:nrow(B), maxN)
B <- B[rcases,]
cc <- cc[rcases,]
}
arguments$nodes <- B
arguments$layout <- cc
arguments$color <- g
arguments$frequencies <- TRUE
arguments$showAxes <- TRUE
arguments$showCoordinates <- TRUE
arguments$degreeFilter <- NULL
if(is.null(arguments$label)) arguments$label <- ""
if(is.null(arguments$controls)) arguments$controls <- c(1,4)
return(do.call(netCoin, arguments))
}
# produces a netCoin object to graph a CoWeb graphic (Upton 2000).
cobCoin<- function(data, variables=names(data), degree=0, significance=.05, ...) {
arguments <- list(...)
if(!exists("color", arguments)) arguments$color<-"variable"
if(!exists("group", arguments)) arguments$group<-"variable"
if(!exists("lwidth", arguments)) arguments$lwidth<-"Haberman"
if(!exists("lcolor", arguments)) arguments$lcolor<-"Haberman"
if(!exists("linkBipolar", arguments)) arguments$linkBipolar <- TRUE
arguments$groupText<-TRUE
arguments$data <- data[,variables]
arguments$maxL <- 1
arguments$commonlabel <- ""
N <- do.call(surCoin, arguments)
control <- N$links[,1:2]
control[[1]] <- gsub(":.*","",control[[1]])
control[[2]] <- gsub(":.*","",control[[2]])
N$links <- N$links[control[[1]]!=control[[2]],]
N$links <- N$links[N$links$`p(Z)` < significance/2 | N$links$`p(Z)` > (1-significance/2),]
N <- addNetCoin(N, layout=layoutCircle(N$nodes, deg=degree, name=N$options$nodeName))
return(N)
}
# Elaborate a netCoin object from a lavaan object.
pathCoin<-function(model, estimates=c("b","se","z","pvalue","beta"), fitMeasures=c("chisq", "cfi", "rmsea"), ...){
arguments <- list(...)
if(!("language" %in% names(arguments))) arguments$language <- "en"
if(!("linkBipolar" %in% names(arguments))) arguments$linkBipolar <- TRUE
if(!("showArrows" %in% names(arguments))) arguments$showArrows <- TRUE
arguments$name <- nameByLanguage(arguments$name,arguments$language,arguments$nodes)
M<-pathParameter(model,estimates=estimates)
names(M$nodes)[1]<-arguments$name
if(!is.null(arguments$note)) arguments$note<-paste0(catFit(model,fitMeasures),arguments$note)
else if(!is.null(fitMeasures)) arguments$note<-catFit(model,fitMeasures)
if("nodes" %in% names(arguments)) {
vvnodes<-setdiff(names(arguments$nodes),arguments$name)
arguments$nodes<-merge(arguments$nodes,M$nodes,by.x=arguments$name,by.y=arguments$name,sort=F)
arguments$nodes<-arguments$nodes[,c(names(M$nodes),vvnodes)]
}
else arguments$nodes<-M$nodes
arguments$links<-M$links
do.call(netCoin,arguments)
}
# dichotomize: Transform character and factor objects into dichotomies.
dichotomize <- function(data,variables, sep="", min=1, length=0, values=NULL,
sparse=FALSE, add=TRUE, sort=TRUE, nas="None") {
if(is.data.frame(data)){
if (min>0 & min<1) min = min*nrow(data)
if(!is.null(nas) & sep=="C") { #Different columns, same factor levels or character.
oldData<-data
if(inherits(data,"tbl_df")) data<-as_factor(data[,variables])
if(sep=="C")sep<-"|"
String<-do.call(paste,c(data,sep=sep))
String<-gsub(paste0("\\",sep,"NA"), "", String)
String<-gsub("NA", nas, String)
data<-data.frame(String=String)
variables<-"String"
}
cn <- colnames(data)
names(cn) <- cn
if (length(sep)!=length(variables)) sep = rep(sep[1],length(variables))
names(sep) <- variables
for(c in variables){
if (sep[c]!="")
L <- lapply(strsplit(as.character(data[[c]]), sep[c],fixed=TRUE), paste, sep[c] ,sep="") # Sep at the end of each element
else {
sep[c]<-"??"
L <- paste(na.omit(data[[c]]),sep[c],sep="")
}
if (is.null(values)) Z <- valuesof(L,length,min,sort,sep)
else Z <- paste(values, sep[c],sep="")
Z <- paste(sep[c],Z,sep="") # To search
C <- paste(sep[c],data[[c]],sep[c],sep="") # Searched
Q <- Matrix(0,nrow=length(Z),ncol=nrow(data),sparse=TRUE) # Result
for(X in 1:length(Z)) {
N <- grep(Z[X],C,fixed=TRUE)
if (length(N) > 0) Q[X,N]<-1
}
Z <- substring(Z,nchar(sep[c])+1,nchar(Z)-nchar(sep[c]))
Z[Z==""]<-".ND."
if (min>0 & !is.null(values)) {
VF <- apply(Q,1,sum)
if (length(VF[VF>=min])>1) {
Q <- Q[VF>=min,]
Z <- Z[VF>=min]
}
else warning("min > empirical observations")
if (sort==TRUE) {
Q<-Q[order(-VF[VF>=min]),]
Z<-Z[order(-VF[VF>=min])]
}
if (length>0) {
Q<-Q[1:min(length,nrow(Q)),]
Z<-Z[1:min(length,nrow(Q))]
}
}
if(length(variables)>1)
Z <- paste(cn[c],Z,sep=":")
if (sparse==TRUE) {
Q<-t(Q)
colnames(Q)<-Z
if (!exists("R_Data")) R_Data<-Q
else R_Data<-cbind(R_Data,Q)
}
else {
W <- as.data.frame(t(as.data.frame(as.matrix(Q), row.names=Z)))
if (!is.null(nas)) W[[paste0(c,":<NA>")]] <- ifelse(apply(W,1,sum)==0,1,0)
W[apply(W,1,sum)==0,] <- NA
if (add==TRUE) {
if (exists("oldData")) data<-oldData
data <- cbind(data,W)
}
else {
if (!exists("R_Data")) R_Data<-W
else R_Data<-cbind(R_Data,W)
}
}
}
if (!exists("R_Data")) R_Data<-data
if (!is.null(nas) & sparse==FALSE) {
names(R_Data)<-gsub("\\:<NA>",paste0("\\:",nas),names(R_Data))
names(R_Data)<-gsub("\\.ND\\.",nas,names(R_Data))
}
return(R_Data)
}
else warning("You must pass a data frame!")
}
valuesof<-function(x,length=0,min=0,sort=TRUE,sep="") {
x <- table(unlist(x))
if (sort) x <- x[order(-x)]
if (min>0) x<-as.matrix(x[x>=min])
if (nrow(x)==0) return(NULL)
if (length>0) x<- x[x>=(x[order(-x)][min(length,nrow(x))])]
return(rownames(x))
}
# Links. See below funcs="shape"
edgeList <- function(data, procedures="Haberman", criteria="Z", level = .95, Bonferroni=FALSE, min=-Inf, max=Inf, support=-Inf,
directed=FALSE, diagonal=FALSE, sort=NULL, decreasing=TRUE, pairwise=FALSE) {
level <- checkLevel(level)
if (tolower(substr(criteria,1,2))%in%c("z","hy") & substr(tolower(procedures[1]),1,2)!="sh") {
if (max==Inf) max<-.50
if (Bonferroni ) max<-max/choose(nrow(data),2) # Changes of Z max criterium (Bonferroni)
}
if (substr(tolower(procedures)[1],1,2)!="sh") { # For coin objects
if (!inherits(data,"coin")) stop("Error: input must be a coin object (see coin function)")
funcs<-c.method(procedures)
if(!is.null(sort)) funcs<-union(c.method(sort),funcs)
criteria<-c.method(criteria)
todas<-union(funcs,criteria)
matrices<-sim(data,todas,level=level, pairwise=pairwise)
funcs<-i.method(funcs)
criteria<-i.method(criteria)
if (length(union(funcs,criteria))==1) {
M<-new.env()
M[[funcs]]<-matrices
matrices<-as.list(M)
}
matrices<-matrices[i.method(todas)]
Mat<-mats2edges(data[,],matrices,criteria=criteria,min=min,max=max,support=support,directed=directed,diagonal=diagonal)
}
else {
if (!inherits(data,"matrix") && !inherits(data,"data.frame"))
stop("Error: input must be a matrix (shape) or a data.frame (tree)")
if (inherits(data,"matrix")){
if(min==-Inf)min<-1
#funcs="value"
#M<-new.env()
#M[[criteria]]<-M[[funcs]]<-data
#matrices<-as.list(M)
#data<-list(f=M[[funcs]],n=NA)
Mat<-mats2edges(data,min=min,max=max,directed=directed,diagonal=diagonal)
}
if (inherits(data,"data.frame")) {
lines<-sapply(data,as.character)
lines<-rbind(c(lines[1,1],rep(NA,ncol(lines)-1)),lines) # Add one blank case in order to avoid mXm problem.
lines<-ifelse(lines=="",NA,lines)
adjlist<-split(lines,seq(nrow(lines))) # splits the character strings into list with different vector for each line
adjlist<-sapply(adjlist,na.omit)
Source=unlist(lapply(adjlist,function(x) rep(x[1],length(x)-1))) # establish first column of edgelist by replicating the 1st element (=ID number) by the length of the line minus 1 (itself)
Target=unlist(lapply(adjlist,"[",-1)) # the second line I actually don't fully understand this command, but it takes the rest of the ID numbers in the character string and transposes it to list vertically
return(as.data.frame(cbind(Source,Target),stringsAsFactors = FALSE,row.names=FALSE))
}
}
# Last transformations: c.Conditional c.Probable and sort
if(length(Mat)>0) {
if (!is.null(Mat$c.conditional))
Mat$c.conditional<-factor(Mat$c.conditional,levels=c(0:8),
labels=c("Null","Mere","Conditional","Significant","Quite significant","Very significant","Subtotal","Suptotal","Total"))
if (!is.null(Mat$c.probable))
Mat$c.probable<-factor(Mat$c.probable,levels=c(0:8),
labels=c("Null","Mere","Probable","Significant","Quite significant","Very significant","Subtotal","Suptotal","Total"))
if (!is.null(sort)) {
if (substr(tolower(procedures)[1],1,2)!="sh") Mat<-Mat[order(Mat[[i.method(c.method(sort))]],decreasing = decreasing),]
else Mat<-Mat[order(Mat$value,decreasing=decreasing),]
}
}
else return(NULL)
names(Mat)[names(Mat)=="Z"]<-"p(Z)"
names(Mat)[names(Mat)=="Fisher"]<-"p(Fisher)"
return(Mat)
}
c.method<-function(method) {
if(is.null(method))return(NULL)
method<-toupper(method)
if ("ALL"==method[1]) method<-c("M","T","G","S","B","J","D","A","O","K","N","Y","P","V","C","R","E","L","H","Z","f","F","X","I","0","Q","U")
method<-sub("FIS","W",method)
method<-sub("HYP","W",method)
method<-sub("HAM","NNH",method)
method<-sub("CC" ,"UC" ,method)
method<-sub("CP" ,"QC" ,method)
method<-sub("OD" ,"CD" ,method)
method<-sub("RO" ,"TA" ,method)
method<-sub("AND","BER",method)
method<-sub("TET","VET",method)
method<-sub("CON","LCO",method)
method<-sub("II" ,"0" ,method)
method<-sub("COI","F", method)
method<-substr(method,1,1)
return(method)
}
i.method<-function(method) {
similarities<-matrix(c("matching","Rogers","Gower","Sneath", "Anderberg",
"Jaccard","dice", "antiDice","Ochiai","Kulczynski",
"Hamann", "Yule", "Pearson", "odds", "Russell", "expected", "Haberman", "confidence", "Z",
"coincidences", "relative", "sConditional","tConditional", "c.conditional","c.probable","tetrachoric","Fisher"),
nrow=1, dimnames=list("Similarity", c("M","T","G","S","B","J","D","A","O","K","N","Y","P","C","R","E","H","L","Z","F","X","I","0","U","Q","V","W")))
similarities<-similarities[,method]
if("L" %in% method) {
cade <- c("conf.L","confidence","conf.U")
posi <- match("confidence",similarities)
if (posi ==1) similarities <- c(cade,similarities[-1])
else if (posi == length(similarities)) similarities<-c(similarities[-length(similarities)], cade)
else similarities <- c(similarities[c(1:max(1,posi-1))], cade,
similarities[(match("confidence",similarities)+1):length(similarities)])
}
return(similarities)
}
r.method <- function(method) {
if(!toupper(substr(method,1,1)) %in% c("P", "K", "S")) method <- "pearson"
rs <- matrix(c("pearson","kendall", "spearman"), nrow=1, dimnames=list("method", c("P", "K", "S")))
method <- rs[,toupper(substr(method, 1, 1))]
return(method)
}
checkLevel <- function(level){
if (level >=1 & level < 100)
level <- level/100
if (level <=0 | level >=100) {
level <- .95
warning("Not valid level")
}
return(level)
}
# Similatities
sim<-function (input, procedures="Jaccard", level=.95, distance=FALSE, minimum=1, maximum=Inf, sort=FALSE, decreasing=FALSE,
weight=NULL, pairwise=FALSE) {
level <- checkLevel(level)
method<-c.method(procedures)
if (is.matrix(input) && !inherits(input,"coin")) {
if (is.null(colnames(input))) dimnames(input)<-list(NULL,paste("X",1:ncol(input),sep=""))
input<-as.data.frame(input)
}
if (is.data.frame(input)) {
C<-coin(input, minimum, maximum, sort, decreasing, weight=weight, pairwise=pairwise)
a<-C[,]
}
else if (inherits(input,"coin")) {
C <- input
a <- input[(diag(input)>=minimum & diag(input)<=maximum),(diag(input)>=minimum &diag(input)<=maximum)]
}
else stop("Error: 1st parameter has to be a data frame or a coin object (see coin function)")
if (pairwise) {
N <- attr(C, "m")
NN <- max(C)
X <- attr(C, "x")
}
else {
N <- attr(C,"n")
X <- matrix(0, nrow=dim(a)[1], ncol = dim(a)[2])
}
if (sort==TRUE | decreasing==TRUE) {
orderM <- order(diag(a),decreasing=decreasing)
a <- a[orderM , orderM]
if (pairwise) N <- N[orderM , orderM]
if (pairwise) X <- X[orderM , orderM]
}
b<--(a-diag(a))-t(X)
c<-t(b)
d=N-a-b-c
if (any(c("G","B","O","K","Y","P","V") %in% method)) m<-ifelse(a+d==N,1,ifelse(b+c==N,-1,0)) #Special values of distances
s<-new.env()
if ("M" %in% method) s$matching <- distant((a + d)/(a + b + c + d),distance)
if ("T" %in% method) s$Rogers <- distant((a + d)/(a + 2 * (b + c) + d), distance)
if ("G" %in% method) {
s$Gower <- distant(a * d/sqrt((a + b) * (a + c) * (d + b) * (d + c)),distance)
s$Gower <- ifelse(is.na(s$Gower),distant(pmax(m,0),distance),s$Gower)
}
if ("S" %in% method) s$Sneath <- distant(2*(a+d)/(2*(a+d)+(b+c)),distance)
if ("B" %in% method) {
s$Anderberg <- distant((a/(a+b)+a/(a+c)+d/(c+d)+d/(b+d))/4,distance)
s$Anderberg <- ifelse(is.na(s$Anderberg),distant(pmax(m,0),distance),s$Anderberg)
}
if ("J" %in% method) s$Jaccard <- distant(a/(a + b + c),distance)
if ("D" %in% method) s$dice <- distant(2 * a/(2 * a + b + c), distance)
if ("A" %in% method) s$antiDice <- distant(a/(a + 2 * (b + c)),distance)
if ("O" %in% method) {
s$Ochiai <- distant(a/sqrt((a + b) * (a + c)),distance)
s$Ochiai <- ifelse(is.na(s$Ochiai),distant(pmax(m,0),distance),s$Ochiai)
}
if ("K" %in% method) {
s$Kulczynski <- distant((a/(a+b)+a/(a+c))/2, distance)
s$Kulczynski <- ifelse(is.na(s$Kulczynski),distant(pmax(m,0),distance),s$Kulczynski)
}
if ("N" %in% method) s$Hamann <- distant((a - (b + c) + d)/(a + b + c + d), distance)
if ("Y" %in% method) {
s$Yule <- distant((a*d-b*c)/(a*d+b*c))
s$Yule <- ifelse(is.na(s$Yule),distant(m,distance),s$Yule)
}
if ("P" %in% method) {
s$Pearson <- distant((a * d - b * c)/sqrt((a + b) * (a + c) * (b + d) * (d + c)),distance)
s$Pearson <- ifelse(is.na(s$Pearson),distant(m,distance),s$Pearson)
}
if ("V" %in% method) {
s$tetrachoric<-((a*d/(b*c))^(pi/4)-1)/((a*d/(b*c))^(pi/4)+1)
s$tetrachoric <- ifelse(is.na(s$tetrachoric),distant(m,distance),s$tetrachoric)
}
if ("C" %in% method) {
s$odds <- (pmax(a,.5)*pmax(d,.5))/(pmax(b,.5)*pmax(c,.5))
if (distance) s$odds<--s$odds
diag(s$odds)<-ifelse(distance,-Inf,Inf)
}
if ("R" %in% method) {
s$Russell <- distant(a/(a + b + c + d),distance)
if (!distance) diag(s$Russell) <- 1
}
if ("E" %in% method) s$expected <- (a+b)*(a+c)/N
if ("L" %in% method) {
s$'conf.L' <- pmax(a-qt(level+(1-level)/2, N-1)*sqrt(((a+b)*(a+c)/N)*((1-(a+b)/N)*(1-(a+c)/N))),0)
signo<-2*(((a+b)*(a+c)/N)<a)-1
s$confidence <- pmax((a+b)*(a+c)/N+signo*qt(level,N-1)*sqrt(((a+b)*(a+c)/N)*((1-(a+b)/N)*(1-(a+c)/N))),0)
diag(s$confidence) <- diag(a)
s$'conf.U' <- pmin(a+qt(level+(1-level)/2, N-1)*sqrt(((a+b)*(a+c)/N)*((1-(a+b)/N)*(1-(a+c)/N))),N)
}
if ("H" %in% method) {
s$Haberman <- sqrt(N) * (a * d - b * c)/sqrt((a + b) * (a + c) * (b + d) * (d + c))
if (pairwise) s$Haberman[is.na(s$Haberman)]<-sqrt(NN)
else s$Haberman[is.na(s$Haberman)]<-sqrt(N)[is.na(s$Haberman)]
if (distance) s$Haberman<-(N+s$Haberman)/(2*N)
}
if ("Z" %in% method) {
s$Z <- 1-pt(sqrt(N) * (a * d - b * c)/sqrt((a + b) * (a + c) * (b + d) * (d + c)),N)
s$Z[is.na(s$Z)]<-0
}
if ("W" %in% method) s$Fisher<-1-phyper(a-1,pmin((a+b),(a+c)),N-pmin((a+b),(a+c)),pmax((a+b),(a+c)))
if ("x" %in% method) {
s$Fisher<-matrix(NA,nrow=nrow(a),ncol=ncol(a))
for (Ro in c(1:nrow(a))) {
for (Co in c(Ro:ncol(a))) {
inMatrix=matrix(c(a[Ro,Co],b[Ro,Co],c[Ro,Co],d[Ro,Co]),nrow=2)
(s$Fisher[Ro,Co]<-fisher.test(inMatrix,alternative="greater")$p.value)
s$Fisher[Co,Ro] = s$Fisher[Ro,Co]
}
}
# diag(s$Fisher)<-0
rownames(s$Fisher)<-colnames(s$Fisher)<-rownames(a)
}
if ("F" %in% method) s$coincidences <- a
if ("X" %in% method) s$relative <- a/N*100
if ("I" %in% method) s$sConditional <-a/(a+c)*100
if ("0" %in% method) s$tConditional <-a/(a+b)*100
if ("U" %in% method) {
Z <- 1-pt(sqrt(N) * (a * d - b * c)/sqrt((a + b) * (a + c) * (b + d) * (d + c)),N)
s$c.conditional<-matrix(ifelse(b+c==0, 8,
ifelse(c==0, 7,
ifelse(b==0, 6,
ifelse(Z<.001,5,
ifelse(Z<.01, 4,
ifelse(Z<.05, 3,
ifelse(Z<.50, 2,
ifelse(a>0, 1, 0)))))))),nrow=nrow(a),dimnames=dimnames(a))
}
if ("Q" %in% method) {
Z <- 1-pt((a/(a+c)-.50)/(1/(2*sqrt(a+c))),(a+c))
s$c.probable<-matrix(ifelse(b+c==0, 8,
ifelse(c==0, 7,
# ifelse(b==0, "Subtotal",
ifelse(Z<.001, 5,
ifelse(Z<.01, 4,
ifelse(Z<.05, 3,
ifelse(Z<.50, 2,
ifelse(a>0, 1, 0))))))),nrow=nrow(a),dimnames=dimnames(a))
}
if (length(method)==1) return(s[[names(s)[1]]])
else return(as.list(s,sorted=TRUE))
}
# Convert similatities into dissimilarities
distant<-function(s,t=FALSE) {
if (t==TRUE) s<-as.dist(1-s)
return(s)
}
# http://pbil.univ-lyon1.fr/ade4/ade4-html/dist.binary.html
# Print lower matrices
lower<-function(matrix,decimals=3) { # Add an option to hiden diagonal
m<-as.matrix(matrix)
form=paste("%1.",decimals,"f",sep="")
lower<-apply(m,1,function(x) sprintf(form,x))
lower[upper.tri(lower)]<-""
lower<-as.data.frame(lower, stringsAsFactors=FALSE)
if (ncol(m)==1) rownames(lower)<-colnames(lower)<-names(matrix)
rownames(lower)<-names(lower)
return(lower)
}
# List of coincidences
coin<-function(incidences,minimum=1, maximum=nrow(incidences), sort=FALSE, decreasing=TRUE,
total=FALSE, subsample=FALSE, weight=NULL, pairwise=FALSE) {
if (pairwise){
n <- sum(rowSums(!is.na(incidences))>0)
nomiss <- ifelse(is.na(incidences),0,1)
incidences <- replace(incidences,is.na(incidences),0)
}
else {
if (!is.null(weight)) weight <- weight[rowSums(is.na(incidences))<1]
incidences<-na.omit(incidences)
}
if (subsample){
vector<-apply(incidences,1,sum)
incidences<-incidences[vector>0,]
}
if (total & is.null(weight)) incidences<-data.frame(Total=1,incidences)
if (all(is.na(incidences) | incidences==0 | incidences==1)) {
if (!pairwise) n <- nrow(incidences)
names(n)<-"n"
if (is.null(weight)) f<-crossprod(as.matrix(incidences))
else {
if (length(weight)!=dim(incidences)[1]) warning("weight has not the appropiate length!")
f<-crossprod(t(crossprod(as.matrix(incidences),diag(weight,length(weight)))),as.matrix(incidences))
n<-maximum<-sum(weight)
}
if (is.null(colnames(f))) dimnames(f)<-list(paste("X",1:ncol(f),sep=""),paste("X",1:ncol(f),sep=""))
d<-diag(f)
if (sort) d<-sort(d,decreasing=decreasing)
S<-names(d[(d>=minimum & d<=maximum)])
if (total & is.null(weight)) S<-c("Total",S)
if (total & !is.null(weight)) warning("total cannot be applied in weighted tables")
if (length(S)>0) {
if (!pairwise) structure(f[S,S], n=n, class=c("coin"))
else {
colnames(nomiss) <- colnames(incidences)
nomiss <- nomiss[,S]
incidences <- incidences[,S]
if (is.null(weight)) {
m<-crossprod(as.matrix(nomiss))
x<-crossprod(1-as.matrix(nomiss),as.matrix(incidences))
}
else {
if (length(weight)!=dim(nomiss)[1]) warning("weight has not the appropiate length!")
m<-crossprod(t(crossprod(as.matrix(nomiss),diag(weight,length(weight)))),as.matrix(nomiss))
x<-crossprod(t(crossprod(1-as.matrix(nomiss),diag(weight,length(weight)))),as.matrix(incidences))
}
structure(f[S,S], n=n, m=m[S,S], x=x[S,S], class=c("coin"))
}
}
else cat("No variables left")
}
else warning("All data in incidence matrix has to be dichotomous.")
}
coocur<-function (ocurrences, minimum = 1, maximum = Inf, sort = FALSE, decreasing=TRUE)
{
result <- matrix(nrow = ncol(ocurrences), ncol = ncol(ocurrences), dimnames = list(colnames(ocurrences),
colnames(ocurrences)))
for (val in c(1:nrow(result))) {
for (cal in c(val:ncol(result))) {
result[cal, val] <- sum(pmin(ocurrences[, cal], ocurrences[, val]))
if (val != cal)
result[val, cal] = result[cal, val]
}
}
d<-diag(result)
if (sort) d<-sort(d,decreasing=decreasing)
S<-names(d[(d>=minimum & d<=maximum)])
if (length(S)>=1) {
result<-result[S,S]
if (length(S)==1) names(result)<-S
n <- sum(ocurrences[,S])
m <- sum(apply(as.matrix(ocurrences[,S]), 1, max))
attr(result, "n") <- n
attr(result, "m") <- m
structure(result, class = "cooc")
}
else cat("No variables left")
}
print.coin<-function(x, ...) {
cat("n= ",attr(x,"n"),"\n",sep="")
print(lower(x,0))
}
print.cooc<-function(x, ...) {
cat("n= ",attr(x,"n"),"; m= ", attr(x,"m"),"\n",sep="")
print(lower(x,0))
}
tempDir <- function(){
dir.create("temp", showWarnings = FALSE)
return(paste("temp",round(as.numeric(Sys.time())),sep="/"))
}
plot.coin <- function(x, dir=tempDir(), language=c("en","es","ca"), ...){
N <- asNodes(x, language = language)
colnames(N)[2] <- "incidences"
E <- edgeList(x,c("Frequencies","Expected"))
bar <- barplot_rd3(N, E, name = names(N)[1], coincidences = "coincidences", incidences = "incidences", expected = "expected", cex = 1, language = language)
pie <- pieCoin(x, cex = 1, language = language)
rd3_multigraph(barplot=bar, piechart=pie, dir = dir)
}
summary.coin <- function(object, ...){
cat(attr(object,"n"),"scenarios and", dim(object[,])[1], "events\n")
diag(object[,])/attr(object,"n")
}
summary.netCoin <- function(object, ...){
summaryNet(object)
}
summary.barCoin <- function(object, ...){
summaryNet(object)
}
summaryNet <- function(x){
cat(dim(x$nodes)[1], "nodes and", dim(x$links)[1], "links.\n")
freq <- frequencyList[frequencyList %in% names(x$nodes)]
if(length(freq)==1) {
cat(freq," distribution of nodes:","\n", sep="")
print(summary(x$nodes[[freq]]))
}
lwidth <- NULL
if(!is.null(x$options$linkWidth))
lwidth <- x$options$linkWidth
else if(length(x$links)>2)
lwidth <- names(x$links)[3]
if(length(lwidth)==1){
cat(lwidth, "'s distribution:","\n",sep="")
print(summary(x$links[[lwidth]]))
}
}
propCoin<-function(x, margin= 0, decimals=1) {
if (!inherits(x,"coin")) stop("Error: input must be a coin object (see coin function)")
if ("m" %in% names(attributes(x))) n <- attr(x, "m")
else n <- attr(x,"n")
x <- x[,]
switch(format(margin),
"0" = round(100 * x / n, decimals),
"1" = round(100 * x / diag(x), decimals),
"2" = round(t(100 * x / diag(x)), decimals))
}
# Transform a coin object into a data frame with name and frequency
asNodes<-function(C, frequency = TRUE, percentages = FALSE, language = c("en","es","ca")){
nodes <- NULL
if (inherits(C,"coin")) {
if ("m" %in% names(attributes(C))) divider <- diag(attr(C,"m"))
else divider <- attr(C,"n")
if (!percentages & frequency) nodes<-data.frame(name=as.character(colnames(C)),frequency=diag(C))
else if (!frequency & percentages) nodes<-data.frame(name=as.character(colnames(C)),"%"=diag(C)/divider*100,check.names=FALSE)
else if (percentages & frequency)nodes<-data.frame(name=as.character(colnames(C)),frequency=diag(C), "%"=diag(C)/divider*100,check.names=FALSE)
else nodes<-data.frame(name=as.character(colnames(C)),check.names=FALSE)
if(language[1]!="en"){
colnames(nodes)[colnames(nodes)=="frequency"] <- getByLanguage(frequencyList,language)
colnames(nodes)[colnames(nodes)=="name"] <- getByLanguage(nameList,language)
}
}
else if (min(c("Source", "Target") %in% names(C))) nodes<-data.frame(name=sort(union(C$Source,C$Target)))
else warning("Is neither a coin object or an edge data frame")
return(nodes)
}
rescale <- function(x) {
to <- c(0, 1)
from <- range(x, na.rm = TRUE, finite = TRUE)
return((x - from[1]) / diff(from) * diff(to) + to[1])
}
toColorScale <- function(items){
if(is.numeric(items)){
return(hsv(1,1,rescale(items)))
}else{
colors <- c(
"#1f77b4", # blue
"#2ca02c", # green
"#d62728", # red
"#9467bd", # purple
"#ff7f0e", # orange
"#8c564b", # brown
"#e377c2", # pink
"#7f7f7f", # grey
"#bcbd22", # lime
"#17becf", # cyan
"#aec7e8", # light blue
"#98df8a", # light green
"#ff9896", # light red
"#c5b0d5", # light purple
"#ffbb78", # light orange
"#c49c94", # light brown
"#f7b6d2", # light pink
"#c7c7c7", # light grey
"#dbdb8d", # light lime
"#9edae5" # light cyan
)
items <- as.numeric(as.factor(items))
items <- ((items-1) %% length(colors))+1
return(colors[items])
}
}
savePajek<-function(net, file="file.net", arcs=NULL, edges=NULL, partitions= NULL, vectors=NULL){
if(length(setdiff(partitions,names(net[["nodes"]])))>0) stop("At least one partition is not amongst ",paste(names(net$nodes),collapse=", "),".")
if(length(setdiff(vectors,names(net[["nodes"]])))>0) stop("At least one vector is not amongst ",paste(names(net$nodes),collapse=", "),".")
if(length(setdiff(arcs,names(net[["links"]])))>0) stop("At least one arc is not amongst ",paste(names(net$links),collapse=", "),".")
if(length(setdiff(edges,names(net[["links"]])))>0) stop("At least one edge is not amongst ",paste(names(net$links),collapse=", "),".")
if(!grepl("\\.",file))file<-paste0(file,".net")
if(!is.null(vectors) | !is.null(partitions)) file<-gsub(".net",".paj",file)
connec<-file(file,"w")
writeLines(paste0("*Network ",net[["options"]]$main),con=connec)
close(connec)
connec<-file(file,"a")
writeLines(paste0("*Vertices ",as.character(nrow(net[["nodes"]]))),con=connec)
writeLines(paste0(seq(1:nrow(net[["nodes"]])),' "',net[["nodes"]]$name,'" '), con=connec)
N<-cbind(n=seq(1:nrow(net[["nodes"]])),net[["nodes"]][1])
L<-cbind(N[unlist(net[["links"]]$Source),1],N[unlist(net[["links"]]$Target),1])
if(!is.null(arcs)) {
cont=1
for(weights in arcs) {
writeLines(paste0("*Arcs : ",cont,' "',weights,'"'), con=connec)
writeLines(paste(L[,1],L[,2],net[["links"]][[weights]]), con=connec)
cont=cont+1
}
}
if(!is.null(edges)) {
ifelse(exists("cont"),cont<-cont,cont<-1)
for(weights in edges) {
writeLines(paste0("*Edges : ",cont,' "',weights,'"'), con=connec)
writeLines(paste(L[,1],L[,2],net[["links"]][[weights]]), con=connec)
cont=cont+1
}
}
if(!is.null(partitions)){
for(partition in partitions) {
writeLines(paste0("*Partition ", partition), con=connec)
writeLines(paste0("*Vertices ", nrow(net$nodes)), con=connec)
writeLines(as.character(as.numeric(as.factor(net[["nodes"]][[partition]]))),con=connec)
}
}
if(!is.null(vectors)){
for(vector in vectors) {
writeLines(paste0("*Vector ", vector), con=connec)
writeLines(paste0("*Vertices ", nrow(net$nodes)), con=connec)
Line<-as.character(net[["nodes"]][[vector]])
Line[is.na(Line)]<-"0"
writeLines(Line, con=connec)
}
}
close(connec)
}
saveGhml <- function(net, file="netCoin.graphml"){
if(!inherits(net, "netCoin")) stop("This program only works with netCoin objects")
if(!grepl("\\.",file))file<-paste0(file,".graphml")
graph <- toIgraph(net)
write_graph(graph, file=file, format="graphml")
}
expectedList<- function(data, names=NULL, min=1, confidence=FALSE) {
if (!inherits(data,"coin")) stop("Error: input must be a coin object")
if (!is.null(names)) colnames(data[,])<-rownames(data[,])<-names
a<-data[,]
b<--(a-diag(a))
c<--t((t(a)-diag(a)))
d=attr(data,"n")-a-b-c
attr(data,"e")<-(a+b)*(a+c)/(a+b+c+d)
E<-edgeList(attr(data,"e"),"shape",min=0,max=Inf)
F<-edgeList(data[,],"shape",min=0,max=Inf)
if (!confidence) {
dataL<-cbind(F,E[,3])[F[,3]>=min,]
colnames(dataL)[3:4]<-c("coincidences","expected")
}
else {
N<-a+b+c+d
signo<-2*(((a+b)*(a+c)/N)<a)-1
attr(data,"l") <- pmax((a+b)*(a+c)/N+signo*1.64*sqrt(((a+b)*(a+c)/N)*((1-(a+b)/N)*(1-(a+c)/N))),0)
diag(attr(data,"l")) <- diag(a)
L<-edgeList(attr(data,"l"),"shape",min=-Inf,max=Inf)
dataL<-cbind(F,E[,3],L[,3])[F[,3]>=min,]
colnames(dataL)[3:5]<-c("coincidences","expected","confidence")
}
return(dataL)
}
mats2edges<-function(data,list=NULL,criteria=1,min=-Inf,max=Inf,support=-Inf,directed=FALSE,diagonal=FALSE){
# Input control
if (!is.null(list)) {
if (!identical(dim(data),dim(list[[1]]))) {
warning("data & list must have the same dimensions")
return()
}
if (criteria!=1 & !(criteria %in% names(list))) warning("criteria are not in the matrices list")
}
if (is.null(rownames(data))) rownames(data)<-as.character(1:nrow(data))
if (is.null(colnames(data))) colnames(data)<-as.character(1:ncol(data))
# type of edgelist
if (nrow(data)<ncol(data)) { # For asymmetric matrices. Improve later for directed cases.
l<-data>-Inf
l[,rownames(l)]<-lower.tri(l[,rownames(l)])
}
else { # For symmetric matrices
if (directed) l <- as.vector(lower.tri(data,diag=diagonal) | upper.tri(data))
else l <- as.vector(lower.tri(data,diag=diagonal))
}
# data.frame building
sources<-rep(colnames(data),each=dim(data)[1])
targets<-rep(rownames(data),dim(data)[2])
Mat <- data.frame(Source=sources,Target=targets)
value <- as.vector(data)
# data.frame for matrices list
if (!is.null(list)) {
c <- as.vector(list[[criteria]])
a<-as.data.frame(lapply(list,as.vector))
Mat<-cbind(Mat,a)[l==TRUE & c<=max & c>=min & value>=support,]
}
# data.frame for alone matrix
else Mat<-{
if (criteria!=1) warning("Criteria don't apply for only one matrix. Use min and max")
cbind(Mat,value)[l==TRUE & value<=max & value>=min,]
}
# return
if (nrow(Mat)>0) {
row.names(Mat)<-NULL
return(Mat)
}
else return(NULL)
}
orderEdges<-function(links,nodes){ #Used in surCoin to order arrows
A<-unlist(sapply(paste0("^",links[,"Source"],"$"),grep,x=nodes))
B<-unlist(sapply(paste0("^",links[,"Target"],"$"),grep,x=nodes))
links[A>B,c("Source","Target")]<-links[A>B,c("Target","Source")]
links<-links[!is.na(links$Source) & !is.na(links$Target),]
return(links)
}
dicho<-function(input,variables,value,newlabel=TRUE) {
if(length(value)!=length(variables)){
value <- rep(value[1],length(variables))
}
for(i in seq_along(variables)){
dicho <- variables[i]
vector0 <- rep(0,nrow(input))
vector0[as.character(input[,dicho])==value[i]] <- 1
input[,dicho] <- vector0
}
datum <- as.data.frame(input[, variables, drop=FALSE])
j=0
for (i in variables) {
j=j+1
if (!is.null(attributes(input[[i]]))) {
if (newlabel) names(datum)[j] <- ifelse(exists("label",attributes(input[[i]])),attr(input[[i]],"label"),i)
}
}
return(datum)
}
languages <- c("en","es","ca")
nameList <- c('name','nombre','nom')
names(nameList) <- languages
frequencyList <- c("frequency","frecuencia","freq\uFC\uE8ncia","%")
names(frequencyList) <- languages
labelList <- c('label','etiqueta','etiqueta')
names(labelList) <- languages
nameByLanguage <- function(name,language,nodes){
if(is.null(name)){
name <- getByLanguage(nameList,language)
}
if(!is.null(nodes)){
if(!(name %in% colnames(nodes)))
warning(paste0("name: '",name,"' column missing in nodes data frame"))
else if(sum(duplicated(nodes[[name]])))
warning(paste0("name: '",name,"' column values must be unique"))
}
return(name)
}
getByLanguage <- function(varlist,language){
if(!is.null(language) && language[1] %in% names(varlist))
language <- language[1]
else
language <- "en"
return(unname(varlist[language]))
}
rbind.all.columns <- function(x, y) {
x.diff <- setdiff(colnames(x), colnames(y))
y.diff <- setdiff(colnames(y), colnames(x))
x[, c(as.character(y.diff))] <- NA
y[, c(as.character(x.diff))] <- NA
return(rbind(x, y))
}
pathParameter<-function(model,estimates=c("b","se","z","pvalue","beta")){
if(inherits(model,"lavaan")){
links<-lavaan::parameterEstimates(model,standardized = T)
names(links)<-gsub("^est$","b",names(links))
names(links)<-gsub("^std.all$","beta",names(links))
links<-links[links$op=="~",c("rhs","lhs",estimates)]
names(links)[1:2]<-c("Source","Target")
if(length(intersect(union(links$Source,links$Target),model@Data@ov$name))>0) {
nodes<-as.data.frame(model@Data@ov,stringsAsFactors=F)[,intersect(names(model@Data@ov),c("name","mean","var"))]
nodes$stdev<-sqrt(nodes$var)
row.names(nodes)<-nodes$name
}
else {
nodes<-data.frame(name=union(links$Source,links$Target))
}
# nodes<-nodes[,c("name","mean","stdev")]
nodes$name<-iconv(nodes$name,"","UTF-8")
links<-links[,c("Source","Target",estimates)]
links$Source<-iconv(links$Source,"","UTF-8")
links$Target<-iconv(links$Target,"","UTF-8")
structure(list(links = links, nodes = nodes))
}
else stop("Model has to be a lavaan object")
}
catFit<-function(model,fitMeasures){
if("chisq" %in% fitMeasures) fitMeasures<-union(fitMeasures,c("df","pvalue"))
fit<-lavaan::fitMeasures(model,fitMeasures)
text<-NULL
if("chisq" %in% names(fit)) text<-paste(text, paste0("Chi2=",format(fit["chisq"],digits=4)," (",fit["df"]," df)",", pvalue=", format(fit["pvalue"],digits=2)), sep=". ")
if("cfi" %in% names(fit)) text<-paste(text, paste0("CFI= ",format(fit["cfi"],digits=3)), sep=". ")
if("rmsea" %in% names(fit)) text<-paste(text, paste0("RMSEA= ",format(fit["rmsea"],digits=3)), sep=". ")
paste0(gsub("^. ","<p>",text),".</p>")
}
corr <- function (a, b = a, weight = NULL )
{
if (is.null(weight)) weight= rep(1/nrow(a), nrow(a))
s<-complete.cases(cbind(a,b,weight))
a<-as.matrix(a[s,]);b<-as.matrix(b[s,])
# normalize weights
weight <- weight[s] / sum(weight[s])
# center matrices
a <- sweep(a, 2, apply((a * weight),2,sum))
b <- sweep(b, 2, apply((b * weight),2,sum))
# compute weighted correlatio
t(a*weight) %*% b / sqrt(apply((a**2 *weight),2,sum) %*% t(apply((b**2 *weight),2,sum)))
}
corrp <- function (a, b = a, weight = NULL ){
count <- 0
CC <- matrix( NA , nrow=(ncol(a)), ncol=ncol(b))
rownames(CC) <- colnames(a)
colnames(CC) <- colnames(b)
for(i in colnames(a)) {
for(j in colnames(b)) {
if(i==j) CC[i,j] <- 1
else {
CC[i,j] <- corr(b[j],a[i],weight=weight)
}
}
count <- count +1
}
return(CC)
}
layoutMCA<-function(matrix) { # Correspondencias simples clasicas aplicadas a dicotomicas.
matrix<-cbind(matrix,1-matrix)
n<-sum(matrix)
P=matrix/n
column.masses<-colSums(P)
row.masses=rowSums(P)
E=row.masses %o% column.masses
R=P-E
I=R/E
Z=R/sqrt(E) # Corrected
SVD=svd(Z)
rownames(SVD$v)=colnames(P)
standard.coordinates.columns = sweep(SVD$v[1:(ncol(Z)/2),1:2], 1, sqrt(column.masses[1:(ncol(Z)/2)]), "/")
principal.coordinates.columns = sweep(standard.coordinates.columns, 2, SVD$d[1:2], "*")
colnames(principal.coordinates.columns)<-c("F1","F2")
return(principal.coordinates.columns)
}
layoutMca<-function(matrix, nfactors=2, rows=FALSE){ # Correspondencias simples clasicas aplicadas a dicotomicas.
P=matrix/nrow(matrix)
column.masses<-colSums(P)
row.masses=rowSums(P)
E=row.masses %o% column.masses
R=P-E
I=R/E
Z=R/sqrt(E) # Corrected
SVD=svd(Z)
rownames(SVD$v)=colnames(P)
CC <-list()
CC$standard.coordinates.rows = sweep(SVD$u, 1, sqrt(row.masses), "/")
CC$principal.coordinates.rows = sweep(CC$standard.coordinates.rows, 2, SVD$d, "*")
CC$standard.coordinates.columns = sweep(SVD$v, 1, sqrt(column.masses), "/")
CC$principal.coordinates.columns = sweep(CC$standard.coordinates.columns, 2, SVD$d, "*")
CC <- lapply(CC, function(X){X<-X[,2:(nfactors+1)];colnames(X) <- paste0("F",1:nfactors); return(X)})
if(rows) return(CC$principal.coordinates.rows)
else return(CC$principal.coordinates.columns)
}
layoutPCA<-function(coin) { # Coordenadas a partir de Pearson: Haberman/raiz(n)
if(inherits(coin, "coin")) A<-eigen(sim(coin,"P"))
else A <- eigen(coin)
C<-sweep(A$vectors[,1:2],2,sqrt(A$values[1:2]),"*")
rownames(C)<-rownames(coin)
colnames(C)<-c("F1","F2")
return(C)
}
# For d_netCorr Factorial coordinates.
layoutFact <- function(data, method=c("pc","nf","vf", "of")) {
difference <- setdiff(method, c("null", "pc","nf","vf","of"))
if(length(difference)>0) warning(paste0(difference, " method is not implemented. "))
R <- list()
if(method[1]=="null") return(NULL)
if("pc" %in% method) R$`Principal Components` <- princomp(na.omit(data), cor=T)$loadings[,1:2]
if("nf" %in% method) R$`Non-rotated Factorial` <- factanal(na.omit(data), factors=2, rotation="none")$loadings[,1:2]
if("vf" %in% method) R$`Varimax factorial` <- factanal(na.omit(data), factors=2, rotation="varimax")$loadings[,1:2]
if("of" %in% method) R$`Oblimin factorial` <- factanal(na.omit(data), factors=2, rotation="oblimin")$loadings[,1:2]
if(length(R)>1) return(R) else return(R[[1]])
}
# For d_netCorr: Labelling the graph components
compon <- function(comps) {
s <- table(comps)
o <- sort(table(comps), decreasing=T)
xdig <- paste0("%",max(c(nchar(trunc(abs(as.numeric(names(s))))))),"d")
j <- 1
for(n in as.numeric(names(o))) {
d <- as.numeric(names(s[j]))
comps <- ifelse(comps==n, paste0("C-",sprintf(xdig, d)), comps)
j <- j+1
}
return(comps)
}
mobileEdges<-function(data, name=1, number=2, difference=0) {
if(!is.numeric(data[[number]])) data[[number]]<-as.numeric(paste(data[[number]]))
DC<-matrix(NA,nrow=nrow(data),ncol=nrow(data))
colnames(DC)<-rownames(DC)<-data[[name]]
for(i in 1:nrow(data))DC[i,]=ifelse(abs(data[[number]][i]-t(data[[number]]))<=difference,(1+difference-abs(data[[number]][i]-t(data[[number]]))),0)
diag(DC)<-0
DCLinks<-edgeList(DC,"shape",min=1)
colnames(DCLinks)[3]<-"sim."
DCLinks$dist.<-(1+difference-DCLinks$sim.)
return(DCLinks)
}
incTime<-function(data, name="name", beginning="birth", end="death") {
if(!is.integer(data[[beginning]])) data[[beginning]]<-as.integer(paste(data[[beginning]]))
if(!is.integer(data[[end]])) data[[end]]<-as.integer(paste(data[[end]]))
anos<-min(na.omit(data[[beginning]])):max(na.omit(data[[end]]))
E<-matrix(NA,nrow=nrow(data),ncol=length(anos))
colnames(E)<-anos
for(i in 1:nrow(data)) E[i,]<-ifelse(anos>=data[[beginning]][i] & (anos<=data[[end]][i] | is.na(data[[end]][i])),1,0)
Datos<-as.data.frame(t(E))
colnames(Datos)<-data[[name]]
return(Datos)
}
shinyCoin <- function(x){
shiny_rd3(x)
}
# Atencion a lwidth y en lenguaje espanol o catalan
glmCoin <- function(formulas, data, weights=NULL, pmax=.05, twotail=FALSE, showArrows=TRUE,
frequency = FALSE, percentage = TRUE,
color="variable", lwidth="z.value", circle= NA, language=c("en","es","ca"),
igraph=FALSE, ...){
if (is.character(weights)) weights<-data[[weights]]
Links <- data.frame(A=NA,B=NA,C=NA,D=NA, E=NA, G=NA, H=NA)[-1,]
names(Links) <- headreg(language)
if (lwidth=="z.value" & language[1]!="en") lwidth<-"val.z"
Formulas <- as.list(gsub("[[:space:]]","",unlist(strsplit(formulas,"\n"))))
formulas <- sapply(Formulas,function(X){substr(X,start=1,stop=as.numeric(gregexpr(",",X))-1)})
familias <- sapply(Formulas,family)
# dependent <- sapply(Formulas,function(X){substr(X,start=1,stop=as.numeric(gregexpr("\u7E",X))-1)})
# dependent <- gsub("\u60","",dependent)
variables<-extract(formulas, data)
for(instance in 1:length(formulas)) {
m <- net_lm(formulas[instance], familias[instance], data, weights, pmax, twotail)
if (nrow(m)>0) {
names(m) <- headreg(language)[-7]
m[[headreg(language)[7]]] <- gsub("`","",formulas[instance])
Links <-rbind(Links,m)
}
}
if (!("nodes" %in% names(list(...)))) {
Nodes<-data.frame(name=iconv(union(Links$Source,Links$Target),to="UTF-8"),
variable=gsub(":.*","",iconv(union(Links$Source,Links$Target),to="UTF-8"))
,stringsAsFactors = FALSE)
row.names(Nodes)<-Nodes$name
arguments<-list(nodes=Nodes, links=Links, showArrows = showArrows, color = color, lwidth = lwidth, language = language, ...)
}
else arguments<-list(links=Links, showArrows = showArrows, color = color, lwidth = lwidth, language = language, ...)
if (frequency | percentage) arguments$nodes<-meanPer(data, variables, arguments$nodes, names(arguments$nodes)[1] , frequency, percentage, weights)
#ADD to N percentages/means. vid extract.R and as.nodes(surCoin)
if (!"name" %in% names(arguments)) {
arguments$name <- nameByLanguage(arguments$name,language,arguments$nodes)
names(arguments$nodes)[1]<-nameList[language[1]]
}
if (!is.na(circle)) arguments$layout <- layoutCircle(arguments$nodes, variables$D, circle)
if(twotail){
arguments$linkBipolar <- TRUE
arguments$lcolor <- "Estimate"
}
if(nrow(arguments$nodes)+nrow(arguments$links)>0) {
xNx <- do.call(netCoin,arguments)
}
else(stop("No nodes, no relations"))
# xNx$nodes$name<-iconv(xNx$nodes$name,to="UTF-8")
# xNx$links[[c("Source","Target")]]<-iconv(xNx$links[[c("Source","Target")]],to="UTF-8")
if (igraph) return(toIgraph(xNx))
else return(xNx)
}
net_lm<-function(formula, family=gaussian, data, weights=NULL, pmax=.05, twotail=FALSE, ...){
arguments <- list(formula=formula, family=family, data= data, weights= weights, ...)
glm <- do.call(glm, arguments)
links <- linkregress(remodel(glm),pmax=pmax,twotails=twotail)
return(links)
}
remodel<-function(model){
DD <- model$model
FO <- model$formula
FF <- model$family$family
if (!is.null(model$model$`(weights)`)) WW <- model$model$`(weights)`
else WW<-NULL
C <- model$coefficients
## VV<-row.names(attr(terms(F),"factors"))[-1]
VV <- gsub("[[:space:]]","",unlist(strsplit(strsplit(FO,"\\~")[[1]][2],"\\+")))
V <- VV[vapply(DD[VV],inherits,TRUE,what="factor")]
if (length(V)==0) return(model)
Vq <- paste0("^",V)
y <- lapply(Vq,grep,names(model$coefficients))
names(y) <- Vq
if (is.null(names(y))) colnames(y)<-gsub("\\^","",colnames(y))
else names(y)<-gsub("\\^","",names(y))
c<-f<-m<-z<-list()
for (i in names(y)) {
z[i]<-list(C[y[[i]]])
m[i]<-which.min(unlist(z[i]))
f[i]<-ifelse(z[[i]][m[[i]]]<0,m[[i]],0)
if(f[i]>0) DD[[i]] <- relevel(DD[[i]],sub(names(z[i]),"",names(z[[i]])[f[[i]]]))
c[[i]]<-matriz(DD[[i]])
}
if (is.null(WW)) neomodel<-glm(FO, contrasts= c, data=DD, family=FF)
else neomodel <- glm(FO, contrasts= c, data=cbind(DD,WW), family=FF, weights = WW)
return(neomodel)
}
linkregress<-function(model, pmax=0.05, twotails=FALSE){
pmax<-abs(ifelse(twotails,qnorm(pmax/2,FALSE),qnorm(pmax,FALSE)))
q<-summary(model)$coefficients
y<-names(model$xlevels) #names(model$model)[-1]
for (i in y) {
rownames(q)<-gsub(paste0("(^",i,")"),"\\1:",rownames(q))
}
m<-data.frame(Source=rownames(q),Target=names(model$model)[1], q, stringsAsFactors = FALSE)
rownames(m)<-NULL
names(m)[-c(1:2)]<-colnames(q)
m<-m[-1,]
crit<-names(m)[grep(" value",names(m))]
if (twotails) m <- m[abs(m[[crit]])>pmax,]
else m <- m[m[[crit]]>pmax,]
row.names(m)<-NULL
return(m)
}
matriz<-function(factor){
T<-t(table(factor)/length(factor))
M<-matrix(rep(T,length(T)),nrow=length(T),
dimnames=list(colnames(T),colnames(T)))
M[,2:dim(M)[1]]<--M[,2:dim(M)[1]]
diag(M)[2:dim(M)[1]]<-1-M[2:dim(M)[1],1]
m<-t(matrix(M[,2:dim(M)[2]],nrow=dim(M)[2],dimnames=list(colnames(T),colnames(T)[-1])))
m.t<-rbind(constant=1/ncol(m),m)
matriz<-matrix(solve(m.t)[,-1], nrow=nrow(m.t), dimnames=list(colnames(m.t),rownames(m.t)[-1]))
return(matriz)
}
family<-function(formula) {
FAM <- ifelse(as.numeric(gregexpr(",",formula))==-1,
"",
substr(formula,start=as.numeric(gregexpr(",",formula))+1,1E6))
family<-c.family(FAM)
return(family)
}
c.family<-function(method=NULL) {
if (is.null(method) | method=="") family<-"GAU"
else {
family<-substr(toupper(method),1,3)
if (family=="QUA" & nchar(method)>5) family<-substr(toupper(method),1,6)
}
families<-matrix(c("gaussian","binomial","Gamma","inverse.gaussian","poisson","quasi","quasibinomial","quasipoisson"),
nrow=1, dimnames=list("Family", c("GAU","BIN","GAM","INV","POI","QUA","QUASIB","QUASIP")))
return(families[,family])
}
headreg<-function(language){
switch (language[1],
en = c('Source','Target','Estimate','Std.error', 'z.value', 'Pr(>|z|)','Equation'),
es = c('Source','Target','Estimador','Err.t\uEDp.', 'val.z', 'Pr(>|z|)','Ecuaci\uF3n'),
ca = c('Source','Target','Estimador','Err.t\uEDp.', 'val.z', 'Pr(>|z|)', 'Ecuaci\uF3')
)
}
extract <- function(formulas, data) {
formulas<-gsub("`","",formulas)
dependent<-variables<-factors<-NULL
for (formula in formulas) {
dependent <- union(dependent, gsub("[[:space:]]","",unlist(strsplit(strsplit(formula,"\\~")[[1]][1],"\\+"))))
variables <- union(variables, gsub("[[:space:]]","",unlist(strsplit(strsplit(formula,"\\~")[[1]][2],"\\+"))))
factors <- union(factors, variables[vapply(data[variables],inherits,TRUE,what="factor")])
}
independent<-setdiff(variables,factors)
variables<-list(D=dependent, I=independent, F=factors)
return(variables)
}
meanPer<-function(data, variables, frame, name=names(frame[1]), frequency= FALSE, percentage= TRUE, weights = NULL){
if (is.null(weights)) weights <- rep(1, nrow(data))
columns<-setdiff(names(frame),c(name, "n.","%"))
l.frame<-length(frame)
row.names(frame)<-frame[[name]]
frame.order<-row.names(frame)
quantitatives<-c(variables$D,variables$I)
data<-na.omit(cbind(as.data.frame(data)[,unlist(variables)],weights))
for (varF in variables$F) {
beginF <- length(data)+1
data<-dichotomize(data, varF, "")
names(data)[beginF:length(data)] <- paste0(varF,":",names(data)[beginF:length(data)])
}
setvariables <- setdiff(names(data),union(weights,variables$F))
sta<-data.frame(names=setvariables)
row.names(sta)<-sta$names
sta$N. <- 0
if (frequency & !is.null(weights)) sta$N. <-round(apply(data[,setvariables]*data$weights, 2, sum),0)
sta$N. <- ifelse(sta$names %in% quantitatives, nrow(data), sta$N.)
if (percentage) {
means<-apply(data[,setvariables], 2, weighted.mean, data$weights)
maxs <-apply(data[,setvariables], 2, max)
mins <-apply(data[,setvariables], 2, min)
sta$M.<-(means-mins)/(maxs-mins)*100
}
frame<-merge(frame,sta[,-1,drop=F], by="row.names", all.x = TRUE)[,-1]
row.names(frame)<-frame[[name]]
adds<-c()
if (frequency) {
frame$n. <- frame$N.
adds <- "n."
}
if (percentage) {
frame$`%`<- frame$M.
adds <- c(adds, "%")
}
out<- grep('^[KLMN]\\.$',names(frame))
frame<-frame[,-c(out)]
return(frame[frame.order,c(name,adds,columns)])
}
logCoin<-function(data, variables=names(data), exogenous=NULL, noFirstCat=NULL, weight=NULL,
order=2, pairwise=FALSE, twotails=FALSE, pmax=.05,
frequency=FALSE, percentage=FALSE,
directed=FALSE, igraph=FALSE, ...) {
arguments <- list(...)
names(data) <- gsub(" ","_", names(data))
variables <- gsub(" ","_", variables)
if(!is.null(exogenous)) exogenous <- gsub(" ","_", exogenous)
if (!is.null(noFirstCat)) noFirstCat <- gsub(" ","_", noFirstCat)
if(!is.null(weight)) weight <- gsub(" ","_", weight)
variables <- union(setdiff(union(variables, noFirstCat), exogenous), exogenous)
if(!is.null(weight)) {
if(inherits(weight,"character")){
variables <- setdiff(variables, weight)
weight<-data.frame(weight=data[[weight]])
data<-cbind(weight,data[,variables])
}
else{
if(length(weight)!=dim(data)[1]) stop("Weights have not the correct dimensions")
if (pairwise) data <- cbind(data[,variables],weight)[,1:length(data[,variables])]
else data <- na.omit(cbind(data[,variables],weight))[,1:length(data[,variables])]
}
}
else data<-cbind(data[,variables], data.frame(weight=rep(1, nrow(data))))
formula <- comb(setdiff(variables,exogenous), exogenous, "weight", order)
varOrder <- variables # To order variables later before coin
#Check methods. No necessary because edgeList call these routines.
#procedures<-i.method(c.method(procedures))
#criteria<-i.method(c.method(criteria))
#procedures<-union(procedures,unlist(arguments[c("lwidth","lweight","lcolor","ltext")]))
#Names
if(!("language" %in% names(arguments))) arguments$language <- "en"
nodes <- arguments$nodes
if (inherits(nodes,"tbl_df")) nodes<-as.data.frame(nodes)
name <- arguments$name <- nameByLanguage(arguments$name,arguments$language,arguments$nodes)
if (!("level" %in% names(arguments))) level<-.95 else level <-arguments$level
#Data.frame
if (all(inherits(data,c("tbl_df","tbl","data.frame"),TRUE))) data<-as.data.frame(data) # convert haven objects
if (inherits(weight,"character")) variables <- setdiff(variables,weight)
pivots <- setdiff(variables, noFirstCat)
if (length(pivots)>5) stop("This function doesn't support more than 5 variables with first categories included")
if (!pairwise & inherits(weight,"character")) {
if (!is.null(weight)) weight <- data[rowSums(is.na(data[,variables]))<1,weight]
data <- data[complete.cases(data[,variables]),]
}
arguments$scenarios <- sum(rowSums(!is.na(data))>0)
data <- as_factor(data)
dt <- aggregate(formula, data=data, FUN="sum")
dt$weight <- round(dt$weight) # To estimate Poisson (without decimals)
sc <- sum(rowSums(!is.na(dt))>0)
dtab <- xtabs(weight ~ ., dt)
fm <- loglm(comb(setdiff(variables,exogenous), exogenous, "", order), dtab) # numerals as names.
nam <- eti <- var <- cat <- varn <- labn <- NULL
nvar <-1
ncell <- 1
for (x in variables) {
dt[[x]] <- as_factor(dt[[x]])
nval <- length(levels(dt[[x]]))
nam <- c(nam, paste0(x, levels(dt[[x]])))
eti <- c(eti, paste0(x, ":", levels(dt[[x]])))
var <- c(var, rep(x, nval))
cat <- c(cat, levels(dt[[x]]))
varn <- c(varn, rep(nvar, nval))
labn <- c(labn, 1:nval)
contrasts(dt[[x]]) <- contr.first(dt[[x]])
nvar <- nvar+1
ncell <- ncell*nval
}
informa <- data.frame(eti, var, cat, varn, labn, row.names = nam)
arguments$note <- paste0("<p>L2= ", sprintf("%3.2f",fm$lrt), "; d.f.= ", fm$df, "; p(>X^2)= ",
sprintf("%3.3g",1-pchisq(fm$lrt, fm$df)), ".</p><p>Covariance structures: ",
sprintf("%1.0f", sc)," (",sprintf("%1.0f", ncell),").</p><p>",
Reduce(paste, deparse(formula)),"</p>")
coefs <- summary(glm(formula, family=poisson, data=dt))$coefficients
if (length(pivots)==0) coefs <- list(data=dt, coefs=coefs)
if (length(pivots) >0) {
coefs <- addCoefs(c(1,1), coefs, formula, dt, pivots)
}
if (length(pivots) >1) {
coefs <- addCoefs(c(1,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
}
if (length(pivots) >2) {
coefs <- addCoefs(c(1,3), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
}
if (length(pivots) >3) {
coefs <- addCoefs(c(1,4), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,3), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
}
if (length(pivots) >4) {
coefs <- addCoefs(c(1,5), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,3), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,4), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,3), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(1,1), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,2), coefs$coef, formula, coefs$data, pivots)
coefs <- addCoefs(c(0,1), coefs$coef, formula, coefs$data, pivots)
}
coefs <- as.data.frame(coefs$coef[sort(row.names(coefs$coef)),])
nodes <- coefs[!grepl(":",rownames(coefs)), ][-1,1, drop=F]
# nodes
nodes$name <- row.names(nodes)
nodes$Estimate <- round(nodes$Estimate, 3)
nodes <- merge(nodes, informa [1:3], by="row.names")[,-1]
if(frequency | percentage) nodes <- merge(nodes, margins(dtab, frequency, percentage), by="name", all.x=T)
ordervars <- factor(nodes$var, levels=varOrder)
ordernodes <- nodes[order(ordervars), "name"]
# coefs
pmax<-abs(ifelse(twotails,qnorm(pmax/2,FALSE),qnorm(pmax,FALSE)))
if(!twotails) {
coefs <- coefs[coefs$`z value`>pmax & grepl(":",rownames(coefs)), ]
coefs$`Pr(>|z|)` <- coefs$`Pr(>|z|)` /2
names(coefs)[names(coefs)=="Pr(>|z|)"] <- "Pr(>z)"
}
else {
coefs <- coefs[abs(coefs$`z value`)>pmax & grepl(":",rownames(coefs)), ]
arguments$linkBipolar <- TRUE
arguments$lcolor <- "Estimate"
}
coefs$order <- nchar(gsub("[^:]","",row.names(coefs)))+1
links <- coefs[coefs$order==2,]
links$Source <- sub(":.*","", row.names(links))
links$Target <- sub(".*:","", row.names(links))
if(nrow(links)){
links[,c("Source", "Target")] <- olinks(links$Source, links$Target, ordernodes)
}
ulinks <- coefs[coefs$order>2,]
if (nrow(ulinks)>0) {
L <- unlist(strsplit(row.names(ulinks),":"))
CC <- paste0("u", sprintf(paste0("%0",ceiling((nrow(ulinks)+1)/10),"d"),1:nrow(ulinks)))
ulinks$Source <- CC
M <- data.frame(Source=rep(CC, times=ulinks$order), Target=L)
unodes <- data.frame(name=CC, Estimate=floor(min(nodes$Estimate)), eti="", var=paste0(".I", ulinks$order), cat="")
if(frequency) unodes$freq <- sapply(unodes$name, FUN=ofreq, x=M, table=dtab, index=informa)
if(percentage) unodes$prop <- sapply(unodes$name, FUN=ofreq, x=M, table=dtab, index=informa)/marginSums(dtab)
nodes <- rbind(nodes, unodes)
ulinks <- merge(M, ulinks)
links <- rbind(links, ulinks)
}
if(exists("language", arguments)) {
names(nodes)[names(nodes)=="name"] <- getByLanguage(nameList, arguments$language)
names(nodes)[names(nodes)=="eti"] <- getByLanguage(labelList, arguments$language)
}
arguments$nodes <- nodes
if(nrow(links)){
arguments$links <- links
if (!"degreeFilter" %in% names(arguments)) arguments$degreeFilter <- 1
}
if (!"size" %in% arguments) {
arguments$size= "Estimate"
if(frequency) arguments$size <- "freq"
if(percentage) arguments$size <- "prop"
}
if (!"color" %in% names(arguments)) arguments$color <- "var"
if (!"lwidth" %in% names(arguments)) arguments$lwidth <- "Estimate"
if (!"label" %in% names(arguments)) arguments$label <- getByLanguage(labelList, arguments$language)
G <- do.call(netCoin, arguments)
return(G)
}
# Funciones previas
contr.first <- function(variable){
if(!inherits(variable, "factor")) stop
nlevels <- length(levels(variable))
f <- matrix(c(rep(-1, nlevels-1),diag(nlevels-1)), nrow=nlevels, byrow=TRUE, dimnames=list(levels(variable),levels(variable)[-1]))
return(f)
}
contr.last <- function(variable) {
if(!inherits(variable, "factor")) stop
nlevels <- length(levels(variable))
l <- matrix(c(diag(nlevels-1), rep(-1, nlevels-1)), nrow=nlevels, byrow=TRUE, dimnames=list(levels(variable),levels(variable)[-nlevels]))
return(l)
}
addCoefs <- function(vector=c(0,1), coefficients, formula, data, variables, family="poisson") {
if(vector[1]==1) contrasts(data[[variables[[vector[2]]]]]) <- contr.last(data[[variables[[vector[2]]]]])
else contrasts(data[[variables[[vector[2]]]]]) <- contr.first(data[[variables[[vector[2]]]]])
coefs <- rbind(coefficients,
summary(glm(formula, family=family, data=data))$coefficients)
return(list(data=data, coef=coefs[unique(rownames(coefs)),]))
}
comb <-function (endogenous, exogenous=NULL, frequency="x", order=2) {
endogenous <- setdiff(endogenous, exogenous)
combin <- function(x, order) {apply(combn(x, order), 2, paste0, collapse="*")}
formula <- paste0(frequency,"~")
ordEndo <- min(order, length(endogenous))
l <- combin(endogenous, ordEndo)
if (is.null(exogenous) | length(endogenous)>=order ) for (nl in l) formula <- paste0(formula,nl,"+")
if (!is.null(exogenous)) {
for (exo in exogenous) {
if(order<length(endogenous)+length(exogenous)) l<- combin(endogenous, min(ordEndo,order-1))
nexo <- ifelse(l[1]=="","","*")
for (nl in l) formula <- paste0(formula,nl,nexo,exo, "+")
}
}
return(as.formula(sub("\\+$","",formula)))
}
margins <- function(table, freq=T, prop=F) {
names <- names(attr(table,"dimnames"))
vect <- c()
n <- ifelse(prop, marginSums(table),1)
for (x in names) {
marg <- marginSums(table, x)
names(marg) <- paste0(x, names(marg))
vect <- c(vect, marg)
}
dtfrm <- data.frame(name=names(vect), freq=vect)
if (prop) dtfrm$prop <- dtfrm$freq/n
if (!freq) dtfrm$freq <- NULL
return(dtfrm)
}
olinks <- function(source, target, order) {
for (i in 1:length(source)) {
if (which(source[[i]]==order) < which(target[[i]]==order)) {
s <- target[[i]]
target[[i]] <- source[[i]]
source[[i]] <- s
}
}
return(cbind(source, target))
}
ofreq <- function (u, x, table, index) {
vector <- x[x$Source==u, "Target"]
return(marginSums(table, index[vector, 4])[t(index[vector,5])])
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.