Skip to content

Commit

Permalink
add read_color_by
Browse files Browse the repository at this point in the history
  • Loading branch information
@danielmsk
danielmsk committed Aug 24, 2020
1 parent e66183b commit e67b9ed
Show file tree
Hide file tree
Showing 91 changed files with 1,051 additions and 449 deletions.
Binary file modified .DS_Store
View file
Binary file not shown.
47 changes: 0 additions & 47 deletions build/lib/bamsnap/_options.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -170,53 +170,6 @@ def get_options():
parser.add_argument('-' + a1['param'], default=a1['default'], help=a1['help'], nargs=a1['nargs'], type=valuetype)

parser.add_argument('-silence', dest='silence', action="store_true", default=False, help='don\'t print any log.')
# parser.add_argument('-debug', dest='debug', action="store_true", default=False, help='turn on the debugging mode')
# parser.add_argument('-bam', dest='bam', default=[], help='bam file', nargs="*")
# parser.add_argument('-bamlist', dest='bamlist', default=None, help='list file with bam file paths')
# parser.add_argument('-title', dest='title', default=[], help='title (name) of bam file(s)', nargs="*")
# parser.add_argument('-pos', dest='pos', default=[], nargs="*",
# help='genomic position (ex. 1:816687-818057, 12:7462545)')
# parser.add_argument('-vcf', dest='vcf', default=None,
# help='list file with genomic positions with VCF format')
# parser.add_argument('-bed', dest='bed', default=None,
# help='list file with genomic positions with BED format')
# parser.add_argument('-out', dest='out', default='', help='title of output file')
# parser.add_argument('-out_type', dest='out_type',
# choices=['png', 'jpg'], default='png', help='output file type')
# parser.add_argument('-conf', dest='conf', default="", help='configuration file')
# parser.add_argument('-ref', dest='ref', default="", help='Reference sequence fasta file (ex. hg19.fa)')
# parser.add_argument('-width', dest='width', default=1000, type=int, help='image size : width (unit:px)')
# parser.add_argument('-height', dest='height', default=None, help='image size : height (unit:px)')
# parser.add_argument('-read_thickness', dest='read_thickness', default=5, type=int, help='read width (unit:px)')
# parser.add_argument('-read_gap_h', dest='read_gap_h', default=2, type=int, help='read gap height (unit:px)')
# parser.add_argument('-read_gap_w', dest='read_gap_w', default=2, type=int, help='read gap width (unit:px)')
# parser.add_argument('-draw', dest='draw', default="coverage,base,heatmap,read",
# help='plot (default: -draw coverage,heatmap,read )')
# parser.add_argument('-margin', dest='margin', default=200, type=int, help='genomic margin size')
# parser.add_argument('-center_line', dest='center_line', default=False, action="store_true", help='draw center line')
# parser.add_argument('-no_target_line', dest='no_target_line', default=False, action="store_true", help='do not draw target line')

# coverage plot option
# parser.add_argument('-draw_coverage_plot', dest='draw_coverage_plot',action="store_true", default=False, help='draw coverage plot')
# parser.add_argument('-base_height', dest='base_height', default=40, type=int, help='base plot height')

# parser.add_argument('-coverage_height', dest='coverage_height', default=40, type=int, help='coverage plot height')
# parser.add_argument('-coverage_vaf', dest='coverage_vaf', default=10, type=float,
# help='coverage variant allele fraction threshold (unit:%%)')

# coverage heatmap height
# parser.add_argument('-draw_heatmap', dest='draw_heatmap', action="store_true", default=False, help='draw heatmap')
# parser.add_argument('-heatmap_height', dest='heatmap_height', default=5, type=int, help='coverage heatmap height')

# parser.add_argument('-no_label', dest='no_label', default=False, action="store_true", help='do not draw label.')

# parser.add_argument('-no_coordinate', dest='no_coordinate', default=False, action="store_true", help='do not draw coordinate.')
# parser.add_argument('-coordinate_height', dest='coordinate_height', default=20, type=int, help='coordinate height')

# parser.add_argument('-draw_read', dest='draw_read_plot',action="store_true", default=False, help='draw read plot')
# parser.add_argument('-no_geneplot', dest='no_geneplot', action="store_true", default=False, help='do not draw gene structure')
# parser.add_argument('-geneplot_height', dest='geneplot_height', default=50, type=int, help='gene plot height')
# parser.add_argument('-merged_image', dest='merged_image', default=False, action="store_true", help='draw a merged plot')

if len(sys.argv) == 1 or (len(sys.argv) == 2 and sys.argv[1][0] != '-'):
sys.argv.append('-h')
Expand Down
25 changes: 19 additions & 6 deletions build/lib/bamsnap/bamsnap.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -252,6 +252,7 @@ def get_bamplot_image(self, bam, pos1, image_w, xscale, refseq):
rset.read_gap_h = self.opt['read_gap_height']
rset.read_thickness = self.opt['read_thickness']
rset.coverage_vaf = self.opt['coverage_vaf']
rset.opt = self.opt
rset.xscale = xscale
rset.calculate_readmap(is_strand_group=True)

Expand All @@ -271,25 +272,29 @@ def get_bamplot_image(self, bam, pos1, image_w, xscale, refseq):

if plot1 == "coverage":
covplot = CoveragePlot(rset, xscale, self.opt['coverage_vaf'])
covplot.coverage_color = self.opt['coverage_color']
covplot.font = self.get_font(self.opt['coverage_fontsize'])
ia_sub = covplot.get_image(image_w, self.opt['coverage_height'], self.opt['coverage_bgcolor'])
im = self.append_image(im, ia_sub)

if plot1 == "read":
if self.opt['read_group'] == "":
h_all = rset.get_estimated_height('all')
ia_sub = rset.get_image(image_w, h_all, 'all', self.opt['read_color'], self.opt['read_bgcolor'])
ia_sub = rset.get_image(
image_w, h_all, 'all', self.opt['read_color'], self.opt['read_bgcolor'], self.opt['read_color_by'])
im = self.append_image(im, ia_sub)

elif self.opt['read_group'] == "strand":
# self.opt['read_bgcolor'] = "F0F000"
h_pos = rset.get_estimated_height('pos_strand')
ia_sub = rset.get_image(image_w, h_pos, 'pos_strand', self.opt['read_pos_color'], self.opt['read_bgcolor'])
ia_sub = rset.get_image(image_w, h_pos, 'pos_strand',
self.opt['read_pos_color'], self.opt['read_bgcolor'], self.opt['read_color_by'])
im = self.append_image(im, ia_sub)

# self.opt['read_bgcolor'] = "00F0F0"
h_neg = rset.get_estimated_height('neg_strand')
ia_sub = rset.get_image(image_w, h_neg, 'neg_strand', self.opt['read_neg_color'], self.opt['read_bgcolor'])
ia_sub = rset.get_image(image_w, h_neg, 'neg_strand',
self.opt['read_neg_color'], self.opt['read_bgcolor'], self.opt['read_color_by'])
im = self.append_image(im, ia_sub)

if plot1 == "coordinates":
Expand Down Expand Up @@ -350,6 +355,8 @@ def drawplot_bamlist(self, pos1, image_w, bamlist, xscale, refseq):

if plot1 == "bamplot":
for bidx, bam in enumerate(bamlist):
if (bidx + 1) % 10 == 0:
self.opt['log'].info("..processing " + bam.filename + " (" + str(bidx + 1) + ")")
ia_sub = self.get_bamplot_image(bam, pos1, image_w, xscale, refseq)
ia = self.append_image(ia, ia_sub)

Expand All @@ -362,11 +369,17 @@ def drawplot_bamlist(self, pos1, image_w, bamlist, xscale, refseq):

if plot1 == "base":
ia = self.append_baseplot_image(ia, pos1, image_w, xscale, refseq)


if self.opt['grid'] > 0:
drawA = ImageDraw.Draw(ia)
for posi in range(pos1['g_spos'], pos1['g_epos'] + 1, self.opt['grid']):
x1 = xscale.xmap[posi]['spos'] - 1
drawA.line([(x1, 0), (x1, ia.height)], fill=COLOR['GRID_COLOR'], width=1)

if self.opt['center_line']:
drawA = ImageDraw.Draw(ia)
x1 = ((self.opt['g_epos']-self.opt['g_spos'])/2)*self.scale_x
for h1 in range(0, image_h, 7):
for h1 in range(0, ia.height, 7):
drawA.line([(x1, h1), (x1, h1+2)], fill=COLOR['CENTER_LINE'], width=1)

if not self.opt['no_target_line']:
Expand Down Expand Up @@ -412,7 +425,7 @@ def generate_zipfile(self):

def start_process_drawplot(self, image_w, bamlist):
for tno in range(self.opt['process']):
self.process[tno] = mp.Process(target=self.start_process_drawplot_bamlist, args=(image_w, bamlist,self.split_poslist[tno]), name='proc ' + str(tno+1))
self.process[tno] = mp.Process(target=self.start_process_drawplot_bamlist, args=(image_w, bamlist, self.split_poslist[tno]), name='proc ' + str(tno+1))
self.process[tno].start()

def start_process_drawplot_bamlist(self, image_w, bamlist, poslist):
Expand Down
5 changes: 4 additions & 1 deletion build/lib/bamsnap/conf.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,7 @@
COLOR['GENE_POS'] = ImageColor.getrgb("#6961ff")
# COLOR['GENE_NEG'] = ImageColor.getrgb("#e87489")
COLOR['GENE_NEG'] = (232, 116, 137, 255)
COLOR['GRID_COLOR'] = ImageColor.getrgb("#EEEEEE")


# A: "rgb( 0, 200, 0)",
Expand All @@ -46,4 +47,6 @@

IMAGE_MARGIN_BOTTOM = 20

GENE_ANNOT_FILE = "data/Homo_sapiens.GRCh38.99.bed.gz"
GENE_ANNOT_FILE = "Homo_sapiens.#REFSEQVERSION#.bed.gz"

REFER_SEQ_VERSION = {'GRCh37': 'GRCh37', 'hg19':'GRCh37', 'GRCh38':'GRCh38', 'hg38':'GRCh38'}
45 changes: 25 additions & 20 deletions build/lib/bamsnap/drawread.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,36 +16,35 @@ class DrawRead():
x1 = 0
x2 = 0
y1 = 0
opt = {}

def __init__(self, a): # a:alignment
#self.a = a
self.mapq = a.mapq
self.id = ''
self.is_proper_pair = a.is_proper_pair
self.is_reverse = a.is_reverse

self.base_qual = a.query_alignment_qualities # two qual scores : a.query_alignment_qualities and a.query_qualities
self.g_positions = self.base_plus_1(a.positions)

self.g_spos = self.g_positions[0]
self.g_epos = self.g_positions[-1]

self.readseq = a.query_alignment_sequence


self.cigar = a.cigartuples
self.set_cigar()

self.readseqpos = {}
self.set_readseqpos()
# print(self.readseq)
# print(self.g_positions)
# print(self.readseqpos)
# print(self.ins_pos_map)

self.g_len = self.g_epos - self.g_spos + 1
self.set_color()
pass
self.reference_name = a.reference_name
self.mate_reference_name = a.next_reference_name
self.has_interchrom_mate = self.reference_name != self.mate_reference_name
# print(a.mate_is_unmapped, a.mate_is_reverse, a.next_reference_name)
# print(a.is_paired)
# print(a.reference_name, a.next_reference_name)
# print(a.mate_is_unmapped, a.mate_is_reverse)
# print(a.mate())
# if a.reference_name != a.next_reference_name:
# print(a.reference_name, a.next_reference_name)

def base_plus_1(self, poslist):
for i in range(len(poslist)):
Expand All @@ -59,11 +58,9 @@ def set_readseqpos(self):
try:
self.ins_pos_map[self.g_positions[i]]
j += self.ins_pos_map[self.g_positions[i]]
# print(self.g_positions[i], self.ins_pos_map[self.g_positions[i]], self.readseq[j])
except KeyError:
pass
self.readseqpos[self.g_positions[i]] = self.readseq[j]
# print(self.g_positions[i], self.readseqpos[self.g_positions[i]], i, j)
j += 1

def set_cigar(self):
Expand All @@ -73,7 +70,6 @@ def set_cigar(self):
gpos = self.g_spos - 1
for cg in self.cigar:
gpos += cg[1]
# print(gpos, cg, self.cigar, self.g_spos)
if cg[0] == 2:
self.has_del = True
if cg[0] == 1:
Expand Down Expand Up @@ -102,14 +98,22 @@ def read_arrowhead_thickness(self, read_thickness):
t1 = 13
return t1

def draw(self, dr, col1="C8C8C8"):
def get_readcolor_by_interchrom(self, default_color):
if self.has_interchrom_mate:
mate_chrom = self.mate_reference_name.replace('chr','')
print(self.opt)
rst_color = self.opt['read_color_interchrom_chr'] + mate_chrom
else:
rst_color = default_color
return rst_color

def draw(self, dr, col1="C8C8C8", readcolorby=""):
if not self.flag_draw:
x1 = self.xscale.xmap[self.g_spos]['spos']
x2 = self.xscale.xmap[self.g_epos]['epos']
y1 = int(self.get_scaled_y(self.yidx))

xy = []

if self.read_thickness > 1:
raht = self.read_arrowhead_thickness(self.read_thickness)
if not self.is_reverse:
Expand All @@ -126,18 +130,20 @@ def draw(self, dr, col1="C8C8C8"):
xy.append((x1, y1 + int(self.read_thickness / 2) ))
xy.append((x1 - raht, y1))

if readcolorby == "insert_size":
col1 = self.get_readcolor_by_interchrom(col1)

if self.mapq == 0:
dr.polygon(xy, fill=getrgb(col1, 60), outline=self.outline_color)
else:
dr.polygon(xy, fill=getrgb(col1))
else:
dr.line([(x1, y1), (x2, y1)], fill=col1, width=self.read_thickness)
dr.line([(x1, y1), (x2, y1)], fill=getrgb(col1), width=self.read_thickness)

self.draw_cigar(dr, y1)
self.draw_variants(dr, y1)
self.flag_draw = True


def draw_cigar(self, dr, y1):
if self.has_ins or self.has_del:
xidx = 0
Expand Down Expand Up @@ -211,7 +217,6 @@ def draw_variants(self, dr, y1):
no_variant += 1
i += 1

# if no_variant < 5:
if True:
i = 0
# if self.has_ins:
Expand Down
Loading

0 comments on commit e67b9ed

Please sign in to comment.