Check that CpG positions are on chromosome

cangermueller · Mar 2, 2017 · 953cf73 · 953cf73
1 parent 70a0ec5
commit 953cf73
Showing 1 changed file with 3 additions and 1 deletion.
diff --git a/scripts/dcpg_data.py b/scripts/dcpg_data.py
@@ -96,8 +96,10 @@ def extract_seq_windows(seq, pos, wlen, seq_index=1, assert_cpg=False):
 
     for i in range(nb_win):
         p = pos[i] - seq_index
+        if p < 0 or p >= len(seq):
+            raise ValueError('Position %d not on chromosome!' % (p + seq_index))
         if seq[p:p + 2] != 'CG':
-            warnings.warn('No CpG site at position %d!' % p)
+            warnings.warn('No CpG site at position %d!' % (p + seq_index))
         win = seq[max(0, p - delta): min(len(seq), p + delta + 1)]
         if len(win) < wlen:
             win = max(0, delta - p) * 'N' + win