Skip to content

Commit e6049b2

Browse files
committed
Merge branch 'pr/100'
# Conflicts: # AppImageBuilder.yml # web/auto-imports.d.ts # web/components.d.ts # web/package.json # web/pnpm-lock.yaml
2 parents 6be5da2 + 270000f commit e6049b2

5 files changed

Lines changed: 211 additions & 74 deletions

File tree

docs/command.md

Lines changed: 60 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -661,7 +661,63 @@ trackplot \
661661

662662
![](imgs/cmd/intron_scale.png)
663663

664-
### 7. Visualize coverage by cpm or rpkm
664+
665+
### 7. Fixed intron size
666+
667+
668+
There are two BAM files. Here is the output with introns fixed at 250 nt, exons unharmed (1:1), without any filtering of low-count junctions:
669+
670+
![](imgs/cmd/PTBP3_mod_specific_transcipts.png)
671+
672+
With --log 10, we can see all the low-count junctions:
673+
674+
![](imgs/cmd/PTBP3_mod_specific_transcipts_log10_edit.png)
675+
676+
If we use Trackplot's scaling features, we can get reasonably close (exons=20, introns=0.10):
677+
678+
![](imgs/cmd/PTBP3_mod_specific_transcipts_introns_and_exons_scaled.png)
679+
680+
Here is the equivalent plot with Trackplot (note the gaps remaining from excluded transcripts):
681+
682+
![](imgs/cmd/PTBP3_trackplot_specific_trans_all_junctions.exons.20.introns.0.10.png)
683+
684+
With --density-by-strand:
685+
686+
![](imgs/cmd/PTBP3_trackplot_all_junctions_specific_transcripts_by_strand_edit.png)
687+
688+
689+
To be clear, the modifications described here were a small step following a giant leap forward. For instance, here is how the data looks in IGV:
690+
691+
![](imgs/cmd/ptbp3_igv_2.png)
692+
693+
694+
The switch to positive-only Y values is hard-coded, along with other preferences (transparency etc.).
695+
696+
To replicate the figures above, use:
697+
698+
```bash
699+
trackplot \
700+
--intron-scale 250 \ # To assigned a fixed length to introns, use an intron scale greater than 1
701+
-e 9:112225716-112333664 \
702+
-r ptbp3_ext.gtf.gz \
703+
--show-junction-num \
704+
--density bam_files.tsv \
705+
--raster \
706+
--height 2 \
707+
--width 20 \
708+
--font-size 14 \
709+
--transcripts-to-show PTBP3-201,PTBP3-202,PTBP3-203,PTBP3-204,PTBP3-205,PTBP3-206,PTBP3-207,PTBP3-208 \
710+
-o PTBP3_mod_specific_transcipts.log10.png
711+
```
712+
713+
To scale using exising Trackplot features, use:
714+
715+
```bash
716+
--exon-scale 20 \
717+
--intron-scale 0.10 \
718+
```
719+
720+
### 8. Visualize coverage by cpm or rpkm
665721

666722
We also support to visualize the coverage by normalized values.
667723

@@ -1119,4 +1175,6 @@ trackplot \
11191175
--line example/line_list.tsv
11201176
```
11211177

1122-
![](imgs/cmd/additional.png)
1178+
![](imgs/cmd/additional.png)
1179+
1180+

trackplot/file/Annotation.py

Lines changed: 23 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -2,8 +2,10 @@
22
# -*- coding:utf-8 -*-
33
u"""
44
Created by ygidtu@gmail.com at 2020.05.07
5-
65
This scripts contains the class handle the reference file
6+
7+
Modified by AD 2025/01/12 to include only specified transcripts, otherwise their exon coordinates still occupy the alignments.
8+
79
"""
810
import glob
911
import gzip
@@ -349,16 +351,15 @@ def index_gtf(cls, input_gtf):
349351

350352
return output_gtf
351353

352-
def __load_gtf__(self):
353-
354+
def __load_gtf__(self, transcripts_to_show: list[str]|None = None):
354355
u"""
355356
Load transcripts inside of region from gtf file
356357
:param region: target region
357358
:return: list of Transcript
358-
"""
359+
"""
360+
# AD - passing transcripts_to_show
359361
transcripts = {}
360362
exons = {}
361-
362363
for rec in Reader.read_gtf(self.path, self.region):
363364
start = max(rec.start, self.region.start)
364365
end = min(rec.end, self.region.end)
@@ -369,6 +370,13 @@ def __load_gtf__(self):
369370
break
370371

371372
if re.search(r"(rna|transcript|cds)", rec.feature, re.I):
373+
374+
if transcripts_to_show:
375+
_name = rec.transcript_name if "transcript_name" in rec.attributes else rec.transcript_id
376+
if _name not in transcripts_to_show:
377+
logger.info(f"Skipping transcript {_name}")
378+
continue
379+
372380
if rec.transcript_id not in transcripts.keys():
373381
transcripts[rec.transcript_id] = Transcript(
374382
chromosome=rec.contig,
@@ -380,7 +388,8 @@ def __load_gtf__(self):
380388
gene=rec.gene_name if "gene_name" in rec.attributes else "",
381389
transcript=rec.transcript_name if "transcript_name" in rec.attributes else "",
382390
exons=[]
383-
)
391+
)
392+
384393
elif re.search(r"(exon)", rec.feature, re.I):
385394
if rec.transcript_id not in exons.keys():
386395
exons[rec.transcript_id] = []
@@ -581,9 +590,16 @@ def load(self,
581590
assert isinstance(region, GenomicLoci), "region should be a GenomicLoci object"
582591
if transcripts is None:
583592
transcripts = []
593+
elif isinstance(transcripts, str):
594+
# AD when no transcripts are specified, transcripts is an empty string, not None
595+
if transcripts == "":
596+
transcripts = []
597+
else:
598+
transcripts = transcripts.split(",")
584599
self.region = region
600+
# AD - pass specified transcripts early
585601
if self.category == "gtf":
586-
self.__load_gtf__()
602+
self.__load_gtf__(transcripts)
587603
if self.add_local_domain:
588604
self.__load_local_domain__(region)
589605

trackplot/file/ReadSegments.py

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -382,7 +382,9 @@ def load_bam(self):
382382
current_ignore_num = min(
383383
[self.deletion_ignore, read.query_alignment_length * self.del_ratio_ignore])
384384
else:
385-
current_ignore_num = np.Inf
385+
current_ignore_num = np.inf
386+
# AD - AttributeError: `np.Inf` was removed in the NumPy 2.0 release. Use `np.inf` instead.
387+
386388

387389
exon_bound = []
388390
intron_bound = []

trackplot/plot.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1437,7 +1437,7 @@ def plot(self,
14371437

14381438
if output:
14391439
logger.info(f"saving fig into {output}")
1440-
fig.savefig(output, transparent=True, bbox_inches='tight')
1440+
fig.savefig(output, transparent=False, bbox_inches='tight') # AD
14411441
elif return_image:
14421442
output = io.BytesIO()
14431443
if return_image == "png":

0 commit comments

Comments
 (0)