Skip to content

Instantly share code, notes, and snippets.

@monogenea
Created October 7, 2019 18:48
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save monogenea/73196c1700aa75e2c830a8f8cfbcd6a8 to your computer and use it in GitHub Desktop.
Save monogenea/73196c1700aa75e2c830a8f8cfbcd6a8 to your computer and use it in GitHub Desktop.
# 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