calHe <- function(fre, filename = "hets.txt"){
pops = row.names(fre)
loci = unique(as.matrix(as.data.frame(strsplit(names(fre), split = ".", fixed = TRUE), stringsAsFactors = FALSE))[1, ])
loci <- loci[loci != "Genomes"]
hets <- array(0, dim = c(length(pops),length(loci)), dimnames = list(pops, loci))
for(P in pops){
freqs <- fre[P, ]
genomes <- freqs$Genomes
names(genomes) <- row.names(freqs)
thesegenomes <- genomes
for (L in loci) {
thesefreqs <- freqs[, grep(paste("^", L, "\\.",
sep = ""), names(freqs)), drop = FALSE]
thesefreqs <- thesefreqs[, names(thesefreqs) != paste(L,
"Genomes", sep = "."), drop = FALSE]
hsByPop <- apply(as.matrix(thesefreqs), 1, function(x) 1 -
sum(x^2))
avgfreq <- unlist(lapply(thesefreqs, weighted.mean,
w = thesegenomes))
hets[P,L] <- 1 - sum(avgfreq^2)
}
}
hets<-data.frame(hets)
hets$mean <- apply(hets,1,mean)
write.table(hets, filename, sep = "\t",quote = FALSE)
return(hets)
}
原文:https://www.cnblogs.com/liulele/p/12466350.html