R/classes_generics_methods.R

In polysat: Tools for Polyploid Microsatellite Analysis

## load library just to make things work properly for now
#library(polysat)

## Set classes for containing ploidies
setClass("ploidysuper", contains="VIRTUAL")
setClass("ploidymatrix", representation(pld="matrix"), contains="ploidysuper",
         validity=function(object){
           failures <- character(0)

           if(is.null(dimnames(object@pld)[[1]]))
             failures <- c(failures, "Sample names missing.")
           if(is.null(dimnames(object@pld)[[2]]))
             failrues <- c(failures, "Locus names misisng.")

           if(length(failures)==0){
             return(TRUE)
           } else {
             return(failures)
           }
           })
setClass("ploidysample", representation(pld="integer"), contains="ploidysuper",
         validity=function(object){
           failures <- character(0)

           if(is.null(names(object@pld)))
              failures <- c(failures, "Sample names missing.")

           if(length(failures==0)){
             return(TRUE)
           } else {
             return(failures)
           }
           })
setClass("ploidylocus", representation(pld="integer"), contains="ploidysuper",
         validity=function(object){
           failures <- character(0)

           if(is.null(names(object@pld)))
             failures <- c(failures, "Locus names missing.")

           if(length(failures==0)){
             return(TRUE)
           } else {
             return(failures)
           }
           })
setClass("ploidyone", representation(pld="integer"), contains="ploidysuper",
         validity=function(object){
           failures <- character(0)

           if(!is.null(names(object@pld)))
             failures <- c(failures, "Ploidy vector should not be named.")
           if(length(object@pld)!=1)
             failures <- c(failures, "Ploidy vector should have length 1.")

           if(length(failures==0)){
             return(TRUE)
           } else {
             return(failures)
           }
         })
setMethod("initialize", signature(.Object="ploidymatrix"),
          function(.Object, samples, loci, ...){
            .Object@pld <- matrix(as.integer(NA), nrow=length(samples),
                                 ncol=length(loci), dimnames=list(samples, loci))
            callNextMethod(.Object, ...)})
setMethod("initialize", signature(.Object="ploidysample"),
          function(.Object, samples, loci, ...){
            ploidies <- as.integer(rep(NA, length(samples)))
            names(ploidies) <- samples
            .Object@pld <- ploidies
            callNextMethod(.Object, ...)})
setMethod("initialize", signature(.Object="ploidylocus"),
          function(.Object, samples, loci, ...){
            ploidies <- as.integer(rep(NA, length(loci)))
            names(ploidies) <- loci
            .Object@pld <- ploidies
            callNextMethod(.Object, ...)})
setMethod("initialize", signature(.Object="ploidyone"),
          function(.Object, samples, loci, ...){
            .Object@pld <- as.integer(NA)
            callNextMethod(.Object, ...)})
# generics to get and replace ploidies
setGeneric("pld", function(object, samples, loci) standardGeneric("pld"))
setGeneric("pld<-", function(object, value) standardGeneric("pld<-"))
# methods to get and replace ploidies
setMethod("pld", signature(object = "ploidyone"),
          function(object,samples,loci) return(object@pld))
setMethod("pld", signature(object = "ploidymatrix"),
          function(object,samples,loci){
            if(missing(samples)) samples <- dimnames(object@pld)[[1]]
            if(missing(loci)) loci <- dimnames(object@pld)[[2]]
            return(object@pld[samples,loci, drop=FALSE])
          })
setMethod("pld", signature(object = "ploidysample"),
          function(object,samples,loci){
            if(missing(samples)) samples <- names(object@pld)
            return(object@pld[samples])
          })
setMethod("pld", signature(object = "ploidylocus"),
          function(object,samples,loci){
            if(missing(loci)) loci <- names(object@pld)
            return(object@pld[loci])
          })
setReplaceMethod("pld", "ploidyone", function(object, value){
  if(length(value) > 1) stop("Only one ploidy allowed for this format.")
  object@pld <- unname(as.integer(value))
  object
})
setReplaceMethod("pld", "ploidymatrix", function(object, value){
  dn <- dimnames(object@pld)
  object@pld[dn[[1]],dn[[2]]] <- as.integer(value)
  object
})
setReplaceMethod("pld", "ploidysample", function(object, value){
  object@pld[names(object@pld)] <- as.integer(value)
  object
})
setReplaceMethod("pld", "ploidylocus", function(object, value){
  object@pld[names(object@pld)] <- as.integer(value)
  object
})
# generic and methods to test collapsibility of ploidies
setGeneric("plCollapse",
           function(object, na.rm, returnvalue) standardGeneric("plCollapse"))
setMethod("plCollapse", signature(object="ploidymatrix", na.rm="logical",
                                  returnvalue="logical"),
          function(object, na.rm, returnvalue){
            if(na.rm){
              removeNA <- function(x) return(x[!is.na(x)]) # pulls NA out of vectors
            }
            done <- FALSE # has it been collapsed yet?
            # test to see if there is one ploidy for the whole set
            u <- unique(as.vector(pld(object)))
            if(na.rm) u <- removeNA(u)
            if(length(u)==1){
              done <- TRUE
              if(returnvalue){
                result <- new("ploidyone")
                pld(result) <- u
                return(result)
              } else {
                return(TRUE)
              }
            }
            if(!done){
              # test to see if there is one ploidy per sample
              u <- sapply(data.frame(t(pld(object))), unique, simplify=FALSE)
              if(na.rm) u <- sapply(u, removeNA, simplify=FALSE)
              if(all(sapply(u, length, simplify=TRUE)==1)){
                done <- TRUE
                if(returnvalue){
                  result <- new("ploidysample", samples=dimnames(pld(object))[[1]])
                  pld(result) <- unlist(u)
                  return(result)
                } else {
                  return(TRUE)
                }
              }
            }
            if(!done){
              # test to see if there is one ploidy per locus
              u <- sapply(data.frame(pld(object)), unique, simplify=FALSE)
              if(na.rm) u <- sapply(u, removeNA, simplify=FALSE)
              if(all(sapply(u, length, simplify=TRUE)==1)){
                done <- TRUE
                if(returnvalue){
                  result <- new("ploidylocus", loci=dimnames(pld(object))[[2]])
                  pld(result) <- unlist(u)
                  return(result)
                } else {
                  return(TRUE)
                }
              }
            }
            if(!done) return(FALSE) # cannot be collapsed
          })
setMethod("plCollapse", signature(object="ploidysample", na.rm="logical",
                                  returnvalue="logical"),
          function(object, na.rm, returnvalue){
            u <- unique(pld(object))
            if(na.rm) u <- u[!is.na(u)]
            if(length(u) != 1){
              return(FALSE)
            } else {
              if(returnvalue){
                result <- new("ploidyone")
                pld(result) <- u
                return(result)
              } else {
                return(TRUE)
              }
            }
          })
setMethod("plCollapse", signature(object="ploidylocus", na.rm="logical",
                                  returnvalue="logical"),
          function(object, na.rm, returnvalue){
            u <- unique(pld(object))
            if(na.rm) u <- u[!is.na(u)]
            if(length(u) != 1){
              return(FALSE)
            } else {
              if(returnvalue){
                result <- new("ploidyone")
                pld(result) <- u
                return(result)
              } else {
                return(TRUE)
              }
            }
          })

##                                 gendata Class                             ##
## Contains info about missing data symbol, microsat repeats, sample ploidy, ##
## and populations.                                                          ##
## Is a superclass for other classes containing genotype data.               ##
setClass("gendata", representation(Description="character", Missing="ANY",
                                   Usatnts="integer", Ploidies="ploidysuper",
                                   PopInfo="integer", PopNames="character"),
         # make a virtual class?  (include "VIRTUAL" in contains)
         prototype(Description="Insert dataset description here.",
                   Missing=as.integer(-9)),
         validity=function(object){
             failures <- character(0)
             # check size and names of Ploidies with Usatnts and PopInfo
             if(is(object@Ploidies, "ploidymatrix")){
               np <- dimnames(pld(object@Ploidies))
               if(!identical(np[[1]],names(object@PopInfo)))
                 failures <- c(failures,
                               "Sample names inconsistent between Plodies and PopInfo")
               if(!identical(np[[2]],names(object@Usatnts)))
                 failures <- c(failures,
                               "Locus names inconsistent between Ploidies and Usatnts")
             }
             if(is(object@Ploidies,"ploidysample")){
               np <- names(pld(object@Ploidies))
               if(!identical(np, names(object@PopInfo)))
                 failures <- c(failures,
                        "Names of Ploidies do not match names of PopInfo.")
             }
             if(is(object@Ploidies,"ploidylocus")){
               np <- names(pld(object@Ploidies))
               if(!identical(np, names(object@Usatnts)))
                 failures <- c(failures,
                               "Names of Ploidies do not match names of Usatnts")
             }

             # check to see that Missing only has one element
             if(length(object@Missing) != 1) failures <- c(failures,
                 "Missing needs to have only one element.")

             # check to see that each number in PopInfo has an element in PopNames
             okpopnums <- 1:length(object@PopNames) #numbers that go with a PopName
             usedpopnums <- unique(object@PopInfo) #numbers in PopInfo
             usedpopnums <- usedpopnums[!is.na(usedpopnums)] # ignore NAs
             if(length(usedpopnums) > length(usedpopnums[usedpopnums %in% okpopnums]))
                failures <- c(failures,
                          "PopInfo contains integers that don't index PopNames."
                              )

             # check to see that all sample and locus names are unique
             if(length(unique(names(object@PopInfo))) < length(object@PopInfo)){
                 failures <- c(failures, "Not all sample names are unique.")
             }
             if(length(unique(names(object@Usatnts))) < length(object@Usatnts)){
                 failures <- c(failures, "Not all locus names are unique.")
             }

             # check to see that locus names do not contain periods
             if(length(grep(".", names(object@Usatnts), fixed=TRUE)) > 0){
                 failures <- c(failures, "Locus names may not contain periods.")
             }

             # return TRUE or failures
             if(length(failures) == 0){
                 return(TRUE)
             } else {
                 return(failures)
             }
         })

##                                  genambig class                                ##
##                               Subclass of gendata                              ##
## Contains two-dimensional list of vectors to contain genotypes where            ##
## allele copy number is not necessarily known.                                   ##
setClass("genambig", representation(Genotypes="array"), contains="gendata",
          validity = function(object){
    failures <- character(0)

    # do validity testing from gendata
    gendatatest <- getValidity(getClassDef("gendata"))(object)
    if(!identical(TRUE, gendatatest)){
        failures <- c(failures, gendatatest)
    }

    # check to see that the Genotypes list is two-dimensional
    if(length(dim(object@Genotypes)) != 2)
       failures <- c(failures, "Genotypes should be a two-dimensional list.")

    # include a check to see that all list elements are vectors?

    # check to see that sample names match in Genotypes, Ploidies, PopInfo
    ng <- dimnames(object@Genotypes)
    if(is(object@Ploidies,"ploidymatrix")){
      np <- dimnames(pld(object@Ploidies))
      if(!identical(np[[1]],ng[[1]]))
        failures <- c(failures,
                      "Sample names inconsistent between Plodies and Genotypes")
      if(!identical(np[[2]],ng[[2]]))
        failures <- c(failures,
                      "Locus names inconsistent between Ploidies and Genotypes")
    }
    if(is(object@Ploidies,"ploidysample")){
      np <- names(pld(object@Ploidies))
      if(!identical(np, ng[[1]]))
        failures <- c(failures,
               "Names of Ploidies do not match Sample names in Genotypes")
    }
    if(is(object@Ploidies,"ploidylocus")){
      np <- names(pld(object@Ploidies))
      if(!identical(np, ng[[2]]))
        failures <- c(failures,
               "Names of Ploidies do not match Locus names in Genotypes")
    }
    if(!identical(ng[[1]], names(object@PopInfo)))
        failures <- c(failures,
                      "Sample names are not the same in Genotypes and PopInfo."
                      )

    # check to see that locus names match in Genotypes, Usatnts
    if(!identical(ng[[2]], names(object@Usatnts)))
        failures <- c(failures,
                      "Locus names are not the same in Genotypes and Usatnts."
                      )

    # return TRUE or failures
    if(length(failures) == 0){
        return(TRUE)
    } else {
        return(failures)
    }
})

## genbinary class
## subclass of Gendata
## Genotypes stored as a matrix of symbols indicating the presence and
## absence of alleles.
setClass("genbinary", representation(Genotypes = "matrix", Present = "ANY",
                                     Absent = "ANY"), contains="gendata",
#         where = where,
         prototype(Present = as.integer(1), Absent = as.integer(0)),
         validity = function(object){
             failures <- character(0)

             # do validity testing from gendata
             gendatatest <- getValidity(getClassDef("gendata"))(object)
             if(!identical(TRUE, gendatatest)){
                 failures <- c(failures, gendatatest)
             }

             # check that there is only one element for Present and Absent
             if(length(object@Present)!=1) failures <- c(failures,
                      "object@Present must be a vector of length 1.")
             if(length(object@Absent)!=1) failures <- c(failures,
                      "object@Absent must be a vector of length 1.")

             # check that Present, Absent, and Missing are different
             pam <- c(object@Present, object@Absent, object@Missing)
             if(length(unique(pam))!=3) failures <- c(failures,
                      "Different values required for Present, Absent, and Missing.")

             # check that all symbols in Genotypes are in PAM
             matsym <- unique(object@Genotypes)
             if(length(matsym) > 0){
                 if(FALSE %in% (matsym %in% pam)) failures <- c(failures,
                    "All genotype symbols must match Present, Absent, or Missing.")
             }

             # return TRUE or failures
             if(length(failures) == 0){
                 return(TRUE)
             } else {
                 return(failures)
             }
         }
     )

### New generic functions for classes in polysat ###
# generic function to get sample names
setGeneric("Samples", function(object, populations, ploidies)
           standardGeneric("Samples"))
# generic to replace sample names
setGeneric("Samples<-", function(object, value) standardGeneric("Samples<-"))
# generic to get locus names
setGeneric("Loci", function(object, usatnts, ploidies) standardGeneric("Loci"))
# generic to replace locus names
setGeneric("Loci<-", function(object, value) standardGeneric("Loci<-"))
# generics to get and replace population identities
setGeneric("PopInfo", function(object) standardGeneric("PopInfo"))
setGeneric("PopInfo<-", function(object, value) standardGeneric("PopInfo<-"))
# generics to get and replace population names
setGeneric("PopNames", function(object) standardGeneric("PopNames"))
setGeneric("PopNames<-", function(object, value) standardGeneric("PopNames<-"))
# generics to get and replace ploidies
setGeneric("Ploidies", function(object, samples, loci) standardGeneric("Ploidies"))
setGeneric("Ploidies<-", function(object, value) standardGeneric("Ploidies<-"))
# generics to get and replace usatnts
setGeneric("Usatnts", function(object) standardGeneric("Usatnts"))
setGeneric("Usatnts<-", function(object, value) standardGeneric("Usatnts<-"))
# generics to delete samples and loci
setGeneric("deleteSamples", function(object, samples)
           standardGeneric("deleteSamples"))
setGeneric("deleteLoci", function(object, loci) standardGeneric("deleteLoci"))
# generics to get and replace missing data code
setGeneric("Missing", function(object) standardGeneric("Missing"))
setGeneric("Missing<-", function(object, value) standardGeneric("Missing<-"))
# generics to get and replace description
setGeneric("Description", function(object) standardGeneric("Description"))
setGeneric("Description<-", function(object,value)
           standardGeneric("Description<-"))
# generics to get and replace population number by name
setGeneric("PopNum", function(object, popname) standardGeneric("PopNum"))
setGeneric("PopNum<-", function(object, popname, value)
           standardGeneric("PopNum<-"))
# generics to get and edit genotype data
setGeneric("Genotype", function(object, sample, locus) standardGeneric("Genotype"))
setGeneric("Genotype<-", function(object, sample, locus, value)
           standardGeneric("Genotype<-"))
setGeneric("Genotypes", function(object, samples = Samples(object), loci = Loci(object))
           standardGeneric("Genotypes"))
setGeneric("Genotypes<-", function(object, samples = Samples(object),
                                   loci = Loci(object), value)
           standardGeneric("Genotypes<-"))
# generic to determine if a genotype is missing
setGeneric("isMissing", function(object, samples = Samples(object),
                                 loci = Loci(object))
           standardGeneric("isMissing"))
# generic to print genotypes to console
setGeneric("viewGenotypes",
           function(object, samples = Samples(object), loci = Loci(object))
           standardGeneric("viewGenotypes"))
# generic to edit genotypes
setGeneric("editGenotypes",
           function(object, maxalleles = max(Ploidies(object)),
                    samples = Samples(object),
                    loci = Loci(object)) standardGeneric("editGenotypes"))
# generic to estimate ploidy
setGeneric("estimatePloidy",
           function(object, extrainfo, samples = Samples(object), loci = Loci(object))
           standardGeneric("estimatePloidy"))
# generics to get and replace symbols for allele present/absent
setGeneric("Present", function(object) standardGeneric("Present"))
setGeneric("Present<-", function(object, value) standardGeneric("Present<-"))
setGeneric("Absent", function(object) standardGeneric("Absent"))
setGeneric("Absent<-", function(object, value) standardGeneric("Absent<-"))
# Generic function and methods for finding all unique genotypes and indexing which samples have them
setGeneric("genIndex", function(object, locus) standardGeneric("genIndex"))

# function for making locus names compatible as column headers
fixloci <- function(loci, warn = TRUE){
    loci2 <- make.names(loci)
    loci2 <- gsub(".", "", loci2, fixed = TRUE)
    if(warn && !identical(loci, loci2)){
      warning("Special characters removed from locus names.")
    }
    return(loci2)
}

#### gendata methods
# initialization for a gendata object
setMethod("initialize",
    signature(.Object = "gendata"),
    function (.Object, samples, loci, Missing, ...)
    {
        if(missing(samples)) samples <- c("ind1","ind2")
        if(missing(loci)) loci <- c("loc1","loc2")

        # fix locus names so that they can be column headers if necessary
        loci <- fixloci(loci)

        # make a vector to contain repeat lengths
        usatnts <- as.integer(rep(NA, length(loci)))
        names(usatnts) <- loci
        .Object@Usatnts <- usatnts

        # make a matrix to contain ploidy
        .Object@Ploidies <- new("ploidymatrix",samples=samples,loci=loci)

        # make a vector to contain population identity
        popinfo <- as.integer(rep(NA, length(samples)))
        names(popinfo) <- samples
        .Object@PopInfo <- popinfo

        # put in the missing data symbol
        if(missing(Missing)) Missing <- as.integer(-9)
        .Object@Missing <- Missing

        # call default method
        callNextMethod(.Object, ...)
    }
)

## methods to get sample names
# All samples if only gendata object is given
setMethod("Samples", signature(object = "gendata", populations = "missing",
                               ploidies = "missing"),
          function(object){return(names(object@PopInfo))})
# Just for a subset of population names
setMethod("Samples", signature(object = "gendata", populations = "character",
                               ploidies = "missing"),
          function(object, populations){
              pops <- match(populations, object@PopNames)
              return(names(object@PopInfo)[object@PopInfo %in% pops])
              })
# Population names and ploidies
setMethod("Samples", signature(object = "gendata", populations = "character",
                               ploidies = "numeric"),
          function(object, populations, ploidies){
            if(!is(object@Ploidies,"ploidysample"))
        stop("Ploidies argument only valid if ploidies in dataset are indexed by sample.")
              ploidies <- as.integer(ploidies)
              pops <- match(populations, object@PopNames)
              popsam <- names(object@PopInfo)[object@PopInfo %in% pops]
              ploisam <- names(pld(object@Ploidies))[pld(object@Ploidies) %in% ploidies]
              return(popsam[popsam %in% ploisam])
              })
# Just a subset of population numbers
setMethod("Samples", signature(object = "gendata", populations = "numeric",
                               ploidies = "missing"),
          function(object, populations){
              populations <- as.integer(populations)
              return(names(object@PopInfo)[object@PopInfo %in% populations])
              })
# Population numbers and ploidies
setMethod("Samples", signature(object = "gendata", populations = "numeric",
                               ploidies = "numeric"),
          function(object, populations, ploidies){
            if(!is(object@Ploidies,"ploidysample"))
        stop("Ploidies argument only valid if ploidies in dataset are indexed by sample.")

              ploidies <- as.integer(ploidies)
              populations <- as.integer(populations)
              popsam <- names(object@PopInfo)[object@PopInfo %in% populations]
              ploisam <- names(pld(object@Ploidies))[pld(object@Ploidies) %in% ploidies]
              return(popsam[popsam %in% ploisam])
              })
# Just ploidies
setMethod("Samples", signature(object = "gendata", populations = "missing",
                               ploidies = "numeric"),
          function(object, ploidies){
            if(!is(object@Ploidies,"ploidysample"))
        stop("Ploidies argument only valid if ploidies in dataset are indexed by sample.")

            ploidies <- as.integer(ploidies)
            return(names(pld(object@Ploidies))[pld(object@Ploidies) %in% ploidies])
              })

##                     Replacement method for sample names                   ##
setReplaceMethod("Samples", "gendata", function(object, value){
    # check that all sample names are unique
    if(length(unique(value)) < length(value)){
        stop("Not all sample names are unique")
    }
  if(is(object@Ploidies,"ploidysample"))
    names(object@Ploidies@pld) <- value
  if(is(object@Ploidies,"ploidymatrix"))
    dimnames(object@Ploidies@pld)[[1]] <- value
  names(object@PopInfo) <- value
  object
})

## Methods to get locus names.
# Return all locus names
setMethod("Loci", signature(object = "gendata", usatnts = "missing", ploidies="missing"),
          function(object){
              return(names(object@Usatnts))
              })
# Return locus names for a certain subset of repeat types
setMethod("Loci", signature(object = "gendata", usatnts = "numeric", ploidies="missing"),
          function(object, usatnts){
              usatnts <- as.integer(usatnts)
              return(names(object@Usatnts)[object@Usatnts %in% usatnts])
          })
# Return locus names for a certain subset of ploidies
setMethod("Loci", signature(object = "gendata", usatnts = "missing", ploidies="numeric"),
          function(object, ploidies){
            if(!is(object@Ploidies,"ploidylocus"))
          stop("Ploidies argument only valid if ploidies in dataset are indexed by locus")

            ploidies <- as.integer(ploidies)
            loci <- names(pld(object@Ploidies))[pld(object@Ploidies) %in% ploidies]

            return(loci)
          })
# Return locus names for a subset of ploidies and repeat types
setMethod("Loci", signature(object = "gendata", usatnts = "numeric", ploidies="numeric"),
          function(object, usatnts, ploidies){
            if(!is(object@Ploidies,"ploidylocus"))
          stop("Ploidies argument only valid if ploidies in dataset are indexed by locus")

            ploidies <- as.integer(ploidies)
            ploci <- names(pld(object@Ploidies))[pld(object@Ploidies) %in% ploidies]

            usatnts <- as.integer(usatnts)
            uloci <- names(object@Usatnts)[object@Usatnts %in% usatnts]

            return(uloci[uloci %in% ploci])
          })

## Replacement method for locus names
setReplaceMethod("Loci", "gendata", function(object, value){
    # fix names so that they can be column headers if necessary
    value <- fixloci(value)
    # check that all locus names are unique
    if(length(unique(value)) < length(value)){
        stop("Not all locus names are unique")
    }
    # check that locus names do not contain periods
#    if(length(grep(".", value, fixed=TRUE)) > 0){
#        stop("Locus names may not contain periods.")
#    }
  # Replace Ploidies names if necessary
  if(is(object@Ploidies,"ploidylocus"))
    names(object@Ploidies@pld) <- value
  if(is(object@Ploidies,"ploidymatrix"))
    dimnames(object@Ploidies@pld)[[2]] <- value

  # Replace Usatnts names
  names(object@Usatnts) <- value
  object
})

## Methods to get and replace population identities
setMethod("PopInfo", signature(object = "gendata"),
          function(object) return(object@PopInfo))
setReplaceMethod("PopInfo", "gendata", function(object, value){
    # add values to PopInfo slot
    object@PopInfo[names(object@PopInfo)] <- as.integer(value)
    # make sure there are enough PopNames
    if(!all(is.na(object@PopInfo)) &&
        length(object@PopNames) < max(object@PopInfo, na.rm=TRUE)){
        missingnames <- (length(object@PopNames) + 1):max(object@PopInfo,
                                                          na.rm=TRUE)
        object@PopNames[missingnames] <- paste("Pop", missingnames,sep="")
    }
    # return object
    object
})

## Methods to get and replace population names
setMethod("PopNames", signature(object = "gendata"),
          function(object) return(object@PopNames))
setReplaceMethod("PopNames", "gendata", function(object, value){
    object@PopNames <- value
    object
})

## Methods to get and replace ploidies
setMethod("Ploidies", signature(object = "gendata", samples="ANY", loci="ANY"),
          function(object, samples, loci){
            return(pld(object@Ploidies, samples=samples, loci=loci))
            })
setReplaceMethod("Ploidies", "gendata", function(object, value){
    if(!all(is.na(as.integer(value)) | as.integer(value) >= 0))
        stop("Ploidy less than zero not allowed.")
    if(!all(is.na(as.integer(value)) | as.integer(value) <= 12))
        warning("Large ploidies (>12) detected.")
    pld(object@Ploidies) <- value
    object
})
# Function to change the ploidy format between the four types
reformatPloidies <- function(object, output="collapse", na.rm=FALSE,
                             erase=FALSE){

  if(!output %in% c("matrix", "sample", "locus", "one", "collapse"))
    stop("Unrecognized output argument.  See reformatPloidies documentation.")

  if(erase){ # just make a new ploidies object filled with NA
    if(output=="collapse") output <- "one"
    object@Ploidies <- new(paste("ploidy",output,sep=""),
                           samples=Samples(object),
                           loci=Loci(object))
  } else {

    ## If you start with ploidymatrix - collapse or do nothing
  if(is(object@Ploidies,"ploidymatrix") && output %in% c("sample","locus","one",
                                                         "collapse")){
    pl <- plCollapse(object@Ploidies, na.rm=na.rm, returnvalue=TRUE)
    if(is(pl, "ploidysuper")){ # if the ploidy was collapsible
      if(output=="collapse" || is(pl, paste("ploidy",output,sep=""))){
        # if collapsed ploidy is the right format, we're done
        object@Ploidies <- pl
      } else { # if collapsed ploidy is not right format
        if(is(pl, "ploidyone")){ # expand from ploidyone
          pl2 <- new(paste("ploidy",output,sep=""),
                     samples=Samples(object),
                     loci=Loci(object))
          pld(pl2) <- pld(pl)
          object@Ploidies <- pl2
        } else {
          stop("Ploidies not collapsible to desired format.")
        }
      }
    } else { # if ploidy was not collapsible at all
      if(output != "collapse") stop("Ploidies not collapsible to desired format.")
    }
  }

  # if you start with ploidysample and want to collapse
  if(is(object@Ploidies,"ploidysample") && output %in% c("one","collapse","locus")){
    pl <- plCollapse(object@Ploidies, na.rm=na.rm, returnvalue=TRUE)
    if(is(pl, "ploidyone")){ # if it is collapsible
      if(output %in% c("collapse","one")){
        object@Ploidies <- pl
      }
      if(output=="locus"){
        pl2 <- new("ploidylocus", loci=Loci(object))
        pld(pl2) <- pld(pl)
        object@Ploidies <- pl2
      }
    } else { # if it is not collapsible
      if(output != "collapse") stop("Ploidies not collapsible to desired format.")
    }
  }

  # if you start with ploidysample and want to expand
  if(is(object@Ploidies,"ploidysample") && output=="matrix"){
    pl <- new("ploidymatrix", samples=Samples(object), loci=Loci(object))
    pld(pl) <- pld(object@Ploidies)
    object@Ploidies <- pl
  }

  # if you start with ploidylocus and want to collapse
  if(is(object@Ploidies,"ploidylocus") && output %in% c("sample","one","collapse")){
    pl <- plCollapse(object@Ploidies, na.rm=na.rm, returnvalue=TRUE)
    if(is(pl, "ploidyone")){ # if it is collapsible
      if(output %in% c("collapse","one")){
        object@Ploidies <- pl
      }
      if(output=="sample"){
        pl2 <- new("ploidysample", samples=Samples(object))
        pld(pl2) <- pld(pl)
        object@Ploidies <- pl2
      }
    } else { # if it is not collapsible
      if(output != "collapse") stop("Ploidies not collapsible to desired format.")
    }
  }

  # if you start with ploidylocus and want to expand
  if(is(object@Ploidies,"ploidylocus") && output == "matrix"){
    pl <- new("ploidymatrix", samples=Samples(object), loci=Loci(object))
    pld(pl) <- rep(pld(object@Ploidies), each=length(Samples(object)))
    object@Ploidies <- pl
  }

  # if you start with ploidyone and want to expand
  if(is(object@Ploidies,"ploidyone") && output %in% c("matrix","sample","locus")){
    pl <- new(paste("ploidy",output,sep=""),
                     samples=Samples(object), loci=Loci(object))
    pld(pl) <- pld(object@Ploidies)
    object@Ploidies <- pl
  }
  }

  return(object)
}


## Methods to get and replace repeat lengths
setMethod("Usatnts", signature(object = "gendata"),
          function(object) return(object@Usatnts))
setReplaceMethod("Usatnts", "gendata", function(object, value){
    object@Usatnts[names(object@Usatnts)] <- as.integer(value)
    object
})

## summary method for gendata
setMethod("summary", "gendata", function(object){
    cat(paste(length(Samples(object)), "samples,", length(Loci(object)), "loci."),
        paste(length(unique(PopInfo(object))), "populations."),
        sep="\n")
    cat("Ploidies:", unique(as.vector(Ploidies(object))))
    cat("\n", sep="")
    cat("Length(s) of microsatellite repeats:", unique(Usatnts(object)))
    cat("\n", sep="")
})

## Methods to delete samples and loci
setMethod("deleteSamples", "gendata", function(object, samples){
    samtouse <- Samples(object)[!Samples(object) %in% samples]
    object@PopInfo <- object@PopInfo[samtouse]
    if(is(object@Ploidies,"ploidysample")){
      object@Ploidies@pld <- object@Ploidies@pld[samtouse]
    }
    if(is(object@Ploidies,"ploidymatrix")){
      object@Ploidies@pld <- object@Ploidies@pld[samtouse,, drop=FALSE]
    }
    return(object)
})

setMethod("deleteLoci", "gendata", function(object, loci){
    loctouse <- Loci(object)[!Loci(object) %in% loci]
    object@Usatnts <- object@Usatnts[loctouse]
    if(is(object@Ploidies,"ploidylocus")){
      object@Ploidies@pld <- object@Ploidies@pld[loctouse]
    }
    if(is(object@Ploidies,"ploidymatrix")){
      object@Ploidies@pld <- object@Ploidies@pld[,loctouse, drop=FALSE]
    }
    return(object)
})

## Subscripting method
setMethod("[", signature(x="gendata", i="ANY", j="ANY"), function(x, i, j){
    # i is samples, j is loci
    x@PopInfo <- x@PopInfo[i]
    if(is(x@Ploidies, "ploidymatrix")){
      x@Ploidies@pld <- x@Ploidies@pld[i,j, drop=FALSE]
    }
    if(is(x@Ploidies, "ploidysample")){
      x@Ploidies@pld <- x@Ploidies@pld[i]
    }
    if(is(x@Ploidies, "ploidylocus")){
      x@Ploidies@pld <- x@Ploidies@pld[j]
    }
    x@Usatnts <- x@Usatnts[j]
    return(x)
})

## Methods to get and replace missing data code
setMethod("Missing", "gendata", function(object) return(object@Missing))
setReplaceMethod("Missing", "gendata", function(object, value){
    if(length(value) != 1) stop("Assigned value should only have one element.")
    object@Missing <- value
    object
})

## Methods to get and replace description
setMethod("Description", "gendata", function(object) return(object@Description))
setReplaceMethod("Description", "gendata", function(object, value){
    object@Description <- value
    object
})

## Methods to get and replace population numbers by name
setMethod("PopNum", signature(object="gendata", popname="character"),
          function(object, popname){
              return(match(popname, PopNames(object)))
          })
setReplaceMethod("PopNum", signature(object="gendata", popname="character"),
                 function(object, popname, value){
                     value <- as.integer(value)
                     if(length(Samples(object, populations=value))>0 &&
                        length(Samples(object, populations=popname))>0 &&
                        !identical(Samples(object, populations=value),
                                   Samples(object, populations=popname)))
                         cat(paste("Samples already present in population ",
                                   value,".  These have been merged into",
                                   popname," .", sep=""), sep="\n")
                     PopInfo(object)[Samples(object, populations=popname)] <-
                         value
                     PopNames(object)[PopNum(object, popname)] <- NA
                     PopNames(object)[value] <- popname
                     object
                     })

## Method for merging objects
setMethod("merge", signature(x="gendata", y="gendata"),
          function(x, y, objectm, samples, loci, overwrite){
            # error if ploidies are not in the same format
            if(!identical(class(x@Ploidies),class(y@Ploidies)))
              stop("x and y must have Ploidies in the same format.")

              # set up new gendata object if this wasn't called from the method
              # of a subclass.
              if(missing(objectm)){
                  if(missing(samples)) samples <- unique(c(Samples(x), Samples(y)))
                  if(missing(loci)) loci <- unique(c(Loci(x), Loci(y)))
                  objectm <- new("gendata", samples=samples, loci=loci)
              }

              # determine what to do with conflicting data: overwrite or give error
              if(missing(overwrite)) overwrite <- "empty"
              owerror <- ifelse(overwrite != "x" && overwrite != "y", TRUE, FALSE)
              if(overwrite=="x"){
                  objectB <- x
                  objectT <- y
              } else {
                  objectB <- y
                  objectT <- x
              }

            # function definition for merging vectors
            # B and T are vectors, m is final index (Samples(objectm) or Loci(objectm)), x is output
            vectmerge <- function(m, B, T, type){
              x <- rep(NA, length(m)) # set up an empty vector
              names(x) <- m

              b <- m[m %in% names(B)] # subset indexing vectors
              t <- m[m %in% names(T)]

              x[b] <- B[b]
              x[t] <- T[t]

              if(owerror && !identical(x[b], B[b]))
                stop(paste("Conflicting", type,"data; set overwrite to \"x\" or \"y\"."))

              return(x)
            }

              # merge Usatnts
            Usatnts(objectm) <- vectmerge(Loci(objectm), Usatnts(objectB), Usatnts(objectT), "Usatnts")

              # merge Ploidies
            if(is(objectT@Ploidies,"ploidysample")){
              objectm@Ploidies <- new("ploidysample", samples=Samples(objectm))
              Ploidies(objectm) <- vectmerge(Samples(objectm), Ploidies(objectB), Ploidies(objectT),
                                             "Ploidies")
            }
            if(is(objectT@Ploidies,"ploidylocus")){
              objectm@Ploidies <- new("ploidylocus", loci=Loci(objectm))
              Ploidies(objectm) <- vectmerge(Loci(objectm), Ploidies(objectB), Ploidies(objectT),
                                             "Ploidies")
            }
            if(is(objectT@Ploidies,"ploidyone")){
              if(owerror && Ploidies(objectB) != Ploidies(objectT))
                stop("Conflicting Ploidies data. Set overwrite to \"x\" or \"y\".")
              objectm@Ploidies <- new("ploidyone")
              Ploidies(objectm) <- Ploidies(objectT)
            }
            if(is(objectT@Ploidies,"ploidymatrix")){
              a <- Samples(objectm)[Samples(objectm) %in% Samples(objectB)]
              b <- Loci(objectm)[Loci(objectm) %in% Loci(objectB)]
              c <- Samples(objectm)[Samples(objectm) %in% Samples(objectT)]
              d <- Loci(objectm)[Loci(objectm) %in% Loci(objectT)]

              Ploidies(objectm)[a,b] <- Ploidies(objectB)[a,b]
              Ploidies(objectm)[c,d] <- Ploidies(objectT)[c,d]

              if(owerror && !identical(Ploidies(objectm)[a,b], Ploidies(objectB)[a,b]))
                stop("Conflicting Ploidies data. Set overwrite to \"x\" or \"y\".")
            }

              # take description from "top" object or both
              if(owerror && Description(objectB)[1] != Description(objectT)[1]){
                  Description(objectm) <- c(Description(objectT), Description(objectB))
              } else {
                  Description(objectm) <- Description(objectT)
              }

              ## merge PopInfo and PopNames
              # Fill vector of object names
              PopNames(objectm) <- unique(c(PopNames(objectT),PopNames(objectB)))
              # Scenerio for identical sample sets
              if(identical(Samples(objectB), Samples(objectT))){
                  if(owerror && !identical(PopInfo(objectB),PopInfo(objectT))){
                stop("Conflicting data in PopInfo. Set overwrite to \"x\" or \"y\".")
                  } else {
                      PopInfo(objectm) <- PopInfo(objectT)[Samples(objectm)]
                  }
              } else {
              # Scenario for non-overlapping or partially overlapping sample sets
                  for(p in PopNames(objectm)){
                      # fill in values from "bottom" object from this population
                      if(p %in% PopNames(objectB)){
                          samtofill <- Samples(objectB,populations=p)
                          samtofill <- samtofill[samtofill %in% Samples(objectm)]
                          # give error if something would be overwritten
                          if(owerror && (FALSE %in% is.na(PopInfo(objectm)[samtofill])))
                      stop("Conflicting data in PopInfo.  Set overwrite to \"x\" or \"y\".")
                          # only fill in values that are currently NA
                          PopInfo(objectm)[samtofill][is.na(PopInfo(objectm)[samtofill])] <-
                              PopNum(objectm, p)
                      }
                      if(p %in% PopNames(objectT)){
                          samtofill <- Samples(objectT, populations=p)
                          samtofill <- samtofill[samtofill %in% Samples(objectm)]
                          # give error if something would be overwritten
                          alreadyfilled <- samtofill[!is.na(PopInfo(objectm)[samtofill])]
                          if(owerror && !identical(PopInfo(objectm)[alreadyfilled],
                                                   PopInfo(objectT)[alreadyfilled])) stop(
                            "Conflicting data in PopInfo.  Set overwrite to \"x\" or \"y\".")
                          # fill in the population number
                          PopInfo(objectm)[samtofill] <-
                              PopNum(objectm, p)
                      }
                  }
              }
              #
              # easiest scenario- PopNames and PopInfo numbers are the same, do like above
              # more complicated- different sets of populations, but numbers are consistent
              # most complicated- numbers have different meanings in the two objects

              # return the merged object.
              return(objectm)
          })

#### genambig methods
# initialization for a genambig object
setMethod("initialize", signature(.Object = "genambig"),
          function(.Object, samples, loci, Missing, ...){
              # check to see that arguments are all there
              if(missing(samples)) samples <- c("ind1","ind2")
              if(missing(loci)) loci <- c("loc1","loc2")
              if(missing(Missing)) Missing <- as.integer(-9)
              
              loci <- fixloci(loci) # fix locus names

              # add empty genotype list
    .Object@Genotypes <- array(list(Missing), dim=c(length(samples),
                                                length(loci)),
                                 dimnames=list(samples, loci))

              callNextMethod(.Object, samples = samples, loci = loci, Missing=Missing,...)
              })

# replacement method for sample names
setReplaceMethod("Samples", "genambig", function(object, value){
    dimnames(object@Genotypes)[[1]] <- value
    callNextMethod(object, value)
})

# replacement method for locus names
setReplaceMethod("Loci", "genambig", function(object, value){
    value <- fixloci(value) # correct locus names if necessary
    dimnames(object@Genotypes)[[2]] <- value
    callNextMethod(object, value)
})

## Methods to get and edit genotype data
setMethod("Genotype", signature(object = "genambig", sample = "ANY",
                                locus = "ANY"),
          function(object, sample, locus){
              if(length(sample) != 1 | length(locus) != 1){
                  stop("sample and locus arguments must have only one element.")
              }
              return(object@Genotypes[[sample, locus]])
          })
setReplaceMethod("Genotype", "genambig", function(object, sample, locus, value){
    if(length(sample) != 1 | length(locus) != 1){
        stop("sample and locus arguments must have only one element.")
    }
    if(!is.vector(value)) stop("Assigned value must be a vector.")
    if(length(value) == 0) stop("Genotype vector must have at least one element.  See ?Missing.")
    if(any(is.na(value))) stop("NA values not allowed in genotypes.  See ?Missing.")
    object@Genotypes[[sample, locus]] <- value
    object
})
setMethod("Genotypes", signature(object = "genambig", samples = "ANY",
                                 loci = "ANY"),
          function(object, samples, loci){
            return(object@Genotypes[samples, loci, drop=FALSE])
          })
setReplaceMethod("Genotypes", "genambig",
                 function(object, samples,
                          loci, value){
                     if(length(dim(value))> 2)
                         stop("Assigned value has too many dimensions.")
                     if(FALSE %in% mapply(is.vector, value))
                         stop("All array elements must be vectors.")
                     if(0 %in% mapply(length, value))
                       stop("All genotype vectors must have at least one element.  See ?Missing.")
                     if(any(mapply(function(x) any(is.na(x)), value)))
                       stop("NA values not allowed in genotypes.  see ?Missing.")
    object@Genotypes[samples, loci] <- value
    object
})

## Method to determine if a genotype is missing
setMethod("isMissing", "genambig", function(object, samples, loci){
    if(length(samples) == 1 && length(loci) == 1){
        return(Genotype(object, samples, loci)[1] == Missing(object))
    } else {
        result <- array(FALSE, dim=c(length(samples), length(loci)),
                        dimnames=list(samples, loci))
        for(s in samples){
            for(L in loci){
                result[s,L] <- isMissing(object, s, L)
            }
        }
        return(result)
    }
})

setMethod(
  f = "show",
  signature = "genambig",
  function(object){
    cat(Description(object), "\n")
    cat("Missing values: ", (Missing(object)), "\n")
    cat("\nGenotypes: \n")
    genSum <- apply(object@Genotypes, 2,
                    function(x)
                      unlist(lapply(x, function(y)
                        paste(sort(y), collapse = "/"))))
    print(genSum, quote = FALSE)

    cat("\nSSR motif lengths: \n")
    if(all(is.na(Usatnts(object)))){
      cat(" none\n")
    } else {
      print(Usatnts(object))
    }

    cat("\nPloidies: \n")
    if(all(is.na(Ploidies(object)))){
      cat(" none\n")
    } else {
      print(Ploidies(object))
    }

    cat("\nPopNames: \n")
    if(length(PopNames(object)) > 0){
      print(PopNames(object))
    } else {
      cat(" none\n")
    }

    cat("\nPopInfo: \n")
    if(all(is.na(PopInfo(object)))){
      cat(" none\n")
    } else {
      print(PopInfo(object))
    }
  }
)

# summary method for genambig
setMethod("summary", "genambig", function(object){
    nummissing <- 0
    for(L in Loci(object)){
        for(s in Samples(object)){
            if(isMissing(object, s, L)) nummissing <- nummissing + 1
        }
    }
    cat("Dataset with allele copy number ambiguity.",
        Description(object),
        paste("Number of missing genotypes:", nummissing), sep="\n")
    callNextMethod(object)
})

# Method to print the genotypes to the console
setMethod("viewGenotypes", signature(object = "genambig", samples = "ANY",
                                     loci = "ANY"),
          function(object, samples, loci){
                  # print a header
                  cat("Sample\tLocus\tAlleles", sep="\n")

                  # loop through genotypes and print them
                  for(L in loci){
                      for(s in samples){
                          cat(s, L, Genotype(object, s, L), sep="\t")
                          cat("\n", sep="")
                      }
                  }
          })

# method for deleting samples
setMethod("deleteSamples", "genambig", function(object, samples){
    samtouse <- Samples(object)[!Samples(object) %in% samples]
    object@Genotypes <- as.array(object@Genotypes[samtouse, Loci(object)])
    callNextMethod(object, samples)
})

# method for deleting loci
setMethod("deleteLoci", "genambig", function(object, loci){
    loctouse <- Loci(object)[!Loci(object) %in% loci]
    object@Genotypes <- as.array(object@Genotypes[Samples(object), loctouse])
    callNextMethod(object, loci)
})

# subscripting method
setMethod("[", signature(x = "genambig", i="ANY", j="ANY"), function(x, i, j){
    x@Genotypes <- x@Genotypes[i,j, drop=FALSE]
    x@PopInfo <- x@PopInfo[i]
    if(is(x@Ploidies, "ploidymatrix")){
      x@Ploidies@pld <- x@Ploidies@pld[i,j, drop=FALSE]
    }
    if(is(x@Ploidies, "ploidysample")){
      x@Ploidies@pld <- x@Ploidies@pld[i]
    }
    if(is(x@Ploidies, "ploidylocus")){
      x@Ploidies@pld <- x@Ploidies@pld[j]
    }
    x@Usatnts <- x@Usatnts[j]
    return(x)
})

# method to change code for missing genotypes
setReplaceMethod("Missing", "genambig", function(object, value){
    if(length(value) != 1) stop("Assigned value should have only one element.")
    for(L in Loci(object)){
        for(s in Samples(object)){
            if(isMissing(object, s, L)) Genotype(object, s, L) <- value
        }
    }
    callNextMethod(object, value)
})

# Method to edit genotypes in a spreadsheet-like format
setMethod("editGenotypes", "genambig",
          function(object, maxalleles, samples, loci){
              if(is.na(maxalleles))
                  stop("Enter ploidies first or give argument for maxalleles.")
              # take the selected genotypes and convert to a data frame
              dummyarray <- matrix(NA, nrow=length(samples)*length(loci), ncol=maxalleles,
                                   dimnames=list(NULL,
                                   paste("Allele", 1:maxalleles, sep="")))
              dummyrow <- 1
              for(L in loci){
                  for(s in samples){
                      numalleles <- length(Genotype(object, s, L))
                      if(numalleles > maxalleles)
                          stop(paste("Set maxalleles to at least", numalleles))
                      dummyarray[dummyrow, 1:numalleles] <- Genotype(object, s, L)
                      dummyrow <- dummyrow + 1
                  }
              }
              samvect <- rep(samples, times=length(loci))
              locvect <- rep(loci, each=length(samples))
              genframe <- data.frame(Samples=samvect, Loci=locvect, dummyarray,
                                     stringsAsFactors=FALSE)
              # open data frame for editing by user
              cat("Edit the alleles, then close the data editor window.", sep="\n")
              genframe <- edit(genframe)
              # overwrite genotypes from those in data frame
              for(r in 1:dim(genframe)[1]){
                  thesealleles <- genframe[r, 3:dim(genframe)[2]]
                  thesealleles <- unique(thesealleles[!is.na(thesealleles)])
                  Genotype(object, genframe[r, "Samples"],
                           genframe[r, "Loci"]) <- thesealleles
              }

              # return object
              return(object)
          })

# Generic function and method to estimate ploidy
setMethod("estimatePloidy", "genambig",
          function(object, extrainfo, samples, loci){
            # get Ploidies into ploidysample format if necessary
              if(!is(object@Ploidies, "ploidysample")){
                object <- reformatPloidies(object, output="sample", erase=TRUE)
              }
              # set up array to contain the maximum and average number of alleles
              ploidyinfo <- array(dim=c(length(samples),2),
                                  dimnames=list(samples,
                                  c("max.alleles","mean.alleles")))
              # fill the array
              for(s in samples){
                  genotypes <- Genotypes(object, s, loci)[!isMissing(object, s, loci)]
                  if(length(genotypes) ==0){
                      ploidyinfo[s,] <- c(NA, NA)
                  } else {
                      numalleles <- mapply(length, genotypes)
                      ploidyinfo[s, "max.alleles"] <- max(numalleles)
                      ploidyinfo[s, "mean.alleles"] <- mean(numalleles)
                  }
              }

              ## Build data frame
              ploidyinfo <- as.data.frame(ploidyinfo)
              # Add in extrainfo, allowing for named or unnamed vectors, arrays, and
              # data frames.
              if(!missing(extrainfo)){
                  if(is.vector(extrainfo)){
                      if(!is.null(names(extrainfo))) extrainfo <- extrainfo[samples]
                      ploidyinfo <- data.frame(extrainfo = extrainfo, ploidyinfo)
                  } else {
                      if(is.array(extrainfo) | is.data.frame(extrainfo)){
                          extrainfo <- as.data.frame(extrainfo)
                          if(!identical(row.names(extrainfo),
                                        as.character(1:dim(extrainfo)[1])))
                              extrainfo <- extrainfo[samples,]
                          ploidyinfo <- data.frame(extrainfo, ploidyinfo)
                      }
                  }
              }
              # Add in PopInfo from the object
              if(length(unique(PopInfo(object))) >1)
                  ploidyinfo <- data.frame(population = PopInfo(object)[samples],
                                           ploidyinfo)
              # Add in Ploidies from the object
              if(length(unique(Ploidies(object))) >1)
                  ploidyinfo <- data.frame(ploidyinfo,
                                           current.ploidy = Ploidies(object)[samples])
              # Make a column for new ploidies
              ploidyinfo <- data.frame(ploidyinfo, new.ploidy = ploidyinfo$max.alleles)

              # let user edit data frame
              cat("Edit the new.ploidy values, then close the data editor window.",
                  sep="\n")
              ploidyinfo <- edit(ploidyinfo)

              # write new.ploidy to object@Ploidies
              Ploidies(object)[samples] <- ploidyinfo$new.ploidy

              # return the object with edited ploidies
              return(object)
          })

setMethod("merge", signature(x="genambig", y="genambig"),
          function(x, y, objectm, samples, loci, overwrite){
              # set up new genambig object if this wasn't called from the method
              # of a subclass.
              if(missing(objectm)){
                  if(missing(samples)) samples <- unique(c(Samples(x), Samples(y)))
                  if(missing(loci)) loci <- unique(c(Loci(x), Loci(y)))
                  objectm <- new("genambig", samples=samples, loci=loci)
              }

              # determine what to do with conflicting data: overwrite or give error
              if(missing(overwrite)) overwrite <- "empty"
              owerror <- ifelse(overwrite != "x" && overwrite != "y", TRUE, FALSE)
              if(overwrite=="x"){
                  objectB <- x
                  objectT <- y
              } else {
                  objectB <- y
                  objectT <- x
              }


              # make sure missing data values are the same, change if appropriate
              if(Missing(x) != Missing(y)){
                  if(owerror) stop(
                   "Different missing data symbols used; set overwrite to \"x\" or \"y\".")
                  Missing(objectB) <- Missing(objectT)
              }
              Missing(objectm) <- Missing(objectT)

              # merge genotypes into objectm
              # put in objectB first
              samtofillB <- Samples(objectm)[Samples(objectm) %in% Samples(objectB)]
              loctofillB <- Loci(objectm)[Loci(objectm) %in% Loci(objectB)]
              Genotypes(objectm, samtofillB, loctofillB) <-
                  Genotypes(objectB, samtofillB, loctofillB)
              # then put in objectT
              samtofillT <- Samples(objectm)[Samples(objectm) %in% Samples(objectT)]
              loctofillT <- Loci(objectm)[Loci(objectm) %in% Loci(objectT)]
              Genotypes(objectm, samtofillT, loctofillT) <-
                  Genotypes(objectT, samtofillT, loctofillT)
              # check to see if objectB overwritten if there wasn't supposed to be
              # conflicting data
              if(owerror && !identical(Genotypes(objectm, samtofillB, loctofillB),
                                       Genotypes(objectB, samtofillB, loctofillB)))
                 stop("Conflicting genotype data.  Set overwrite to \"x\" or \"y\".")

              # call gendata method to merge popinfo etc
              callNextMethod(x, y, objectm, overwrite=overwrite)
          })

#### genbinary methods
setMethod("initialize", signature(.Object = "genbinary"),
          function(.Object, samples, loci, ...){
              # fill in empty arguments if necessary
              if(missing(samples)) samples <- c("ind1", "ind2")
              if(missing(loci)) loci <- c("loc1", "loc2")
#              if(missing(Missing)) Missing <- as.integer(-9)
#              if(missing(Present)) Present <- as.integer(1)
#              if(missing(Absent)) Absent <- as.integer(0)

              loci <- fixloci(loci)

              # add empty genotype matrix
              .Object@Genotypes <- matrix(nrow=length(samples), ncol=0,
                                          dimnames=list(samples, NULL))
              # go to the gendata method
              callNextMethod(.Object, samples=samples, loci=loci, ...)
          })

# access and replacement methods for Present and Absent values
setMethod("Present", "genbinary", function(object) return(object@Present))
setReplaceMethod("Present", signature(object="genbinary"), function(object, value){
    # potential errors
    if(length(value) != 1) stop("Value must be vector of length 1.")
    if(value == object@Absent) stop("Value must be different from Absent.")
    if(value == object@Missing) stop("Value must be different from Missing.")

    # replace the value in the genotype matrix
    object@Genotypes[object@Genotypes == object@Present] <- value

    # replace the value in the Present slot
    object@Present <- value

    return(object)
})
setMethod("Absent", "genbinary", function(object) return(object@Absent))
setReplaceMethod("Absent", signature(object="genbinary"), function(object, value){
    # potential errors
    if(length(value) != 1) stop("Value must be vector of length 1")
    if(value == object@Present) stop("Value must be different from Present.")
    if(value == object@Missing) stop("Value must be different from Missing.")

    # replace the value in the genotype matrix
    object@Genotypes[object@Genotypes == object@Absent] <- value

    # replace the value in the Absent slot
    object@Absent <- value

    return(object)
})

# replacement method for Missing
setReplaceMethod("Missing", signature(object="genbinary"), function(object, value){
    # potential errors
    if(length(value) != 1) stop("Value must be vector of length 1.")
    if(value == object@Present) stop("Value must be different from Present")
    if(value == object@Absent) stop("Value must be different from Absent")

    # replace the value in the genotype matrix
    object@Genotypes[object@Genotypes == object@Missing] <- value

    # go to the gendata method to fill in Missing slot
    callNextMethod(object, value)
})

# methods for retrieving and replacing genotypes
setMethod("Genotypes", "genbinary", function(object, samples, loci){
    # set up vector to index columns for these loci
    loccolumns <- integer(0)
    for(L in loci){
        loccolumns <- c(loccolumns, grep(paste("^", L,"\\.",sep=""),
                                         dimnames(object@Genotypes)[[2]]))
    }

    # return subset of the matrix
    return(object@Genotypes[samples, loccolumns, drop=FALSE])
})
setMethod("Genotype", "genbinary", function(object, sample, locus){
    return(Genotypes(object, sample, locus))
})

setReplaceMethod("Genotypes", "genbinary",
                 function(object, samples, loci, value){
                     # make data frame of loci and alleles by column
                     tempcol <- strsplit(dimnames(value)[[2]], split=".", fixed=TRUE)
                     colinfo <- data.frame(Locus=rep("blank", length(tempcol)),
                                           Allele=rep("blank", length(tempcol)),
                                           stringsAsFactors=FALSE)
                     for(i in 1:length(tempcol)){
                         if(length(tempcol[[i]]) != 2)
                             stop("Column names should be of form \"locus.allele\".")
                         colinfo[i,]<-tempcol[[i]]
                     }
                     # check that loci are in the object
                     if(FALSE %in% (unique(colinfo$Locus) %in% loci))
                         stop("Locus names in column names should match loci in object.")
                     # check that length of first dimension matches # samples
                     if(dim(value)[1] != length(samples))
                         stop("Number of rows should match number of samples to replace.")
                     # check that only valid symbols are used
                     if(FALSE %in% (unique(value) %in% c(Absent(object), Present(object),
                                                         Missing(object))))
                         stop("Symbols in matrix must match Absent(object), Present(object), or Missing(object).")

                     # get samples that aren't being filled in (for indexing blank
                     # spaces below)
                     samblank <- Samples(object)[!Samples(object) %in% samples]

                     # fill in the values
                     for(i in 1:dim(value)[2]){
                         # get the name of this column
                         colname <- dimnames(value)[[2]][i]
                     # for allele columns that already exist, replace
                         if(colname %in% dimnames(object@Genotypes)[[2]]){
                             object@Genotypes[samples, colname] <-
                                 value[,i]
                         } else {
                     # for allele columns that don't exist, add them
                             object@Genotypes <- cbind(object@Genotypes,
                                  matrix(nrow=length(Samples(object)), ncol=1,
                                         dimnames=list(NULL, colname)))
                             object@Genotypes[samples, colname] <-
                                 value[,i]
                     # fill any new blank spaces with absent or missing
                             for(s in samblank){
                                 g <- Genotype(object, s, colinfo[i,"Locus"])
                                 if(all(is.na(g)) || all(g[!is.na(g)]==Missing(object))){
                                     object@Genotypes[s,colname] <- Missing(object)
                                 } else {
                                     object@Genotypes[s,colname] <- Absent(object)
                                 }
                             }
                         }
                     }

                     return(object)
                 })

# replacement methods for samples and loci
setReplaceMethod("Samples", "genbinary", function(object, value){
    # change the names in the Genotypes slot
    dimnames(object@Genotypes)[[1]] <- value

    # go to the gendata method to change sample names in other slots.
    callNextMethod(object, value)
})
setReplaceMethod("Loci", "genbinary", function(object, value){
    # get the original locus names and systematically replace them
    oldloci <- Loci(object)
    value <- fixloci(value) # correct locus names if necessary
    for(i in 1:length(oldloci)){
        dimnames(object@Genotypes)[[2]] <- gsub(paste(oldloci[i],".",sep=""),
                                                paste(value[i],".",sep=""),
                                                dimnames(object@Genotypes)[[2]],
                                                fixed=TRUE)
    }
    # go to gendata method to change locus names in Usatnts
    callNextMethod(object, value)
})

# method to determine if a genotype is missing
setMethod("isMissing", "genbinary", function(object, samples, loci){
    if(length(samples) == 1 && length(loci) == 1){
        return(TRUE %in% (Genotypes(object, samples, loci) == Missing(object)))
    } else {
        result <- array(FALSE, dim=c(length(samples), length(loci)),
                        dimnames=list(samples, loci))
        for(s in samples){
            for(L in loci){
                result[s,L] <- isMissing(object, s, L)
            }
        }
        return(result)
    }
})

# summary method for genbinary
setMethod("summary", "genbinary", function(object){
    # tally missing genotypes
    nummissing <- 0
    for(L in Loci(object)){
        for(s in Samples(object)){
            if(isMissing(object, s, L)) nummissing <- nummissing + 1
        }
    }
    # print stuff about dataset
    cat("Binary presence/absence dataset.",
        Description(object),
        paste("Number of missing genotypes:", nummissing), sep="\n")
    # go to the summary method for gendata
    callNextMethod(object)
})

# method for viewing genotypes
setMethod("viewGenotypes", "genbinary", function(object, samples, loci){
    for(L in loci){
        print(Genotypes(object, samples, L))
        cat("\n", sep="")
    }
})

# method for editing genotypes
setMethod("editGenotypes", "genbinary", function(object, maxalleles, samples, loci){
    Genotypes(object, samples) <- edit(Genotypes(object, samples, loci))
    return(object)
})

# methods for deleting samples and loci
setMethod("deleteSamples", "genbinary", function(object, samples){
    samtouse <- Samples(object)[!Samples(object) %in% samples]
    object@Genotypes <- object@Genotypes[samtouse,]
    callNextMethod(object, samples)
})
setMethod("deleteLoci", "genbinary", function(object, loci){
    loccolumns <- integer(0)
    for(L in loci){
        loccolumns <- c(loccolumns, grep(paste("^", L,"\\.",sep=""),
                                         dimnames(object@Genotypes)[[2]]))
    }
    object@Genotypes <- object@Genotypes[,-loccolumns]
    callNextMethod(object, loci)
})

# subscripting methods
setMethod("[", signature(x = "genbinary", i = "ANY", j= "ANY"), function(x, i, j){
    x@Genotypes <- Genotypes(x, i, j)
    x@PopInfo <- x@PopInfo[i]
    if(is(x@Ploidies, "ploidymatrix")){
      x@Ploidies@pld <- x@Ploidies@pld[i,j, drop=FALSE]
    }
    if(is(x@Ploidies, "ploidysample")){
      x@Ploidies@pld <- x@Ploidies@pld[i]
    }
    if(is(x@Ploidies, "ploidylocus")){
      x@Ploidies@pld <- x@Ploidies@pld[j]
    }
    x@Usatnts <- x@Usatnts[j]
    return(x)
})

# method for estimating ploidy
setMethod("estimatePloidy", "genbinary",
          function(object, extrainfo, samples, loci){
            # get Ploidies into ploidysample format if necessary
              if(!is(object@Ploidies, "ploidysample")){
                object <- reformatPloidies(object, output="sample", erase=TRUE)
              }
              # set up array to contain the maximum and average number of alleles
              ploidyinfo <- array(dim=c(length(samples),2),
                                  dimnames=list(samples,
                                  c("max.alleles","mean.alleles")))
              # Convert Present and Absent to 1's and 0's so
              # simple addition can be done
              objectA <- object
              Missing(objectA) <- as.integer(-9)
              Absent(objectA) <- as.integer(0)
              Present(objectA) <- as.integer(1)

              # fill the array
              for(s in samples){
                  loctouse <- loci[!isMissing(objectA, s, loci)]
                  if(length(loctouse) ==0){
                      ploidyinfo[s,] <- c(NA, NA)
                  } else {
                      numalleles <- rep(0, times=length(loctouse))
                      names(numalleles) <- loctouse
                      for(L in loctouse){
                          numalleles[L] <- sum(Genotypes(objectA,s,L))
                      }
                      ploidyinfo[s, "max.alleles"] <- max(numalleles)
                      ploidyinfo[s, "mean.alleles"] <- mean(numalleles)
                  }
              }

              ## Build data frame
              ploidyinfo <- as.data.frame(ploidyinfo)
              # Add in extrainfo, allowing for named or unnamed vectors, arrays, and
              # data frames.
              if(!missing(extrainfo)){
                  if(is.vector(extrainfo)){
                      if(!is.null(names(extrainfo))) extrainfo <- extrainfo[samples]
                      ploidyinfo <- data.frame(extrainfo = extrainfo, ploidyinfo)
                  } else {
                      if(is.array(extrainfo) | is.data.frame(extrainfo)){
                          extrainfo <- as.data.frame(extrainfo)
                          if(!identical(row.names(extrainfo),
                                        as.character(1:dim(extrainfo)[1])))
                              extrainfo <- extrainfo[samples,]
                          ploidyinfo <- data.frame(extrainfo, ploidyinfo)
                      }
                  }
              }
              # Add in PopInfo from the object
              if(length(unique(PopInfo(object))) >1)
                  ploidyinfo <- data.frame(population = PopInfo(object)[samples],
                                           ploidyinfo)
              # Add in Ploidies from the object
              if(length(unique(Ploidies(object))) >1)
                  ploidyinfo <- data.frame(ploidyinfo,
                                           current.ploidy = Ploidies(object)[samples])
              # Make a column for new ploidies
              ploidyinfo <- data.frame(ploidyinfo, new.ploidy = ploidyinfo$max.alleles)

              # let user edit data frame
              cat("Edit the new.ploidy values, then close the data editor window.",
                  sep="\n")
              ploidyinfo <- edit(ploidyinfo)

              # write new.ploidy to object@Ploidies
              Ploidies(object)[samples] <- ploidyinfo$new.ploidy

              # return the object with edited ploidies
              return(object)
          })

# method for merging objects
setMethod("merge", signature(x = "genbinary", y="genbinary"),
          function(x, y, objectm, samples, loci, overwrite){
    # set up new genbinary object if this wasn't called from the method of
    # a subclass
    if(missing(objectm)){
        if(missing(samples)) samples <- unique(c(Samples(x), Samples(y)))
        if(missing(loci)) loci <- unique(c(Loci(x), Loci(y)))
        objectm <- new("genbinary", samples, loci)
    }

    # determine what to do with conflicting data
    if(missing(overwrite)) overwrite <- "empty"
    owerror <- (!overwrite %in% c("x", "y"))
    if(overwrite=="x"){
        objectB <- x
        objectT <- y
    } else {
        objectB <- y
        objectT <- x
    }

    # make sure Missing, Present, and Absent are the same
    if(Missing(x) != Missing(y)){
        if(owerror) stop("Different missing data symbols used; set overwrite to \"x\" or \"y\".")
        Missing(objectB) <- Missing(objectT)
    }
    Missing(objectm) <- Missing(objectT)
    if(Present(x) != Present(y)){
        if(owerror) stop("Different values for Present; set overwrite to \"x\" or \"y\".")
        Present(objectB) <- Present(objectT)
    }
    Present(objectm) <- Present(objectT)
    if(Absent(x) != Absent(y)){
        if(owerror) stop("Different values for Absent; set overwrite to \"x\" or \"y\".")
        Absent(objectB) <- Absent(objectT)
    }
    Absent(objectm) <- Absent(objectT)

    # merge genotypes into objectm
    # put in objectB first
    samtofillB <- Samples(objectm)[Samples(objectm) %in% Samples(objectB)]
    loctofillB <- Loci(objectm)[Loci(objectm) %in% Loci(objectB)]
    Genotypes(objectm, samtofillB, loctofillB) <-
    Genotypes(objectB, samtofillB, loctofillB)
    # then put in objectT
    samtofillT <- Samples(objectm)[Samples(objectm) %in% Samples(objectT)]
    loctofillT <- Loci(objectm)[Loci(objectm) %in% Loci(objectT)]
    Genotypes(objectm, samtofillT, loctofillT) <-
    Genotypes(objectT, samtofillT, loctofillT)
    # some missing data symbols may have been introduced; get rid of them.
    # (this can occur if B has alleles that T doesn't)
    for(s in samtofillT){
        for(L in loctofillT){
            if(isMissing(objectm, s, L) && !isMissing(objectT, s, L)){
                Ballele <- names(Genotypes(objectm,s,L))[Genotypes(objectm,s,L)
                                                         == Missing(objectm)]
                objectm@Genotypes[s,Ballele] <- Absent(objectm)
            }
        }
    }
    # check to see if objectB overwritten if there wasn't supposed to be
    # conflicting data
    if(owerror && !identical(Genotypes(objectm, samtofillB, loctofillB),
                             Genotypes(objectB, samtofillB, loctofillB)))
       stop("Conflicting genotype data.  Set overwrite to \"x\" or \"y\".")

    # call gendata method to merge popinfo etc
    callNextMethod(x, y, objectm, overwrite=overwrite)
})

# method for finding all unique genotypes for a locus, and indexing the genotype for each individual
setMethod("genIndex", signature(object = "genambig", locus = "ANY"), 
          function(object, locus){
            # list of genotypes for this locus, put alleles in ascending order for each genotype
            thesegen <- lapply(Genotypes(object, loci = locus), sort)
            uniquegen <- list() # to hold each unique genotype
            nGeno <- 0 # number of unique genotypes found
            nind <- length(thesegen) # number of individuals
            genindex <- integer(nind) # to hold index identifying genotype for each individual
            names(genindex) <- Samples(object)
            
            for(i in 1:nind){ # loop through individuals
              thisgen <- thesegen[[i]] # genotype of this sample
              genfind <- sapply(uniquegen, function(x) identical(thisgen, x)) # this this genotype in the list already, and where
              if(!any(genfind)){ # new genotype
                nGeno <- nGeno + 1
                uniquegen[[nGeno]] <- thisgen
                genindex[i] <- nGeno
              } else {           # existing genotype
                stopifnot(sum(genfind) == 1) # should only match once
                genindex[i] <- which(genfind)
              }
            }
            
            return(list(uniquegen = uniquegen, genindex = genindex))
          })

# method for the array-list of genotypeProbs results produced by meandistance.matrix2
setMethod("genIndex", signature(object = "array", locus = "ANY"),
          function(object, locus){
            stopifnot(!missing(locus))
            thesegen <- lapply(object[,locus], function(x) x$genotypes)
            uniquegen <- list()
            nGeno <- as.integer(0)
            nind <- length(thesegen)
            genindex <- list() # each item in genindex will be a vector of genotype indices for all possible genotypes
            length(genindex) <- nind
            names(genindex) <- dimnames(object)[[1]]
            
            for(i in 1:nind){
              ng <- dim(thesegen[[i]])[1] # number of genotypes
              genindex[[i]] <- integer(ng)
              for(j in 1:ng){
                thisgen <- sort(thesegen[[i]][j,])
                genfind <- sapply(uniquegen, function(x) identical(thisgen, x))
                if(!any(genfind)){ # new genotype
                  nGeno <- nGeno + 1
                  uniquegen[[nGeno]] <- thisgen
                  genindex[[i]][j] <- nGeno
                } else {           # existing genotype
                  stopifnot(sum(genfind) == 1) # should only match once
                  genindex[[i]][j] <- which(genfind)
                }
              }
            }
            
            return(list(uniquegen = uniquegen, genindex = genindex))
          })