R/export_map_list.R

Defines functions export_map_list

Documented in export_map_list

#' Export a genetic map to a CSV file
#' 
#' Function to export genetic linkage map(s) generated by \code{MAPpoly}.
#' The map(s) should be passed as a single object or a list of objects of class \code{mappoly.map}.
#'
#' @param  map.list A list of objects or a single object of class \code{mappoly.map}
#' 
#' @param file either a character string naming a file or a connection open for writing. 
#'    "" indicates output to the console.
#' 
#' @examples
#'  export_map_list(solcap.err.map[[1]], file = "")
#'
#' @author Marcelo Mollinari, \email{mmollin@ncsu.edu}
#'
#' @references
#'     Mollinari, M., and Garcia, A.  A. F. (2019) Linkage
#'     analysis and haplotype phasing in experimental autopolyploid
#'     populations with high ploidy level using hidden Markov
#'     models, _G3: Genes, Genomes, Genetics_. 
#'     \doi{10.1534/g3.119.400378}
#'
#' @export export_map_list
#' @importFrom dplyr tibble
export_map_list <- function(map.list, file = "map_output.csv"){
  if (inherits(map.list, "mappoly.map")) 
    map.list <- list(map.list)
  if(is.list(map.list)){
    if (any(!sapply(map.list, inherits, "mappoly.map"))) 
      stop("All elemnts in 'map.list' should be of class 'mappoly.map'")
  }
  R <- NULL
  for(i in 1:length(map.list)){
    ploidy <- map.list[[i]]$info$ploidy
    ph.P <- ph_list_to_matrix(map.list[[i]]$maps[[1]]$seq.ph$P, ploidy)
    ph.Q <- ph_list_to_matrix(map.list[[i]]$maps[[1]]$seq.ph$Q, ploidy)
    if(length(map.list[[i]]$info$seq.ref) == nrow(ph.P))
      for(j in 1:nrow(ph.P)){
        ph.P[j, ph.P[j, ] == 0] <- map.list[[i]]$info$seq.ref[j]
        ph.P[j, ph.P[j, ] == 1] <- map.list[[i]]$info$seq.alt[j]
        ph.Q[j, ph.Q[j, ] == 0] <- map.list[[i]]$info$seq.ref[j]
        ph.Q[j, ph.Q[j, ] == 1] <- map.list[[i]]$info$seq.alt[j]
      }
    ph.P <- as.data.frame(ph.P)
    ph.Q <- as.data.frame(ph.Q)
    colnames(ph.P) <- letters[1:ploidy]
    colnames(ph.Q) <- letters[(1+ploidy):(2*ploidy)]
    if(is.null(map.list[[i]]$info$chrom))
      map.list[[i]]$info$chrom <- rep(NA, nrow(ph.P))
    if(is.null(map.list[[i]]$info$genome.pos))
      map.list[[i]]$info$genome.pos <- rep(NA, nrow(ph.P))
    if(is.null(map.list[[i]]$info$seq.ref))
      map.list[[i]]$info$seq.ref <- rep(NA, nrow(ph.P))
    if(is.null(map.list[[i]]$info$seq.alt))
      map.list[[i]]$info$seq.alt <- rep(NA, nrow(ph.P))
    x <- dplyr::tibble("Marker Name" = map.list[[i]]$info$mrk.names,
                    "LG" = rep(i, nrow(ph.P)),
                    "Ref Chrom" = map.list[[i]]$info$chrom,
                    "Ref Position" = map.list[[i]]$info$genome.pos,
                    "Ref Allele" = map.list[[i]]$info$seq.ref,
                    "Alt Allele" = map.list[[i]]$info$seq.alt,
                    "Map Position" = round(cumsum(c(0, imf_h(map.list[[i]]$maps[[1]]$seq.rf))),2),
                    "Dose in P1" = map.list[[i]]$info$seq.dose.p1,
                    "Dose in P2" = map.list[[i]]$info$seq.dose.p2,
                    ph.P, ph.Q)
    R <- rbind(R, x)
  }
  write.csv(R , file = file, row.names = FALSE)
  invisible(R)
}
mmollina/MAPPoly documentation built on March 8, 2024, 2:04 a.m.