| # Sample call rate & heterozygosity | |
| callMat <- !is.na(genData$SNP) | |
| Sampstats <- row.summary(genData$SNP) | |
| hetExp <- callMat %*% (2 * SNPstats$MAF * (1 - SNPstats$MAF)) # Hardy-Weinberg heterozygosity (expected) | |
| hetObs <- with(Sampstats, Heterozygosity * (ncol(genData$SNP)) * Call.rate) | |
| Sampstats$hetF <- 1-(hetObs/hetExp) | |
| # Use sample call rate of 100%, het threshold of 0.1 (very stringent) | |
| het <- 0.1 # Set cutoff for inbreeding coefficient; | |
| het_call <- with(Sampstats, abs(hetF) < het & Call.rate == 1) | |
| genData$SNP <- genData$SNP[het_call,] | |
| genData$LIP <- genData$LIP[het_call,] |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment