-
Notifications
You must be signed in to change notification settings - Fork 7
/
genind2structure.R
56 lines (50 loc) · 2.06 KB
/
genind2structure.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
# Function to export to STRUCTURE format from genind object.
# genind objects are created in the R package adegenet. The function below is an R function.
# Lindsay V. Clark, 26 July 2015
# obj: genind object
# file: file name to write
# pops: whether to include population info in the file
# Function is flexible with regards to ploidy, although genotypes are
# considered to be unambiguous.
# Missing data must be recorded as NA in obj@tab.
# example use:
# data(nancycats)
# genind2structure(nancycats, file="nancy_structure.txt", pops=TRUE)
genind2structure <- function(obj, file="", pops=FALSE){
if(!"genind" %in% class(obj)){
warning("Function was designed for genind objects.")
}
# get the max ploidy of the dataset
pl <- max(obj@ploidy)
# get the number of individuals
S <- adegenet::nInd(obj)
# column of individual names to write; set up data.frame
tab <- data.frame(ind=rep(adegenet::indNames(obj), each=pl))
# column of pop ids to write
if(pops){
popnums <- 1:adegenet::nPop(obj)
names(popnums) <- as.character(unique(adegenet::pop(obj)))
popcol <- rep(popnums[as.character(adegenet::pop(obj))], each=pl)
tab <- cbind(tab, data.frame(pop=popcol))
}
loci <- adegenet::locNames(obj)
# add columns for genotypes
tab <- cbind(tab, matrix(-9, nrow=dim(tab)[1], ncol=adegenet::nLoc(obj),
dimnames=list(NULL,loci)))
# begin going through loci
for(L in loci){
thesegen <- obj@tab[,grep(paste("^", L, "\\.", sep=""),
dimnames(obj@tab)[[2]]),
drop = FALSE] # genotypes by locus
al <- 1:dim(thesegen)[2] # numbered alleles
for(s in 1:S){
if(all(!is.na(thesegen[s,]))){
tabrows <- (1:dim(tab)[1])[tab[[1]] == adegenet::indNames(obj)[s]] # index of rows in output to write to
tabrows <- tabrows[1:sum(thesegen[s,])] # subset if this is lower ploidy than max ploidy
tab[tabrows,L] <- rep(al, times = thesegen[s,])
}
}
}
# export table
write.table(tab, file=file, sep="\t", quote=FALSE, row.names=FALSE)
}