From da4ae79ddce7c36b48d4bb4ba314503ef9743cd7 Mon Sep 17 00:00:00 2001 From: Alexander_Kabui Date: Tue, 29 Oct 2024 15:30:21 +0300 Subject: Refactor: refactor functions to compute genetic probabilities: --- scripts/rqtl2_wrapper.R | 45 +++++++++++++++++++++++++++++---------------- 1 file changed, 29 insertions(+), 16 deletions(-) (limited to 'scripts') diff --git a/scripts/rqtl2_wrapper.R b/scripts/rqtl2_wrapper.R index 1a56207..5decd85 100644 --- a/scripts/rqtl2_wrapper.R +++ b/scripts/rqtl2_wrapper.R @@ -131,32 +131,45 @@ ind_ids_pheno(dataset) cat("Name of the founder Strains/n") founders(dataset) -# Work on computing the genetic probabilities - -analysis_type <- "single" +# Function for computing the genetic probabilities perform_genetic_pr <- function(cross, cores = 1, - error_prob = 0.002, - analysis_type = "single") { - # improve on this - if (analysis_type == "single") { - pr <- calc_genoprob(cross, - error_prob = error_prob, - quiet = FALSE, - cores = cores) - return (pr) + step=1, + map=NULL, + use_pseudomarkers=FALSE, + map_function=c("haldane", "kosambi", "c-f", "morgan"), + error_prob = 0.002 + ) { + #' Function to calculate the genetic probabilities + #' @description function to perform genetic probabilities + #' @param cores number no of cores to use Defaults to "1" + #' @param map Genetic map of markers. defaults to "NONE" + #' @param use_pseudomarkers option to insert pseudo markers in the gmap default "FALSE" + #' @param error_prob + #' @param map_function Character string indicating the map function to use to convert genetic + #' @param step for default "1" + #' @return a list of three-dimensional arrays of probabilities, individuals x genotypes x pst + + cat("Finding the genetic Probabilities\n") + if (use_pseudomarkers){ + map <- insert_pseudomarkers(cross$gmap, step=step) + return(calc_genoprob(cross, map=map, + error_prob=error_prob, map_function=map_function, + quiet=FALSE, cores=cores)) } -} - -# get the genetic probability + else { + return (calc_genoprob(cross, map=map, error_prob=error_prob, + quiet = FALSE, map_function =map_function, + cores = cores)) + }} Pr = perform_genetic_pr(dataset) cat("Summaries on the genetic probabilites \n") print(Pr) summary(Pr) -#calculate genotyping error LOD scores +#calculate genotyping error LOD scores error_lod <- calc_errorlod(dataset, Pr, quiet = FALSE, cores = 4) print(error_lod) -- cgit 1.4.1