This is an automated email from the git hooks/post-receive script. tille pushed a commit to branch master in repository mapdamage.
commit d046a9d5a9deb7c2a41abd7caade686e8dd1e8cc Author: Andreas Tille <ti...@debian.org> Date: Mon Jun 19 09:36:00 2017 +0200 New upstream version 2.0.8+dfsg --- README.md | 3 +++ mapdamage/parseoptions.py | 32 +++++++++++++++++++++++++------- mapdamage/rescale.py | 23 ++++++++++++++++------- mapdamage/version.py | 3 +-- setup.py | 2 +- 5 files changed, 46 insertions(+), 17 deletions(-) diff --git a/README.md b/README.md index f908dba..a0915b8 100644 --- a/README.md +++ b/README.md @@ -9,7 +9,10 @@ among ancient DNA sequencing reads generated by Next-Generation Sequencing platf Mapdamage is developed at the [Centre for GeoGenetics](http://geogenetics.ku.dk/) by the [Orlando Group ](http://geogenetics.ku.dk/research/research_groups/palaeomix_group/). +### Web interface +Anna Kostikova from [insideDNA](https://insidedna.me) implemented a web interface for mapDamage. +Users can adjust the cloud ressources in terms of RAM/CPU to perform their analysis. She wrote a [tutorial](https://insidedna.me/tutorials/view/Analysis-ancient-DNA-samples-using-mapDamage) explaining how to use this [web interface (sign up required)](https://insidedna.me/app#/tools/100648/) ### Citation diff --git a/mapdamage/parseoptions.py b/mapdamage/parseoptions.py index 7102cf2..bf5020f 100644 --- a/mapdamage/parseoptions.py +++ b/mapdamage/parseoptions.py @@ -128,20 +128,27 @@ def options(): help="How long sequence to use from each side [%default]", type = int, default=12, action="store") group3.add_option("--stats-only", dest="stats_only", help="Run only statistical estimation from a valid result folder", \ default=False, action="store_true") - group3.add_option("--rescale", dest="rescale", help="Rescale the quality scores in the BAM file using the output from the statistical estimation", \ + group3.add_option("--no-stats", help="Disabled statistical estimation, active by default", default=False, action="store_true") + group3.add_option("--check-R-packages", help="Check if the R modules are working", default=False, action="store_true") + parser.add_option_group(group3) + + group4 = OptionGroup(parser,"Options for rescaling of BAM files") + group4.add_option("--rescale", dest="rescale", help="Rescale the quality scores in the BAM file using the output from the statistical estimation", default=False, action="store_true") - group3.add_option("--rescale-only", dest="rescale_only", help="Run only rescaling from a valid result folder", \ + group4.add_option("--rescale-only", dest="rescale_only", help="Run only rescaling from a valid result folder", default=False, action="store_true") - group3.add_option("--rescale-out", dest="rescale_out", help="Write the rescaled BAM to this file", \ + group4.add_option("--rescale-out", dest="rescale_out", help="Write the rescaled BAM to this file", default=None, action="store") - group3.add_option("--no-stats", help="Disabled statistical estimation, active by default", default=False, action="store_true") - group3.add_option("--check-R-packages", help="Check if the R modules are working", default=False, action="store_true") + group4.add_option("--rescale-length-5p", dest="rescale_length_5p", + help="How many bases to rescale at the 5' termini; defaults to --seq-length.", type=int, action="store") + group4.add_option("--rescale-length-3p", dest="rescale_length_3p", + help="How many bases to rescale at the 5' termini; defaults to --seq-length.", type=int, action="store") + parser.add_option_group(group4) - parser.add_option_group(group3) #Parse the arguments (options, args) = parser.parse_args() - + # check python version if not check_py_version(): return None @@ -238,6 +245,17 @@ def options(): sys.stderr.write("Error, %s does not exist while plot/stats/rescale only was used\n" % options.folder) return None + if options.rescale_length_3p is None: + options.rescale_length_3p = options.seq_length + elif not (0 <= options.rescale_length_3p <= options.seq_length): + parser.error("--rescale-length-3p must be less than or equal to " + "--seq-length and greater than zero") + + if options.rescale_length_5p is None: + options.rescale_length_5p = options.seq_length + elif not (0 <= options.rescale_length_5p <= options.seq_length): + parser.error("--rescale-length-5p must be less than or equal to " + "--seq-length and greater than zero") # check if the Rscript executable is present on the system if not whereis('Rscript'): diff --git a/mapdamage/rescale.py b/mapdamage/rescale.py index 789d21e..ca71c3e 100644 --- a/mapdamage/rescale.py +++ b/mapdamage/rescale.py @@ -17,7 +17,7 @@ def phred_char_to_pval(ch): """Transforming ASCII character in the Phred scale to the error rate""" return 10**(-(float(ord(ch))-float(33))/10) -def get_corr_prob(folder): +def get_corr_prob(folder, rescale_length_5p, rescale_length_3p): """ Reads the damage probability correction file, returns a dictionary with this structure @@ -36,6 +36,12 @@ def get_corr_prob(folder): (folder, fi_handle.line_num, corr_prob[line["Position"]])) else: corr_prob[int(line["Position"])] = {'C.T':float(line["C.T"]), 'G.A':float(line["G.A"])} + + # Exclude probabilities for positions outside of user-specified region + for key in corr_prob.keys(): + if key < -rescale_length_3p or key > rescale_length_5p: + corr_prob.pop(key) + return corr_prob except csv.Error as e: sys.exit('File %s, line %d: %s' % (os.path.join(folder,"Stats_out_MCMC_correct_prob.csv"), \ @@ -269,11 +275,12 @@ def rescale_qual_read(bam, read, ref, corr_prob,subs, debug = False,direction="b read.qual = new_qual # truncate this to 5 digits number_of_rescaled_bases = float("%.5f" % number_of_rescaled_bases) - # check if the read has a MR tag - for tag in read.tags: - if tag[0] == "MR": - raise SystemExit("Read: %s already has a MR tag, can't rescale" % (read)) - read.tags = read.tags + [("MR:f",number_of_rescaled_bases)] + + if read.has_tag("MR"): + raise SystemExit("Read: %s already has a MR tag, can't rescale" % read) + + read.set_tag("MR", number_of_rescaled_bases, 'f') + return read @@ -300,7 +307,9 @@ def rescale_qual(ref, options,debug=False): else: write_mode = "wb" bam_out = pysam.Samfile(options.rescale_out,write_mode, template = bam) - corr_prob = get_corr_prob(options.folder) + corr_prob = get_corr_prob(options.folder, + rescale_length_5p=options.rescale_length_5p, + rescale_length_3p=options.rescale_length_3p) subs = initialize_subs() first_pair = True number_of_non_proper_pairs = 0 diff --git a/mapdamage/version.py b/mapdamage/version.py index 4e95d9f..5d2b8ff 100644 --- a/mapdamage/version.py +++ b/mapdamage/version.py @@ -2,5 +2,4 @@ try: from mapdamage._version import __version__ except ImportError: - __version__ = "2.0.6" - + __version__ = "2.0.8" diff --git a/setup.py b/setup.py index b166396..d0c70ca 100644 --- a/setup.py +++ b/setup.py @@ -55,7 +55,7 @@ class compileInstall(DistutilsInstall): setup( cmdclass={'install': compileInstall}, name='mapdamage', - version='2.0', + version='2.0.8', author='Aurélien Ginolhac, Mikkel Schubert, Hákon Jónsson', author_email='mschub...@snm.ku.dk, jonsson.ha...@gmail.com', packages=['mapdamage'], -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/mapdamage.git _______________________________________________ debian-med-commit mailing list debian-med-commit@lists.alioth.debian.org http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/debian-med-commit