[med-svn] [mapdamage] 01/06: New upstream version 2.0.8+dfsg

Andreas Tille tille at debian.org
Mon Jun 19 07:55:51 UTC 2017


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 <tille at 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='MSchubert at snm.ku.dk, jonsson.hakon at gmail.com',
     packages=['mapdamage'],

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/mapdamage.git



More information about the debian-med-commit mailing list