R/gl.plot.network.r

In dartR: Importing and Analysing 'SNP' and 'Silicodart' Data Generated by Genome-Wide Restriction Fragment Analysis

#' Represents a distance or dissimilarity matrix as a network
#'
#' This script takes a distance matrix generated by dist() and represents the
#' relationship among the specimens as a network diagram. In order to use this
#' script, a decision is required on a threshold for relatedness to be
#' represented as link in the network, and on the layout used to create the
#' diagram.
#'
#' The threshold for relatedness to be represented as a link in the network is
#' specified as a quantile. Those relatedness measures above the quantile are
#' plotted as links, those below the quantile are not. Often you are looking for
#' relatedness outliers in comparison with the overall relatedness among
#' individuals, so a very conservative quantile is used (e.g. 0.004), but
#' ultimately, this decision is made as a matter of trial and error. One way to
#' approach this trial and error is to try to achieve a sparse set of links
#' between unrelated 'background' individuals so that the stronger links are
#' preferentially shown.
#'
#' There are several layouts from which to choose. The most popular are given as
#' options in this script.
#' \itemize{
#' \item fr -- Fruchterman, T.M.J. and Reingold, E.M. (1991). Graph Drawing by
#' Force-directed Placement. Software -- Practice and Experience 21:1129-1164.
#' \item kk -- Kamada, T. and Kawai, S.: An Algorithm for Drawing General
#' Undirected Graphs. Information Processing Letters 31:7-15, 1989.
#' \item drl -- Martin, S., Brown, W.M., Klavans, R., Boyack, K.W., DrL:
#' Distributed Recursive (Graph) Layout. SAND Reports 2936:1-10, 2008.
#' }
#'
#' Colors of node symbols are those of the rainbow.
#'
#'@param D A distance or dissimilarity matrix generated by dist() or gl.dist()
#' [required].
#'@param x A genlight object from which the D matrix was generated
#'[default NULL].
#'@param method One of "fr", "kk" or "drl" [default "fr"].
#'@param node.size Size of the symbols for the network nodes [default 3].
#'@param node.label TRUE to display node labels [default FALSE].
#'@param node.label.size Size of the node labels [default 0.7].
#'@param node.label.color Color of the text of the node labels
#'[default 'black'].
#'@param alpha Upper threshold to determine which links between nodes to display
#'[default 0.005].
#'@param title Title for the plot
#'[default "Network based on genetic distance"].
#' @param verbose Verbosity: 0, silent or fatal errors; 1, begin and end; 2,
#' progress log; 3, progress and results summary; 5, full report
#' [default 2, unless specified using gl.set.verbosity].
#'@return returns no value (i.e. NULL)
#'@importFrom grDevices rgb
#'@importFrom graphics legend
#'@export
#' @author Custodian: Arthur Georges -- Post to
#' \url{https://groups.google.com/d/forum/dartr}
#'@examples
#' if ((requireNamespace("rrBLUP", quietly = TRUE)) & (requireNamespace("gplots", quietly = TRUE))) {
#' test <- gl.subsample.loci(platypus.gl, n = 100)
#' test <- gl.keep.ind(test,ind.list = indNames(test)[1:10])
#' D <- gl.grm(test, legendx=0.04)
#' gl.plot.network(D,test)
#' }

gl.plot.network <- function(D,
                            x = NULL,
                            method = "fr",
                            node.size = 3,
                            node.label = FALSE,
                            node.label.size = 0.7,
                            node.label.color = "black",
                            alpha = 0.005,
                            title = "Network based on genetic distance",
                            verbose = NULL) {
    # CHECK IF PACKAGES ARE INSTALLED
    pkg <- "igraph"
    if (!(requireNamespace(pkg, quietly = TRUE))) {
      cat(error(
        "Package",
        pkg,
        " needed for this function to work. Please install it.\n"
      ))
      return(-1)
    }
    
    # SET VERBOSITY
    verbose <- gl.check.verbosity(verbose)
    
    # FLAG SCRIPT START
    funname <- match.call()[[1]]
    utils.flag.start(func = funname,
                     build = "Jody",
                     verbosity = verbose)
    
    # CHECK DATATYPE
    datatype <- utils.check.datatype(x, verbose = verbose)
    
    # FUNCTION SPECIFIC ERROR CHECKING
    
    if (!is(D, "dist") & !is(D, "matrix")) {
        stop(error(
            "Fatal Error: distance matrix required for gl.dist.network!\n"
        ))
    }
    
    if (!is.null(x)) {
        if (!is(x, "genlight")) {
            stop(
                error(
                    "Fatal Error: if specified, genlight object required for gl.dist.network!\n"
                )
            )
        }
    } else {
        if (verbose >= 2) {
            cat(
                warn(
                    "Note: genlight object not specified, population assignments not available for plotting\n"
                )
            )
        }
    }
    
    if (!(method == "fr" ||
          method == "kk" || method == "drl")) {
        cat(warn(
            "Warning: Layout method must be one of fr, or kk, or drl, set to fr\n"
        ))
        method <- "fr"
    }
    
    # DO THE JOB
    
    m <- as.matrix(D)
    len <- length(m[, 1])
    links <-
        data.frame(array(NA, dim = c((len * len - len) / 2, 3)))
    count <- 1
    for (i in 1:(len - 1)) {
        for (j in (i + 1):len) {
            links[count, 1] <- row.names(m)[i]
            links[count, 2] <- row.names(m)[j]
            links[count, 3] <- m[i, j]
            count <- count + 1
        }
    }
    colnames(links) <- c("from", "to", "weight")
    
    if (!is.null(x)) {
        nodes <- data.frame(cbind(x$ind.names, as.character(pop(x))))
        colnames(nodes) <- c("name", "pop")
    } else {
        nodes <- data.frame(row.names(m))
        colnames(nodes) <- "name"
    }
    
    network <-
        igraph::graph_from_data_frame(d = links,
                                      vertices = nodes,
                                      directed = FALSE)
    
    if (!is.null(x)) {
        colors <-rainbow(nlevels(pop(x)))
        my_colors <- colors[pop(x)]
    } else {
        my_colors <- "red"
    }
    
    q <- stats::quantile(links$weight, p = 1 - alpha)
    network.FS <-
        igraph::delete_edges(network, igraph::E(network)[links$weight < q])
    
    if (method == "fr") {
        layout.name <- "Fruchterman-Reingold layout"
        l <- igraph::layout_with_fr(network.FS)
    }
    if (method == "kk") {
        layout.name <- "Kamada-Kawai layout"
        l <- igraph::layout_with_kk(network.FS)
    }
    if (method == "drl") {
        layout.name <- "DrL Graph layout"
        l <- igraph::layout_with_drl(network.FS)
    }
    title <- paste(title, "\n[", layout.name, "]")
    
    if (node.label) {
        node.label <- igraph::V(network)$name
    } else {
        node.label <- NA
        node.label.size <- NA
        node.label.color <- NA
    }
    
    plot(
        network.FS,
        edge.arrow.size = 0,
        edge.curved = 0,
        edge.width = links$weight,
        vertex.size = node.size,
        vertex.color = my_colors,
        vertex.frame.color = "#555555",
        vertex.label = node.label,
        vertex.label.color = node.label.color,
        vertex.label.cex = node.label.size,
        layout = l,
        main = title
    )
    
    if (!is.null(x)) {
        legend(
            "bottomleft",
            legend = levels(pop(x)),
            col = colors,
            bty = "n",
            pch = 20,
            pt.cex = 3,
            cex = 1,
            text.col = colors,
            horiz = FALSE,
            inset = c(0.1, 0.1)
        )
    }
    
    # FLAG SCRIPT END
    
    if (verbose > 0) {
        cat(report(report("Completed:", funname, "\n")))
    }
    
    return(invisible())
    
}