Merge branch 'transcript_table_threshold' into 'dev'
Switch isoform table thresholding to just displaying to n rows based on coverage See merge request epi2melabs/workflow-containers/wf-transcriptomes!73
This commit is contained in:
commit
b140d3791e
@ -6,7 +6,7 @@ from collections import Counter, defaultdict, OrderedDict
|
||||
import math
|
||||
from pathlib import Path
|
||||
|
||||
from aplanat import bars, hist, lines
|
||||
from aplanat import bars, hist
|
||||
from aplanat.components import simple as scomponents
|
||||
from aplanat.components.fastcat import read_length_plot, read_quality_plot
|
||||
from aplanat.report import WFReport
|
||||
@ -594,21 +594,25 @@ def cluster_quality(cluster_qc_dir, report, sample_ids):
|
||||
section.plot(cover_panel)
|
||||
|
||||
|
||||
def transcript_table(report, df_tmaps, covr_threshold):
|
||||
"""Create searchable table of transcripts."""
|
||||
def transcript_table(report, df_tmaps, max_rows):
|
||||
"""Create searchable table of transcripts.
|
||||
|
||||
:param df_tmaps: pd.DataFrame of all gffcomapre `.tmap` files from all
|
||||
samples
|
||||
"""
|
||||
section = report.add_section()
|
||||
|
||||
# Should we put data from each sample into it's own table or have it
|
||||
# Should wße put data from each sample into it's own table or have it
|
||||
# all in single table and sample_id column? Currently it's the latter
|
||||
|
||||
# drop some columns for the big table and do some filtering
|
||||
section.markdown('''
|
||||
### Isoforms table
|
||||
|
||||
Low coverage transcripts are removed to speed up the table viewing. <br>
|
||||
Coverage threshold can be set with the parameter
|
||||
`transcript_table_cov_thresh`.
|
||||
''')
|
||||
Table interactivity can be slow if too many isoforms are loaded. <br>
|
||||
The number of isoform rows to load in this table can be set with
|
||||
`isoform_table_nrows`. It is currently set to `{}`
|
||||
'''.format(max_rows))
|
||||
|
||||
df = df_tmaps.drop(
|
||||
columns=[
|
||||
@ -619,31 +623,6 @@ def transcript_table(report, df_tmaps, covr_threshold):
|
||||
section.markdown("No transcripts found")
|
||||
return
|
||||
|
||||
df.sort_values('cov', ascending=True, inplace=True)
|
||||
counts = list(range(len(df)))
|
||||
|
||||
# Filter on coverage threshold
|
||||
df = df[df['cov'] >= covr_threshold]
|
||||
if len(df) < 200: # Min size of table should be 200. Don't filter
|
||||
df = df.sort_values('cov', ascending=False).iloc[:, 0:200]
|
||||
covr_threshold = 0
|
||||
|
||||
# Keep Isoforms with coverage > threshold
|
||||
vline_x = np.argmax(df['cov'] > covr_threshold)
|
||||
vline_y = [0, df['cov'].max()]
|
||||
|
||||
cov_plt = lines.line(
|
||||
[counts, [vline_x, vline_x]], # x-values
|
||||
[df['cov'].values.tolist(), vline_y], # y-values
|
||||
title=(
|
||||
"Read Coverage. Threshold = {}x coverage".format(
|
||||
covr_threshold)
|
||||
), x_axis_label='Num Isoforms',
|
||||
y_axis_label='Coverage',
|
||||
colors=['blue', 'red'])
|
||||
|
||||
section.plot(cov_plt)
|
||||
|
||||
# Make a column of number of isoforms in parent gene
|
||||
gb = df.groupby(['ref_gene_id', 'sample_id']).count()
|
||||
# gb = gb.set_index(['ref_gene_id', 'sample_id'])
|
||||
@ -654,7 +633,13 @@ def transcript_table(report, df_tmaps, covr_threshold):
|
||||
# Uncalssified transcritps should not be lumped togetehr
|
||||
df.loc[df.class_code == 'u', 'parent gene iso num'] = None
|
||||
|
||||
df.sort_values('parent gene iso num', inplace=True, ascending=True)
|
||||
# Keep top n rows with most coverage
|
||||
df.sort_values('cov', ascending=False, inplace=True)
|
||||
df['cov'] = df['cov'].astype(int)
|
||||
df = df.iloc[0: max_rows, :]
|
||||
|
||||
# Sort by transcripts with highest isoform diversity
|
||||
df.sort_values('parent gene iso num', inplace=True, ascending=False)
|
||||
|
||||
section.table(df, index=False)
|
||||
|
||||
@ -664,7 +649,7 @@ def transcriptome_summary(report, gffs, sample_ids, denovo=False):
|
||||
Plot transcriptome summaries.
|
||||
|
||||
Some of this data is available via gffcompare output, but the de novo
|
||||
pipeline skips that, so we do it al here.
|
||||
pipeline skips that, so we do it all here.
|
||||
|
||||
We do not report exon number for the denovo assembly yet. This is because
|
||||
in this case, the gff annotation is generated by aligning to the CDS not
|
||||
@ -718,7 +703,7 @@ def transcriptome_summary(report, gffs, sample_ids, denovo=False):
|
||||
|
||||
bar_isos = hist.histogram(
|
||||
[isoforms_per_gene], colors=[Colors.cerulean],
|
||||
title="isoforms per gene")
|
||||
title="isoforms per gene", binwidth=1)
|
||||
bar_isos.xaxis.axis_label = "Num. isoforms"
|
||||
bar_isos.yaxis.axis_label = "Num. genes"
|
||||
|
||||
@ -864,8 +849,8 @@ def main():
|
||||
"--sample_ids", required=True, nargs='+',
|
||||
help="List of sample ids")
|
||||
parser.add_argument(
|
||||
"--transcript_table_cov_thresh", required=False, type=int, default=50,
|
||||
help="Isoforms without this support will be excluded from the table")
|
||||
"--isoform_table_nrows", required=False, type=int, default=5000,
|
||||
help="Maximum rows to display in isoforms table")
|
||||
parser.add_argument(
|
||||
"--cluster_qc_dirs", required=False, type=str, default=None, nargs='*',
|
||||
help="Directory with various cluster quality csvs")
|
||||
@ -912,7 +897,7 @@ def main():
|
||||
report.write(args.report)
|
||||
|
||||
if df_tmaps is not None:
|
||||
transcript_table(report, df_tmaps, args.transcript_table_cov_thresh)
|
||||
transcript_table(report, df_tmaps, args.isoform_table_nrows)
|
||||
|
||||
if args.cluster_qc_dirs is not None:
|
||||
cluster_quality(args.cluster_qc_dirs, report, sample_ids)
|
||||
|
||||
2
main.nf
2
main.nf
@ -309,7 +309,7 @@ process makeReport {
|
||||
--summaries $seq_summaries \
|
||||
--gffcompare_dir $gffcmp_dir \
|
||||
--gff_annotation $gff_annotation \
|
||||
--transcript_table_cov_thresh $params.transcript_table_cov_thresh \
|
||||
--isoform_table_nrows $params.isoform_table_nrows \
|
||||
$OPT_JAFFAL_CSV \
|
||||
$OPT_DENOVO
|
||||
"""
|
||||
|
||||
@ -17,7 +17,7 @@ params {
|
||||
ref_annotation = null
|
||||
threads = 4
|
||||
// Thresholds for viewing isoforms in report table
|
||||
transcript_table_cov_thresh = 50
|
||||
isoform_table_nrows = 5000
|
||||
|
||||
out_dir = "output"
|
||||
sample = null
|
||||
|
||||
@ -84,10 +84,10 @@
|
||||
"description": "Extra options for stringtie transcript assembly.",
|
||||
"default": " --conservative "
|
||||
},
|
||||
"transcript_table_cov_thresh": {
|
||||
"isoform_table_nrows": {
|
||||
"type": "integer",
|
||||
"description": "Minimum coverage for a transcript to appear in the report table",
|
||||
"default": 50
|
||||
"description": "Maximum rows to dispay in the isoform report table",
|
||||
"default": 5000
|
||||
},
|
||||
"denovo": {
|
||||
"type": "boolean",
|
||||
|
||||
Loading…
Reference in New Issue
Block a user